
(FPCore (wj x) :precision binary64 (let* ((t_0 (* wj (exp wj)))) (- wj (/ (- t_0 x) (+ (exp wj) t_0)))))
double code(double wj, double x) {
double t_0 = wj * exp(wj);
return wj - ((t_0 - x) / (exp(wj) + t_0));
}
real(8) function code(wj, x)
real(8), intent (in) :: wj
real(8), intent (in) :: x
real(8) :: t_0
t_0 = wj * exp(wj)
code = wj - ((t_0 - x) / (exp(wj) + t_0))
end function
public static double code(double wj, double x) {
double t_0 = wj * Math.exp(wj);
return wj - ((t_0 - x) / (Math.exp(wj) + t_0));
}
def code(wj, x): t_0 = wj * math.exp(wj) return wj - ((t_0 - x) / (math.exp(wj) + t_0))
function code(wj, x) t_0 = Float64(wj * exp(wj)) return Float64(wj - Float64(Float64(t_0 - x) / Float64(exp(wj) + t_0))) end
function tmp = code(wj, x) t_0 = wj * exp(wj); tmp = wj - ((t_0 - x) / (exp(wj) + t_0)); end
code[wj_, x_] := Block[{t$95$0 = N[(wj * N[Exp[wj], $MachinePrecision]), $MachinePrecision]}, N[(wj - N[(N[(t$95$0 - x), $MachinePrecision] / N[(N[Exp[wj], $MachinePrecision] + t$95$0), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]]
\begin{array}{l}
\\
\begin{array}{l}
t_0 := wj \cdot e^{wj}\\
wj - \frac{t_0 - x}{e^{wj} + t_0}
\end{array}
\end{array}
Sampling outcomes in binary64 precision:
Herbie found 5 alternatives:
| Alternative | Accuracy | Speedup |
|---|
(FPCore (wj x) :precision binary64 (let* ((t_0 (* wj (exp wj)))) (- wj (/ (- t_0 x) (+ (exp wj) t_0)))))
double code(double wj, double x) {
double t_0 = wj * exp(wj);
return wj - ((t_0 - x) / (exp(wj) + t_0));
}
real(8) function code(wj, x)
real(8), intent (in) :: wj
real(8), intent (in) :: x
real(8) :: t_0
t_0 = wj * exp(wj)
code = wj - ((t_0 - x) / (exp(wj) + t_0))
end function
public static double code(double wj, double x) {
double t_0 = wj * Math.exp(wj);
return wj - ((t_0 - x) / (Math.exp(wj) + t_0));
}
def code(wj, x): t_0 = wj * math.exp(wj) return wj - ((t_0 - x) / (math.exp(wj) + t_0))
function code(wj, x) t_0 = Float64(wj * exp(wj)) return Float64(wj - Float64(Float64(t_0 - x) / Float64(exp(wj) + t_0))) end
function tmp = code(wj, x) t_0 = wj * exp(wj); tmp = wj - ((t_0 - x) / (exp(wj) + t_0)); end
code[wj_, x_] := Block[{t$95$0 = N[(wj * N[Exp[wj], $MachinePrecision]), $MachinePrecision]}, N[(wj - N[(N[(t$95$0 - x), $MachinePrecision] / N[(N[Exp[wj], $MachinePrecision] + t$95$0), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]]
\begin{array}{l}
\\
\begin{array}{l}
t_0 := wj \cdot e^{wj}\\
wj - \frac{t_0 - x}{e^{wj} + t_0}
\end{array}
\end{array}
(FPCore (wj x)
:precision binary64
(let* ((t_0 (* wj (exp wj))))
(if (<= (- wj (/ (- t_0 x) (+ (exp wj) t_0))) 2e-12)
(- (+ (* wj wj) (+ x (* -2.0 (* wj x)))) (* wj (* wj wj)))
(+ wj (/ (- (/ x (exp wj)) wj) (+ wj 1.0))))))
double code(double wj, double x) {
double t_0 = wj * exp(wj);
double tmp;
if ((wj - ((t_0 - x) / (exp(wj) + t_0))) <= 2e-12) {
tmp = ((wj * wj) + (x + (-2.0 * (wj * x)))) - (wj * (wj * wj));
} else {
tmp = wj + (((x / exp(wj)) - wj) / (wj + 1.0));
}
return tmp;
}
real(8) function code(wj, x)
real(8), intent (in) :: wj
real(8), intent (in) :: x
real(8) :: t_0
real(8) :: tmp
t_0 = wj * exp(wj)
if ((wj - ((t_0 - x) / (exp(wj) + t_0))) <= 2d-12) then
tmp = ((wj * wj) + (x + ((-2.0d0) * (wj * x)))) - (wj * (wj * wj))
else
tmp = wj + (((x / exp(wj)) - wj) / (wj + 1.0d0))
end if
code = tmp
end function
public static double code(double wj, double x) {
double t_0 = wj * Math.exp(wj);
double tmp;
if ((wj - ((t_0 - x) / (Math.exp(wj) + t_0))) <= 2e-12) {
tmp = ((wj * wj) + (x + (-2.0 * (wj * x)))) - (wj * (wj * wj));
} else {
tmp = wj + (((x / Math.exp(wj)) - wj) / (wj + 1.0));
}
return tmp;
}
def code(wj, x): t_0 = wj * math.exp(wj) tmp = 0 if (wj - ((t_0 - x) / (math.exp(wj) + t_0))) <= 2e-12: tmp = ((wj * wj) + (x + (-2.0 * (wj * x)))) - (wj * (wj * wj)) else: tmp = wj + (((x / math.exp(wj)) - wj) / (wj + 1.0)) return tmp
function code(wj, x) t_0 = Float64(wj * exp(wj)) tmp = 0.0 if (Float64(wj - Float64(Float64(t_0 - x) / Float64(exp(wj) + t_0))) <= 2e-12) tmp = Float64(Float64(Float64(wj * wj) + Float64(x + Float64(-2.0 * Float64(wj * x)))) - Float64(wj * Float64(wj * wj))); else tmp = Float64(wj + Float64(Float64(Float64(x / exp(wj)) - wj) / Float64(wj + 1.0))); end return tmp end
function tmp_2 = code(wj, x) t_0 = wj * exp(wj); tmp = 0.0; if ((wj - ((t_0 - x) / (exp(wj) + t_0))) <= 2e-12) tmp = ((wj * wj) + (x + (-2.0 * (wj * x)))) - (wj * (wj * wj)); else tmp = wj + (((x / exp(wj)) - wj) / (wj + 1.0)); end tmp_2 = tmp; end
code[wj_, x_] := Block[{t$95$0 = N[(wj * N[Exp[wj], $MachinePrecision]), $MachinePrecision]}, If[LessEqual[N[(wj - N[(N[(t$95$0 - x), $MachinePrecision] / N[(N[Exp[wj], $MachinePrecision] + t$95$0), $MachinePrecision]), $MachinePrecision]), $MachinePrecision], 2e-12], N[(N[(N[(wj * wj), $MachinePrecision] + N[(x + N[(-2.0 * N[(wj * x), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision] - N[(wj * N[(wj * wj), $MachinePrecision]), $MachinePrecision]), $MachinePrecision], N[(wj + N[(N[(N[(x / N[Exp[wj], $MachinePrecision]), $MachinePrecision] - wj), $MachinePrecision] / N[(wj + 1.0), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]]]
\begin{array}{l}
\\
\begin{array}{l}
t_0 := wj \cdot e^{wj}\\
\mathbf{if}\;wj - \frac{t_0 - x}{e^{wj} + t_0} \leq 2 \cdot 10^{-12}:\\
\;\;\;\;\left(wj \cdot wj + \left(x + -2 \cdot \left(wj \cdot x\right)\right)\right) - wj \cdot \left(wj \cdot wj\right)\\
\mathbf{else}:\\
\;\;\;\;wj + \frac{\frac{x}{e^{wj}} - wj}{wj + 1}\\
\end{array}
\end{array}
if (-.f64 wj (/.f64 (-.f64 (*.f64 wj (exp.f64 wj)) x) (+.f64 (exp.f64 wj) (*.f64 wj (exp.f64 wj))))) < 1.99999999999999996e-12Initial program 67.6%
sub-neg67.6%
div-sub67.6%
sub-neg67.6%
+-commutative67.6%
distribute-neg-in67.6%
remove-double-neg67.6%
sub-neg67.6%
div-sub67.6%
distribute-rgt1-in68.2%
associate-/l/68.1%
Simplified68.1%
Taylor expanded in wj around 0 98.2%
Taylor expanded in x around 0 98.2%
unpow298.2%
Simplified98.2%
Taylor expanded in x around 0 98.2%
unpow398.2%
Applied egg-rr98.2%
if 1.99999999999999996e-12 < (-.f64 wj (/.f64 (-.f64 (*.f64 wj (exp.f64 wj)) x) (+.f64 (exp.f64 wj) (*.f64 wj (exp.f64 wj))))) Initial program 92.8%
sub-neg92.8%
div-sub92.8%
sub-neg92.8%
+-commutative92.8%
distribute-neg-in92.8%
remove-double-neg92.8%
sub-neg92.8%
div-sub92.8%
distribute-rgt1-in94.2%
associate-/l/94.2%
Simplified99.9%
Final simplification98.7%
(FPCore (wj x) :precision binary64 (if (<= wj 0.0014) (- (+ (* wj wj) (+ x (* -2.0 (* wj x)))) (* wj (* wj wj))) (- wj (/ wj (+ wj 1.0)))))
double code(double wj, double x) {
double tmp;
if (wj <= 0.0014) {
tmp = ((wj * wj) + (x + (-2.0 * (wj * x)))) - (wj * (wj * wj));
} else {
tmp = wj - (wj / (wj + 1.0));
}
return tmp;
}
real(8) function code(wj, x)
real(8), intent (in) :: wj
real(8), intent (in) :: x
real(8) :: tmp
if (wj <= 0.0014d0) then
tmp = ((wj * wj) + (x + ((-2.0d0) * (wj * x)))) - (wj * (wj * wj))
else
tmp = wj - (wj / (wj + 1.0d0))
end if
code = tmp
end function
public static double code(double wj, double x) {
double tmp;
if (wj <= 0.0014) {
tmp = ((wj * wj) + (x + (-2.0 * (wj * x)))) - (wj * (wj * wj));
} else {
tmp = wj - (wj / (wj + 1.0));
}
return tmp;
}
def code(wj, x): tmp = 0 if wj <= 0.0014: tmp = ((wj * wj) + (x + (-2.0 * (wj * x)))) - (wj * (wj * wj)) else: tmp = wj - (wj / (wj + 1.0)) return tmp
function code(wj, x) tmp = 0.0 if (wj <= 0.0014) tmp = Float64(Float64(Float64(wj * wj) + Float64(x + Float64(-2.0 * Float64(wj * x)))) - Float64(wj * Float64(wj * wj))); else tmp = Float64(wj - Float64(wj / Float64(wj + 1.0))); end return tmp end
function tmp_2 = code(wj, x) tmp = 0.0; if (wj <= 0.0014) tmp = ((wj * wj) + (x + (-2.0 * (wj * x)))) - (wj * (wj * wj)); else tmp = wj - (wj / (wj + 1.0)); end tmp_2 = tmp; end
code[wj_, x_] := If[LessEqual[wj, 0.0014], N[(N[(N[(wj * wj), $MachinePrecision] + N[(x + N[(-2.0 * N[(wj * x), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision] - N[(wj * N[(wj * wj), $MachinePrecision]), $MachinePrecision]), $MachinePrecision], N[(wj - N[(wj / N[(wj + 1.0), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]]
\begin{array}{l}
\\
\begin{array}{l}
\mathbf{if}\;wj \leq 0.0014:\\
\;\;\;\;\left(wj \cdot wj + \left(x + -2 \cdot \left(wj \cdot x\right)\right)\right) - wj \cdot \left(wj \cdot wj\right)\\
\mathbf{else}:\\
\;\;\;\;wj - \frac{wj}{wj + 1}\\
\end{array}
\end{array}
if wj < 0.00139999999999999999Initial program 75.6%
sub-neg75.6%
div-sub75.6%
sub-neg75.6%
+-commutative75.6%
distribute-neg-in75.6%
remove-double-neg75.6%
sub-neg75.6%
div-sub75.6%
distribute-rgt1-in76.4%
associate-/l/76.4%
Simplified76.4%
Taylor expanded in wj around 0 98.2%
Taylor expanded in x around 0 97.9%
unpow297.9%
Simplified97.9%
Taylor expanded in x around 0 97.9%
unpow397.9%
Applied egg-rr97.9%
if 0.00139999999999999999 < wj Initial program 20.0%
sub-neg20.0%
div-sub20.0%
sub-neg20.0%
+-commutative20.0%
distribute-neg-in20.0%
remove-double-neg20.0%
sub-neg20.0%
div-sub20.0%
distribute-rgt1-in20.0%
associate-/l/20.0%
Simplified100.0%
Taylor expanded in x around 0 100.0%
+-commutative100.0%
Simplified100.0%
Final simplification98.0%
(FPCore (wj x) :precision binary64 (+ x (* -2.0 (* wj x))))
double code(double wj, double x) {
return x + (-2.0 * (wj * x));
}
real(8) function code(wj, x)
real(8), intent (in) :: wj
real(8), intent (in) :: x
code = x + ((-2.0d0) * (wj * x))
end function
public static double code(double wj, double x) {
return x + (-2.0 * (wj * x));
}
def code(wj, x): return x + (-2.0 * (wj * x))
function code(wj, x) return Float64(x + Float64(-2.0 * Float64(wj * x))) end
function tmp = code(wj, x) tmp = x + (-2.0 * (wj * x)); end
code[wj_, x_] := N[(x + N[(-2.0 * N[(wj * x), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
x + -2 \cdot \left(wj \cdot x\right)
\end{array}
Initial program 74.5%
sub-neg74.5%
div-sub74.5%
sub-neg74.5%
+-commutative74.5%
distribute-neg-in74.5%
remove-double-neg74.5%
sub-neg74.5%
div-sub74.5%
distribute-rgt1-in75.3%
associate-/l/75.3%
Simplified76.8%
Taylor expanded in wj around 0 83.3%
Final simplification83.3%
(FPCore (wj x) :precision binary64 wj)
double code(double wj, double x) {
return wj;
}
real(8) function code(wj, x)
real(8), intent (in) :: wj
real(8), intent (in) :: x
code = wj
end function
public static double code(double wj, double x) {
return wj;
}
def code(wj, x): return wj
function code(wj, x) return wj end
function tmp = code(wj, x) tmp = wj; end
code[wj_, x_] := wj
\begin{array}{l}
\\
wj
\end{array}
Initial program 74.5%
sub-neg74.5%
div-sub74.5%
sub-neg74.5%
+-commutative74.5%
distribute-neg-in74.5%
remove-double-neg74.5%
sub-neg74.5%
div-sub74.5%
distribute-rgt1-in75.3%
associate-/l/75.3%
Simplified76.8%
Taylor expanded in wj around inf 4.9%
Final simplification4.9%
(FPCore (wj x) :precision binary64 x)
double code(double wj, double x) {
return x;
}
real(8) function code(wj, x)
real(8), intent (in) :: wj
real(8), intent (in) :: x
code = x
end function
public static double code(double wj, double x) {
return x;
}
def code(wj, x): return x
function code(wj, x) return x end
function tmp = code(wj, x) tmp = x; end
code[wj_, x_] := x
\begin{array}{l}
\\
x
\end{array}
Initial program 74.5%
sub-neg74.5%
div-sub74.5%
sub-neg74.5%
+-commutative74.5%
distribute-neg-in74.5%
remove-double-neg74.5%
sub-neg74.5%
div-sub74.5%
distribute-rgt1-in75.3%
associate-/l/75.3%
Simplified76.8%
Taylor expanded in wj around 0 82.2%
Final simplification82.2%
(FPCore (wj x) :precision binary64 (- wj (- (/ wj (+ wj 1.0)) (/ x (+ (exp wj) (* wj (exp wj)))))))
double code(double wj, double x) {
return wj - ((wj / (wj + 1.0)) - (x / (exp(wj) + (wj * exp(wj)))));
}
real(8) function code(wj, x)
real(8), intent (in) :: wj
real(8), intent (in) :: x
code = wj - ((wj / (wj + 1.0d0)) - (x / (exp(wj) + (wj * exp(wj)))))
end function
public static double code(double wj, double x) {
return wj - ((wj / (wj + 1.0)) - (x / (Math.exp(wj) + (wj * Math.exp(wj)))));
}
def code(wj, x): return wj - ((wj / (wj + 1.0)) - (x / (math.exp(wj) + (wj * math.exp(wj)))))
function code(wj, x) return Float64(wj - Float64(Float64(wj / Float64(wj + 1.0)) - Float64(x / Float64(exp(wj) + Float64(wj * exp(wj)))))) end
function tmp = code(wj, x) tmp = wj - ((wj / (wj + 1.0)) - (x / (exp(wj) + (wj * exp(wj))))); end
code[wj_, x_] := N[(wj - N[(N[(wj / N[(wj + 1.0), $MachinePrecision]), $MachinePrecision] - N[(x / N[(N[Exp[wj], $MachinePrecision] + N[(wj * N[Exp[wj], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
wj - \left(\frac{wj}{wj + 1} - \frac{x}{e^{wj} + wj \cdot e^{wj}}\right)
\end{array}
herbie shell --seed 2023228
(FPCore (wj x)
:name "Jmat.Real.lambertw, newton loop step"
:precision binary64
:herbie-target
(- wj (- (/ wj (+ wj 1.0)) (/ x (+ (exp wj) (* wj (exp wj))))))
(- wj (/ (- (* wj (exp wj)) x) (+ (exp wj) (* wj (exp wj))))))