
(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 7 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
(if (<= wj -5.2e-6)
(- wj (/ (- wj (/ x (exp wj))) (+ wj 1.0)))
(+
x
(*
wj
(-
(* wj (+ (* x (+ (* wj -2.6666666666666665) 2.5)) (- 1.0 wj)))
(* x 2.0))))))
double code(double wj, double x) {
double tmp;
if (wj <= -5.2e-6) {
tmp = wj - ((wj - (x / exp(wj))) / (wj + 1.0));
} else {
tmp = x + (wj * ((wj * ((x * ((wj * -2.6666666666666665) + 2.5)) + (1.0 - wj))) - (x * 2.0)));
}
return tmp;
}
real(8) function code(wj, x)
real(8), intent (in) :: wj
real(8), intent (in) :: x
real(8) :: tmp
if (wj <= (-5.2d-6)) then
tmp = wj - ((wj - (x / exp(wj))) / (wj + 1.0d0))
else
tmp = x + (wj * ((wj * ((x * ((wj * (-2.6666666666666665d0)) + 2.5d0)) + (1.0d0 - wj))) - (x * 2.0d0)))
end if
code = tmp
end function
public static double code(double wj, double x) {
double tmp;
if (wj <= -5.2e-6) {
tmp = wj - ((wj - (x / Math.exp(wj))) / (wj + 1.0));
} else {
tmp = x + (wj * ((wj * ((x * ((wj * -2.6666666666666665) + 2.5)) + (1.0 - wj))) - (x * 2.0)));
}
return tmp;
}
def code(wj, x): tmp = 0 if wj <= -5.2e-6: tmp = wj - ((wj - (x / math.exp(wj))) / (wj + 1.0)) else: tmp = x + (wj * ((wj * ((x * ((wj * -2.6666666666666665) + 2.5)) + (1.0 - wj))) - (x * 2.0))) return tmp
function code(wj, x) tmp = 0.0 if (wj <= -5.2e-6) tmp = Float64(wj - Float64(Float64(wj - Float64(x / exp(wj))) / Float64(wj + 1.0))); else tmp = Float64(x + Float64(wj * Float64(Float64(wj * Float64(Float64(x * Float64(Float64(wj * -2.6666666666666665) + 2.5)) + Float64(1.0 - wj))) - Float64(x * 2.0)))); end return tmp end
function tmp_2 = code(wj, x) tmp = 0.0; if (wj <= -5.2e-6) tmp = wj - ((wj - (x / exp(wj))) / (wj + 1.0)); else tmp = x + (wj * ((wj * ((x * ((wj * -2.6666666666666665) + 2.5)) + (1.0 - wj))) - (x * 2.0))); end tmp_2 = tmp; end
code[wj_, x_] := If[LessEqual[wj, -5.2e-6], N[(wj - N[(N[(wj - N[(x / N[Exp[wj], $MachinePrecision]), $MachinePrecision]), $MachinePrecision] / N[(wj + 1.0), $MachinePrecision]), $MachinePrecision]), $MachinePrecision], N[(x + N[(wj * N[(N[(wj * N[(N[(x * N[(N[(wj * -2.6666666666666665), $MachinePrecision] + 2.5), $MachinePrecision]), $MachinePrecision] + N[(1.0 - wj), $MachinePrecision]), $MachinePrecision]), $MachinePrecision] - N[(x * 2.0), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]]
\begin{array}{l}
\\
\begin{array}{l}
\mathbf{if}\;wj \leq -5.2 \cdot 10^{-6}:\\
\;\;\;\;wj - \frac{wj - \frac{x}{e^{wj}}}{wj + 1}\\
\mathbf{else}:\\
\;\;\;\;x + wj \cdot \left(wj \cdot \left(x \cdot \left(wj \cdot -2.6666666666666665 + 2.5\right) + \left(1 - wj\right)\right) - x \cdot 2\right)\\
\end{array}
\end{array}
if wj < -5.20000000000000019e-6Initial program 49.2%
distribute-rgt1-in99.5%
associate-/l/99.7%
div-sub49.7%
associate-/l*49.7%
*-inverses99.7%
*-rgt-identity99.7%
Simplified99.7%
if -5.20000000000000019e-6 < wj Initial program 76.2%
distribute-rgt1-in76.2%
associate-/l/76.2%
div-sub76.2%
associate-/l*76.2%
*-inverses77.4%
*-rgt-identity77.4%
Simplified77.4%
Taylor expanded in wj around 0 98.8%
Taylor expanded in x around 0 98.8%
distribute-lft-out98.8%
+-commutative98.8%
*-commutative98.8%
mul-1-neg98.8%
sub-neg98.8%
Simplified98.8%
Final simplification98.8%
(FPCore (wj x)
:precision binary64
(+
x
(*
wj
(-
(* wj (+ (* x (+ (* wj -2.6666666666666665) 2.5)) (- 1.0 wj)))
(* x 2.0)))))
double code(double wj, double x) {
return x + (wj * ((wj * ((x * ((wj * -2.6666666666666665) + 2.5)) + (1.0 - wj))) - (x * 2.0)));
}
real(8) function code(wj, x)
real(8), intent (in) :: wj
real(8), intent (in) :: x
code = x + (wj * ((wj * ((x * ((wj * (-2.6666666666666665d0)) + 2.5d0)) + (1.0d0 - wj))) - (x * 2.0d0)))
end function
public static double code(double wj, double x) {
return x + (wj * ((wj * ((x * ((wj * -2.6666666666666665) + 2.5)) + (1.0 - wj))) - (x * 2.0)));
}
def code(wj, x): return x + (wj * ((wj * ((x * ((wj * -2.6666666666666665) + 2.5)) + (1.0 - wj))) - (x * 2.0)))
function code(wj, x) return Float64(x + Float64(wj * Float64(Float64(wj * Float64(Float64(x * Float64(Float64(wj * -2.6666666666666665) + 2.5)) + Float64(1.0 - wj))) - Float64(x * 2.0)))) end
function tmp = code(wj, x) tmp = x + (wj * ((wj * ((x * ((wj * -2.6666666666666665) + 2.5)) + (1.0 - wj))) - (x * 2.0))); end
code[wj_, x_] := N[(x + N[(wj * N[(N[(wj * N[(N[(x * N[(N[(wj * -2.6666666666666665), $MachinePrecision] + 2.5), $MachinePrecision]), $MachinePrecision] + N[(1.0 - wj), $MachinePrecision]), $MachinePrecision]), $MachinePrecision] - N[(x * 2.0), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
x + wj \cdot \left(wj \cdot \left(x \cdot \left(wj \cdot -2.6666666666666665 + 2.5\right) + \left(1 - wj\right)\right) - x \cdot 2\right)
\end{array}
Initial program 75.5%
distribute-rgt1-in76.7%
associate-/l/76.7%
div-sub75.5%
associate-/l*75.5%
*-inverses77.9%
*-rgt-identity77.9%
Simplified77.9%
Taylor expanded in wj around 0 96.7%
Taylor expanded in x around 0 96.7%
distribute-lft-out96.7%
+-commutative96.7%
*-commutative96.7%
mul-1-neg96.7%
sub-neg96.7%
Simplified96.7%
Final simplification96.7%
(FPCore (wj x) :precision binary64 (- x (* wj (+ (* x 2.0) (* wj (+ wj -1.0))))))
double code(double wj, double x) {
return x - (wj * ((x * 2.0) + (wj * (wj + -1.0))));
}
real(8) function code(wj, x)
real(8), intent (in) :: wj
real(8), intent (in) :: x
code = x - (wj * ((x * 2.0d0) + (wj * (wj + (-1.0d0)))))
end function
public static double code(double wj, double x) {
return x - (wj * ((x * 2.0) + (wj * (wj + -1.0))));
}
def code(wj, x): return x - (wj * ((x * 2.0) + (wj * (wj + -1.0))))
function code(wj, x) return Float64(x - Float64(wj * Float64(Float64(x * 2.0) + Float64(wj * Float64(wj + -1.0))))) end
function tmp = code(wj, x) tmp = x - (wj * ((x * 2.0) + (wj * (wj + -1.0)))); end
code[wj_, x_] := N[(x - N[(wj * N[(N[(x * 2.0), $MachinePrecision] + N[(wj * N[(wj + -1.0), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
x - wj \cdot \left(x \cdot 2 + wj \cdot \left(wj + -1\right)\right)
\end{array}
Initial program 75.5%
distribute-rgt1-in76.7%
associate-/l/76.7%
div-sub75.5%
associate-/l*75.5%
*-inverses77.9%
*-rgt-identity77.9%
Simplified77.9%
Taylor expanded in wj around 0 96.7%
Taylor expanded in x around 0 96.7%
distribute-lft-out96.7%
+-commutative96.7%
*-commutative96.7%
mul-1-neg96.7%
sub-neg96.7%
Simplified96.7%
Taylor expanded in x around 0 96.5%
Final simplification96.5%
(FPCore (wj x) :precision binary64 (+ x (* wj wj)))
double code(double wj, double x) {
return x + (wj * wj);
}
real(8) function code(wj, x)
real(8), intent (in) :: wj
real(8), intent (in) :: x
code = x + (wj * wj)
end function
public static double code(double wj, double x) {
return x + (wj * wj);
}
def code(wj, x): return x + (wj * wj)
function code(wj, x) return Float64(x + Float64(wj * wj)) end
function tmp = code(wj, x) tmp = x + (wj * wj); end
code[wj_, x_] := N[(x + N[(wj * wj), $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
x + wj \cdot wj
\end{array}
Initial program 75.5%
distribute-rgt1-in76.7%
associate-/l/76.7%
div-sub75.5%
associate-/l*75.5%
*-inverses77.9%
*-rgt-identity77.9%
Simplified77.9%
Taylor expanded in wj around 0 96.7%
Taylor expanded in wj around 0 96.5%
distribute-rgt-out96.5%
metadata-eval96.5%
Simplified96.5%
Taylor expanded in x around 0 96.0%
(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 75.5%
distribute-rgt1-in76.7%
associate-/l/76.7%
div-sub75.5%
associate-/l*75.5%
*-inverses77.9%
*-rgt-identity77.9%
Simplified77.9%
Taylor expanded in wj around 0 85.8%
(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 75.5%
distribute-rgt1-in76.7%
associate-/l/76.7%
div-sub75.5%
associate-/l*75.5%
*-inverses77.9%
*-rgt-identity77.9%
Simplified77.9%
Taylor expanded in wj around inf 4.5%
(FPCore (wj x) :precision binary64 -1.0)
double code(double wj, double x) {
return -1.0;
}
real(8) function code(wj, x)
real(8), intent (in) :: wj
real(8), intent (in) :: x
code = -1.0d0
end function
public static double code(double wj, double x) {
return -1.0;
}
def code(wj, x): return -1.0
function code(wj, x) return -1.0 end
function tmp = code(wj, x) tmp = -1.0; end
code[wj_, x_] := -1.0
\begin{array}{l}
\\
-1
\end{array}
Initial program 75.5%
distribute-rgt1-in76.7%
associate-/l/76.7%
div-sub75.5%
associate-/l*75.5%
*-inverses77.9%
*-rgt-identity77.9%
Simplified77.9%
Taylor expanded in wj around inf 4.4%
Taylor expanded in wj around 0 3.3%
(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 2024146
(FPCore (wj x)
:name "Jmat.Real.lambertw, newton loop step"
:precision binary64
:alt
(! :herbie-platform default (let ((ew (exp wj))) (- wj (- (/ wj (+ wj 1)) (/ x (+ ew (* wj ew)))))))
(- wj (/ (- (* wj (exp wj)) x) (+ (exp wj) (* wj (exp wj))))))