
(FPCore (x y) :precision binary64 (- 1.0 (log (- 1.0 (/ (- x y) (- 1.0 y))))))
double code(double x, double y) {
return 1.0 - log((1.0 - ((x - y) / (1.0 - y))));
}
real(8) function code(x, y)
real(8), intent (in) :: x
real(8), intent (in) :: y
code = 1.0d0 - log((1.0d0 - ((x - y) / (1.0d0 - y))))
end function
public static double code(double x, double y) {
return 1.0 - Math.log((1.0 - ((x - y) / (1.0 - y))));
}
def code(x, y): return 1.0 - math.log((1.0 - ((x - y) / (1.0 - y))))
function code(x, y) return Float64(1.0 - log(Float64(1.0 - Float64(Float64(x - y) / Float64(1.0 - y))))) end
function tmp = code(x, y) tmp = 1.0 - log((1.0 - ((x - y) / (1.0 - y)))); end
code[x_, y_] := N[(1.0 - N[Log[N[(1.0 - N[(N[(x - y), $MachinePrecision] / N[(1.0 - y), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
1 - \log \left(1 - \frac{x - y}{1 - y}\right)
\end{array}
Sampling outcomes in binary64 precision:
Herbie found 7 alternatives:
| Alternative | Accuracy | Speedup |
|---|
(FPCore (x y) :precision binary64 (- 1.0 (log (- 1.0 (/ (- x y) (- 1.0 y))))))
double code(double x, double y) {
return 1.0 - log((1.0 - ((x - y) / (1.0 - y))));
}
real(8) function code(x, y)
real(8), intent (in) :: x
real(8), intent (in) :: y
code = 1.0d0 - log((1.0d0 - ((x - y) / (1.0d0 - y))))
end function
public static double code(double x, double y) {
return 1.0 - Math.log((1.0 - ((x - y) / (1.0 - y))));
}
def code(x, y): return 1.0 - math.log((1.0 - ((x - y) / (1.0 - y))))
function code(x, y) return Float64(1.0 - log(Float64(1.0 - Float64(Float64(x - y) / Float64(1.0 - y))))) end
function tmp = code(x, y) tmp = 1.0 - log((1.0 - ((x - y) / (1.0 - y)))); end
code[x_, y_] := N[(1.0 - N[Log[N[(1.0 - N[(N[(x - y), $MachinePrecision] / N[(1.0 - y), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
1 - \log \left(1 - \frac{x - y}{1 - y}\right)
\end{array}
(FPCore (x y) :precision binary64 (if (<= (/ (- y x) (+ y -1.0)) 0.99998) (- 1.0 (log1p (/ (- x y) (+ y -1.0)))) (+ 1.0 (- (/ -1.0 y) (log (* (/ 1.0 y) (+ x -1.0)))))))
double code(double x, double y) {
double tmp;
if (((y - x) / (y + -1.0)) <= 0.99998) {
tmp = 1.0 - log1p(((x - y) / (y + -1.0)));
} else {
tmp = 1.0 + ((-1.0 / y) - log(((1.0 / y) * (x + -1.0))));
}
return tmp;
}
public static double code(double x, double y) {
double tmp;
if (((y - x) / (y + -1.0)) <= 0.99998) {
tmp = 1.0 - Math.log1p(((x - y) / (y + -1.0)));
} else {
tmp = 1.0 + ((-1.0 / y) - Math.log(((1.0 / y) * (x + -1.0))));
}
return tmp;
}
def code(x, y): tmp = 0 if ((y - x) / (y + -1.0)) <= 0.99998: tmp = 1.0 - math.log1p(((x - y) / (y + -1.0))) else: tmp = 1.0 + ((-1.0 / y) - math.log(((1.0 / y) * (x + -1.0)))) return tmp
function code(x, y) tmp = 0.0 if (Float64(Float64(y - x) / Float64(y + -1.0)) <= 0.99998) tmp = Float64(1.0 - log1p(Float64(Float64(x - y) / Float64(y + -1.0)))); else tmp = Float64(1.0 + Float64(Float64(-1.0 / y) - log(Float64(Float64(1.0 / y) * Float64(x + -1.0))))); end return tmp end
code[x_, y_] := If[LessEqual[N[(N[(y - x), $MachinePrecision] / N[(y + -1.0), $MachinePrecision]), $MachinePrecision], 0.99998], N[(1.0 - N[Log[1 + N[(N[(x - y), $MachinePrecision] / N[(y + -1.0), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision], N[(1.0 + N[(N[(-1.0 / y), $MachinePrecision] - N[Log[N[(N[(1.0 / y), $MachinePrecision] * N[(x + -1.0), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]]
\begin{array}{l}
\\
\begin{array}{l}
\mathbf{if}\;\frac{y - x}{y + -1} \leq 0.99998:\\
\;\;\;\;1 - \mathsf{log1p}\left(\frac{x - y}{y + -1}\right)\\
\mathbf{else}:\\
\;\;\;\;1 + \left(\frac{-1}{y} - \log \left(\frac{1}{y} \cdot \left(x + -1\right)\right)\right)\\
\end{array}
\end{array}
if (/.f64 (-.f64 x y) (-.f64 #s(literal 1 binary64) y)) < 0.99997999999999998Initial program 99.8%
sub-neg99.8%
log1p-define99.9%
distribute-neg-frac299.9%
neg-sub099.9%
associate--r-99.9%
metadata-eval99.9%
+-commutative99.9%
Simplified99.9%
if 0.99997999999999998 < (/.f64 (-.f64 x y) (-.f64 #s(literal 1 binary64) y)) Initial program 4.0%
sub-neg4.0%
log1p-define4.0%
distribute-neg-frac24.0%
neg-sub04.0%
associate--r-4.0%
metadata-eval4.0%
+-commutative4.0%
Simplified4.0%
Taylor expanded in y around -inf 89.3%
Simplified89.3%
*-un-lft-identity89.3%
associate-/l/89.3%
log1p-undefine89.3%
sub-neg89.3%
sum-log100.0%
Applied egg-rr100.0%
*-lft-identity100.0%
*-lft-identity100.0%
times-frac100.0%
associate-*r/100.0%
*-commutative100.0%
rgt-mult-inverse100.0%
*-commutative100.0%
Simplified100.0%
Final simplification99.9%
(FPCore (x y) :precision binary64 (if (<= (/ (- y x) (+ y -1.0)) 0.99998) (- 1.0 (log1p (/ (- x y) (+ y -1.0)))) (- 1.0 (+ (log (/ -1.0 y)) (- (/ 1.0 y) x)))))
double code(double x, double y) {
double tmp;
if (((y - x) / (y + -1.0)) <= 0.99998) {
tmp = 1.0 - log1p(((x - y) / (y + -1.0)));
} else {
tmp = 1.0 - (log((-1.0 / y)) + ((1.0 / y) - x));
}
return tmp;
}
public static double code(double x, double y) {
double tmp;
if (((y - x) / (y + -1.0)) <= 0.99998) {
tmp = 1.0 - Math.log1p(((x - y) / (y + -1.0)));
} else {
tmp = 1.0 - (Math.log((-1.0 / y)) + ((1.0 / y) - x));
}
return tmp;
}
def code(x, y): tmp = 0 if ((y - x) / (y + -1.0)) <= 0.99998: tmp = 1.0 - math.log1p(((x - y) / (y + -1.0))) else: tmp = 1.0 - (math.log((-1.0 / y)) + ((1.0 / y) - x)) return tmp
function code(x, y) tmp = 0.0 if (Float64(Float64(y - x) / Float64(y + -1.0)) <= 0.99998) tmp = Float64(1.0 - log1p(Float64(Float64(x - y) / Float64(y + -1.0)))); else tmp = Float64(1.0 - Float64(log(Float64(-1.0 / y)) + Float64(Float64(1.0 / y) - x))); end return tmp end
code[x_, y_] := If[LessEqual[N[(N[(y - x), $MachinePrecision] / N[(y + -1.0), $MachinePrecision]), $MachinePrecision], 0.99998], N[(1.0 - N[Log[1 + N[(N[(x - y), $MachinePrecision] / N[(y + -1.0), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision], N[(1.0 - N[(N[Log[N[(-1.0 / y), $MachinePrecision]], $MachinePrecision] + N[(N[(1.0 / y), $MachinePrecision] - x), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]]
\begin{array}{l}
\\
\begin{array}{l}
\mathbf{if}\;\frac{y - x}{y + -1} \leq 0.99998:\\
\;\;\;\;1 - \mathsf{log1p}\left(\frac{x - y}{y + -1}\right)\\
\mathbf{else}:\\
\;\;\;\;1 - \left(\log \left(\frac{-1}{y}\right) + \left(\frac{1}{y} - x\right)\right)\\
\end{array}
\end{array}
if (/.f64 (-.f64 x y) (-.f64 #s(literal 1 binary64) y)) < 0.99997999999999998Initial program 99.8%
sub-neg99.8%
log1p-define99.9%
distribute-neg-frac299.9%
neg-sub099.9%
associate--r-99.9%
metadata-eval99.9%
+-commutative99.9%
Simplified99.9%
if 0.99997999999999998 < (/.f64 (-.f64 x y) (-.f64 #s(literal 1 binary64) y)) Initial program 4.0%
sub-neg4.0%
log1p-define4.0%
distribute-neg-frac24.0%
neg-sub04.0%
associate--r-4.0%
metadata-eval4.0%
+-commutative4.0%
Simplified4.0%
Taylor expanded in x around 0 1.4%
Taylor expanded in y around -inf 70.7%
+-commutative70.7%
neg-mul-170.7%
unsub-neg70.7%
Simplified70.7%
Final simplification91.2%
(FPCore (x y) :precision binary64 (if (<= (/ (- y x) (+ y -1.0)) 0.99998) (- 1.0 (log1p (/ (- x y) (+ y -1.0)))) (+ 1.0 (- x (log (/ -1.0 y))))))
double code(double x, double y) {
double tmp;
if (((y - x) / (y + -1.0)) <= 0.99998) {
tmp = 1.0 - log1p(((x - y) / (y + -1.0)));
} else {
tmp = 1.0 + (x - log((-1.0 / y)));
}
return tmp;
}
public static double code(double x, double y) {
double tmp;
if (((y - x) / (y + -1.0)) <= 0.99998) {
tmp = 1.0 - Math.log1p(((x - y) / (y + -1.0)));
} else {
tmp = 1.0 + (x - Math.log((-1.0 / y)));
}
return tmp;
}
def code(x, y): tmp = 0 if ((y - x) / (y + -1.0)) <= 0.99998: tmp = 1.0 - math.log1p(((x - y) / (y + -1.0))) else: tmp = 1.0 + (x - math.log((-1.0 / y))) return tmp
function code(x, y) tmp = 0.0 if (Float64(Float64(y - x) / Float64(y + -1.0)) <= 0.99998) tmp = Float64(1.0 - log1p(Float64(Float64(x - y) / Float64(y + -1.0)))); else tmp = Float64(1.0 + Float64(x - log(Float64(-1.0 / y)))); end return tmp end
code[x_, y_] := If[LessEqual[N[(N[(y - x), $MachinePrecision] / N[(y + -1.0), $MachinePrecision]), $MachinePrecision], 0.99998], N[(1.0 - N[Log[1 + N[(N[(x - y), $MachinePrecision] / N[(y + -1.0), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision], N[(1.0 + N[(x - N[Log[N[(-1.0 / y), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]]
\begin{array}{l}
\\
\begin{array}{l}
\mathbf{if}\;\frac{y - x}{y + -1} \leq 0.99998:\\
\;\;\;\;1 - \mathsf{log1p}\left(\frac{x - y}{y + -1}\right)\\
\mathbf{else}:\\
\;\;\;\;1 + \left(x - \log \left(\frac{-1}{y}\right)\right)\\
\end{array}
\end{array}
if (/.f64 (-.f64 x y) (-.f64 #s(literal 1 binary64) y)) < 0.99997999999999998Initial program 99.8%
sub-neg99.8%
log1p-define99.9%
distribute-neg-frac299.9%
neg-sub099.9%
associate--r-99.9%
metadata-eval99.9%
+-commutative99.9%
Simplified99.9%
if 0.99997999999999998 < (/.f64 (-.f64 x y) (-.f64 #s(literal 1 binary64) y)) Initial program 4.0%
sub-neg4.0%
log1p-define4.0%
distribute-neg-frac24.0%
neg-sub04.0%
associate--r-4.0%
metadata-eval4.0%
+-commutative4.0%
Simplified4.0%
Taylor expanded in x around 0 1.4%
Taylor expanded in y around -inf 70.7%
neg-mul-170.7%
unsub-neg70.7%
Simplified70.7%
Final simplification91.2%
(FPCore (x y) :precision binary64 (if (<= y -6.4e+23) (+ 1.0 (- x (log (/ -1.0 y)))) (- 1.0 (log1p (/ x (+ y -1.0))))))
double code(double x, double y) {
double tmp;
if (y <= -6.4e+23) {
tmp = 1.0 + (x - log((-1.0 / y)));
} else {
tmp = 1.0 - log1p((x / (y + -1.0)));
}
return tmp;
}
public static double code(double x, double y) {
double tmp;
if (y <= -6.4e+23) {
tmp = 1.0 + (x - Math.log((-1.0 / y)));
} else {
tmp = 1.0 - Math.log1p((x / (y + -1.0)));
}
return tmp;
}
def code(x, y): tmp = 0 if y <= -6.4e+23: tmp = 1.0 + (x - math.log((-1.0 / y))) else: tmp = 1.0 - math.log1p((x / (y + -1.0))) return tmp
function code(x, y) tmp = 0.0 if (y <= -6.4e+23) tmp = Float64(1.0 + Float64(x - log(Float64(-1.0 / y)))); else tmp = Float64(1.0 - log1p(Float64(x / Float64(y + -1.0)))); end return tmp end
code[x_, y_] := If[LessEqual[y, -6.4e+23], N[(1.0 + N[(x - N[Log[N[(-1.0 / y), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]), $MachinePrecision], N[(1.0 - N[Log[1 + N[(x / N[(y + -1.0), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]]
\begin{array}{l}
\\
\begin{array}{l}
\mathbf{if}\;y \leq -6.4 \cdot 10^{+23}:\\
\;\;\;\;1 + \left(x - \log \left(\frac{-1}{y}\right)\right)\\
\mathbf{else}:\\
\;\;\;\;1 - \mathsf{log1p}\left(\frac{x}{y + -1}\right)\\
\end{array}
\end{array}
if y < -6.4e23Initial program 12.7%
sub-neg12.7%
log1p-define12.7%
distribute-neg-frac212.7%
neg-sub012.7%
associate--r-12.7%
metadata-eval12.7%
+-commutative12.7%
Simplified12.7%
Taylor expanded in x around 0 0.9%
Taylor expanded in y around -inf 71.6%
neg-mul-171.6%
unsub-neg71.6%
Simplified71.6%
if -6.4e23 < y Initial program 95.2%
sub-neg95.2%
log1p-define95.3%
distribute-neg-frac295.3%
neg-sub095.3%
associate--r-95.3%
metadata-eval95.3%
+-commutative95.3%
Simplified95.3%
Taylor expanded in x around inf 91.6%
Final simplification85.8%
(FPCore (x y) :precision binary64 (- 1.0 (log1p (/ x (+ y -1.0)))))
double code(double x, double y) {
return 1.0 - log1p((x / (y + -1.0)));
}
public static double code(double x, double y) {
return 1.0 - Math.log1p((x / (y + -1.0)));
}
def code(x, y): return 1.0 - math.log1p((x / (y + -1.0)))
function code(x, y) return Float64(1.0 - log1p(Float64(x / Float64(y + -1.0)))) end
code[x_, y_] := N[(1.0 - N[Log[1 + N[(x / N[(y + -1.0), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
1 - \mathsf{log1p}\left(\frac{x}{y + -1}\right)
\end{array}
Initial program 71.4%
sub-neg71.4%
log1p-define71.4%
distribute-neg-frac271.4%
neg-sub071.4%
associate--r-71.4%
metadata-eval71.4%
+-commutative71.4%
Simplified71.4%
Taylor expanded in x around inf 71.2%
Final simplification71.2%
(FPCore (x y) :precision binary64 (- 1.0 (log1p (- x))))
double code(double x, double y) {
return 1.0 - log1p(-x);
}
public static double code(double x, double y) {
return 1.0 - Math.log1p(-x);
}
def code(x, y): return 1.0 - math.log1p(-x)
function code(x, y) return Float64(1.0 - log1p(Float64(-x))) end
code[x_, y_] := N[(1.0 - N[Log[1 + (-x)], $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
1 - \mathsf{log1p}\left(-x\right)
\end{array}
Initial program 71.4%
sub-neg71.4%
log1p-define71.4%
distribute-neg-frac271.4%
neg-sub071.4%
associate--r-71.4%
metadata-eval71.4%
+-commutative71.4%
Simplified71.4%
Taylor expanded in y around 0 61.6%
log1p-define61.6%
mul-1-neg61.6%
Simplified61.6%
(FPCore (x y) :precision binary64 (+ x 1.0))
double code(double x, double y) {
return x + 1.0;
}
real(8) function code(x, y)
real(8), intent (in) :: x
real(8), intent (in) :: y
code = x + 1.0d0
end function
public static double code(double x, double y) {
return x + 1.0;
}
def code(x, y): return x + 1.0
function code(x, y) return Float64(x + 1.0) end
function tmp = code(x, y) tmp = x + 1.0; end
code[x_, y_] := N[(x + 1.0), $MachinePrecision]
\begin{array}{l}
\\
x + 1
\end{array}
Initial program 71.4%
sub-neg71.4%
log1p-define71.4%
distribute-neg-frac271.4%
neg-sub071.4%
associate--r-71.4%
metadata-eval71.4%
+-commutative71.4%
Simplified71.4%
Taylor expanded in x around 0 44.6%
Taylor expanded in y around 0 45.4%
neg-mul-145.4%
Simplified45.4%
Final simplification45.4%
(FPCore (x y)
:precision binary64
(let* ((t_0 (- 1.0 (log (- (/ x (* y y)) (- (/ 1.0 y) (/ x y)))))))
(if (< y -81284752.61947241)
t_0
(if (< y 3.0094271212461764e+25)
(log (/ (exp 1.0) (- 1.0 (/ (- x y) (- 1.0 y)))))
t_0))))
double code(double x, double y) {
double t_0 = 1.0 - log(((x / (y * y)) - ((1.0 / y) - (x / y))));
double tmp;
if (y < -81284752.61947241) {
tmp = t_0;
} else if (y < 3.0094271212461764e+25) {
tmp = log((exp(1.0) / (1.0 - ((x - y) / (1.0 - y)))));
} else {
tmp = t_0;
}
return tmp;
}
real(8) function code(x, y)
real(8), intent (in) :: x
real(8), intent (in) :: y
real(8) :: t_0
real(8) :: tmp
t_0 = 1.0d0 - log(((x / (y * y)) - ((1.0d0 / y) - (x / y))))
if (y < (-81284752.61947241d0)) then
tmp = t_0
else if (y < 3.0094271212461764d+25) then
tmp = log((exp(1.0d0) / (1.0d0 - ((x - y) / (1.0d0 - y)))))
else
tmp = t_0
end if
code = tmp
end function
public static double code(double x, double y) {
double t_0 = 1.0 - Math.log(((x / (y * y)) - ((1.0 / y) - (x / y))));
double tmp;
if (y < -81284752.61947241) {
tmp = t_0;
} else if (y < 3.0094271212461764e+25) {
tmp = Math.log((Math.exp(1.0) / (1.0 - ((x - y) / (1.0 - y)))));
} else {
tmp = t_0;
}
return tmp;
}
def code(x, y): t_0 = 1.0 - math.log(((x / (y * y)) - ((1.0 / y) - (x / y)))) tmp = 0 if y < -81284752.61947241: tmp = t_0 elif y < 3.0094271212461764e+25: tmp = math.log((math.exp(1.0) / (1.0 - ((x - y) / (1.0 - y))))) else: tmp = t_0 return tmp
function code(x, y) t_0 = Float64(1.0 - log(Float64(Float64(x / Float64(y * y)) - Float64(Float64(1.0 / y) - Float64(x / y))))) tmp = 0.0 if (y < -81284752.61947241) tmp = t_0; elseif (y < 3.0094271212461764e+25) tmp = log(Float64(exp(1.0) / Float64(1.0 - Float64(Float64(x - y) / Float64(1.0 - y))))); else tmp = t_0; end return tmp end
function tmp_2 = code(x, y) t_0 = 1.0 - log(((x / (y * y)) - ((1.0 / y) - (x / y)))); tmp = 0.0; if (y < -81284752.61947241) tmp = t_0; elseif (y < 3.0094271212461764e+25) tmp = log((exp(1.0) / (1.0 - ((x - y) / (1.0 - y))))); else tmp = t_0; end tmp_2 = tmp; end
code[x_, y_] := Block[{t$95$0 = N[(1.0 - N[Log[N[(N[(x / N[(y * y), $MachinePrecision]), $MachinePrecision] - N[(N[(1.0 / y), $MachinePrecision] - N[(x / y), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]}, If[Less[y, -81284752.61947241], t$95$0, If[Less[y, 3.0094271212461764e+25], N[Log[N[(N[Exp[1.0], $MachinePrecision] / N[(1.0 - N[(N[(x - y), $MachinePrecision] / N[(1.0 - y), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision], t$95$0]]]
\begin{array}{l}
\\
\begin{array}{l}
t_0 := 1 - \log \left(\frac{x}{y \cdot y} - \left(\frac{1}{y} - \frac{x}{y}\right)\right)\\
\mathbf{if}\;y < -81284752.61947241:\\
\;\;\;\;t\_0\\
\mathbf{elif}\;y < 3.0094271212461764 \cdot 10^{+25}:\\
\;\;\;\;\log \left(\frac{e^{1}}{1 - \frac{x - y}{1 - y}}\right)\\
\mathbf{else}:\\
\;\;\;\;t\_0\\
\end{array}
\end{array}
herbie shell --seed 2024100
(FPCore (x y)
:name "Numeric.SpecFunctions:invIncompleteGamma from math-functions-0.1.5.2, B"
:precision binary64
:alt
(if (< y -81284752.61947241) (- 1.0 (log (- (/ x (* y y)) (- (/ 1.0 y) (/ x y))))) (if (< y 3.0094271212461764e+25) (log (/ (exp 1.0) (- 1.0 (/ (- x y) (- 1.0 y))))) (- 1.0 (log (- (/ x (* y y)) (- (/ 1.0 y) (/ x y)))))))
(- 1.0 (log (- 1.0 (/ (- x y) (- 1.0 y))))))