
(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 5 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 -1.85e+21)
(+ 1.0 (- (log (- y)) (log1p (- x))))
(if (<= y 1.3e+79)
(- 1.0 (log1p (/ (- y x) (- 1.0 y))))
(+ 1.0 (- (log y) (log x))))))
double code(double x, double y) {
double tmp;
if (y <= -1.85e+21) {
tmp = 1.0 + (log(-y) - log1p(-x));
} else if (y <= 1.3e+79) {
tmp = 1.0 - log1p(((y - x) / (1.0 - y)));
} else {
tmp = 1.0 + (log(y) - log(x));
}
return tmp;
}
public static double code(double x, double y) {
double tmp;
if (y <= -1.85e+21) {
tmp = 1.0 + (Math.log(-y) - Math.log1p(-x));
} else if (y <= 1.3e+79) {
tmp = 1.0 - Math.log1p(((y - x) / (1.0 - y)));
} else {
tmp = 1.0 + (Math.log(y) - Math.log(x));
}
return tmp;
}
def code(x, y): tmp = 0 if y <= -1.85e+21: tmp = 1.0 + (math.log(-y) - math.log1p(-x)) elif y <= 1.3e+79: tmp = 1.0 - math.log1p(((y - x) / (1.0 - y))) else: tmp = 1.0 + (math.log(y) - math.log(x)) return tmp
function code(x, y) tmp = 0.0 if (y <= -1.85e+21) tmp = Float64(1.0 + Float64(log(Float64(-y)) - log1p(Float64(-x)))); elseif (y <= 1.3e+79) tmp = Float64(1.0 - log1p(Float64(Float64(y - x) / Float64(1.0 - y)))); else tmp = Float64(1.0 + Float64(log(y) - log(x))); end return tmp end
code[x_, y_] := If[LessEqual[y, -1.85e+21], N[(1.0 + N[(N[Log[(-y)], $MachinePrecision] - N[Log[1 + (-x)], $MachinePrecision]), $MachinePrecision]), $MachinePrecision], If[LessEqual[y, 1.3e+79], N[(1.0 - N[Log[1 + N[(N[(y - x), $MachinePrecision] / N[(1.0 - y), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision], N[(1.0 + N[(N[Log[y], $MachinePrecision] - N[Log[x], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]]]
\begin{array}{l}
\\
\begin{array}{l}
\mathbf{if}\;y \leq -1.85 \cdot 10^{+21}:\\
\;\;\;\;1 + \left(\log \left(-y\right) - \mathsf{log1p}\left(-x\right)\right)\\
\mathbf{elif}\;y \leq 1.3 \cdot 10^{+79}:\\
\;\;\;\;1 - \mathsf{log1p}\left(\frac{y - x}{1 - y}\right)\\
\mathbf{else}:\\
\;\;\;\;1 + \left(\log y - \log x\right)\\
\end{array}
\end{array}
if y < -1.85e21Initial program 16.4%
sub-neg16.4%
log1p-def16.4%
distribute-neg-frac16.4%
sub-neg16.4%
distribute-neg-in16.4%
remove-double-neg16.4%
+-commutative16.4%
sub-neg16.4%
Simplified16.4%
Taylor expanded in y around -inf 99.4%
sub-neg99.4%
metadata-eval99.4%
distribute-lft-in99.4%
metadata-eval99.4%
+-commutative99.4%
log1p-def99.4%
mul-1-neg99.4%
Simplified99.4%
Taylor expanded in y around 0 0.0%
+-commutative0.0%
neg-mul-10.0%
unsub-neg0.0%
associate-+l-0.0%
sub-neg0.0%
mul-1-neg0.0%
log1p-def0.0%
mul-1-neg0.0%
log-div99.4%
metadata-eval99.4%
associate-/l*99.4%
/-rgt-identity99.4%
*-commutative99.4%
neg-mul-199.4%
Simplified99.4%
if -1.85e21 < y < 1.30000000000000007e79Initial program 100.0%
sub-neg100.0%
log1p-def100.0%
distribute-neg-frac100.0%
sub-neg100.0%
distribute-neg-in100.0%
remove-double-neg100.0%
+-commutative100.0%
sub-neg100.0%
Simplified100.0%
if 1.30000000000000007e79 < y Initial program 51.1%
sub-neg51.1%
log1p-def51.1%
distribute-neg-frac51.1%
sub-neg51.1%
distribute-neg-in51.1%
remove-double-neg51.1%
+-commutative51.1%
sub-neg51.1%
Simplified51.1%
Taylor expanded in y around -inf 0.0%
sub-neg0.0%
metadata-eval0.0%
distribute-lft-in0.0%
metadata-eval0.0%
+-commutative0.0%
log1p-def0.0%
mul-1-neg0.0%
Simplified0.0%
Taylor expanded in x around inf 0.0%
associate-+r+0.0%
+-commutative0.0%
+-commutative0.0%
log-div0.0%
associate-+l-0.0%
log-div0.0%
metadata-eval0.0%
associate-/l*0.0%
/-rgt-identity0.0%
log-div98.2%
metadata-eval98.2%
associate-/r*98.2%
*-commutative98.2%
neg-mul-198.2%
neg-mul-198.2%
remove-double-neg98.2%
log-rec98.2%
unsub-neg98.2%
Simplified98.2%
Final simplification99.7%
(FPCore (x y) :precision binary64 (if (<= (/ (- x y) (- 1.0 y)) 1.0) (- 1.0 (log1p (/ (- y x) (- 1.0 y)))) (- 1.0 (log (/ -1.0 y)))))
double code(double x, double y) {
double tmp;
if (((x - y) / (1.0 - y)) <= 1.0) {
tmp = 1.0 - log1p(((y - x) / (1.0 - y)));
} else {
tmp = 1.0 - log((-1.0 / y));
}
return tmp;
}
public static double code(double x, double y) {
double tmp;
if (((x - y) / (1.0 - y)) <= 1.0) {
tmp = 1.0 - Math.log1p(((y - x) / (1.0 - y)));
} else {
tmp = 1.0 - Math.log((-1.0 / y));
}
return tmp;
}
def code(x, y): tmp = 0 if ((x - y) / (1.0 - y)) <= 1.0: tmp = 1.0 - math.log1p(((y - x) / (1.0 - y))) else: tmp = 1.0 - math.log((-1.0 / y)) return tmp
function code(x, y) tmp = 0.0 if (Float64(Float64(x - y) / Float64(1.0 - y)) <= 1.0) tmp = Float64(1.0 - log1p(Float64(Float64(y - x) / Float64(1.0 - y)))); else tmp = Float64(1.0 - log(Float64(-1.0 / y))); end return tmp end
code[x_, y_] := If[LessEqual[N[(N[(x - y), $MachinePrecision] / N[(1.0 - y), $MachinePrecision]), $MachinePrecision], 1.0], N[(1.0 - N[Log[1 + N[(N[(y - x), $MachinePrecision] / N[(1.0 - y), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision], N[(1.0 - N[Log[N[(-1.0 / y), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]]
\begin{array}{l}
\\
\begin{array}{l}
\mathbf{if}\;\frac{x - y}{1 - y} \leq 1:\\
\;\;\;\;1 - \mathsf{log1p}\left(\frac{y - x}{1 - y}\right)\\
\mathbf{else}:\\
\;\;\;\;1 - \log \left(\frac{-1}{y}\right)\\
\end{array}
\end{array}
if (/.f64 (-.f64 x y) (-.f64 1 y)) < 1Initial program 72.5%
sub-neg72.5%
log1p-def72.5%
distribute-neg-frac72.5%
sub-neg72.5%
distribute-neg-in72.5%
remove-double-neg72.5%
+-commutative72.5%
sub-neg72.5%
Simplified72.5%
if 1 < (/.f64 (-.f64 x y) (-.f64 1 y)) Initial program 72.5%
sub-neg72.5%
log1p-def72.5%
distribute-neg-frac72.5%
sub-neg72.5%
distribute-neg-in72.5%
remove-double-neg72.5%
+-commutative72.5%
sub-neg72.5%
Simplified72.5%
Taylor expanded in y around -inf 31.1%
sub-neg31.1%
metadata-eval31.1%
distribute-lft-in31.1%
metadata-eval31.1%
+-commutative31.1%
log1p-def31.1%
mul-1-neg31.1%
Simplified31.1%
Taylor expanded in x around 0 23.0%
Final simplification72.5%
(FPCore (x y) :precision binary64 (if (<= y -1.85) (- 1.0 (log (/ -1.0 y))) (- 1.0 (+ y (log1p (- x))))))
double code(double x, double y) {
double tmp;
if (y <= -1.85) {
tmp = 1.0 - log((-1.0 / y));
} else {
tmp = 1.0 - (y + log1p(-x));
}
return tmp;
}
public static double code(double x, double y) {
double tmp;
if (y <= -1.85) {
tmp = 1.0 - Math.log((-1.0 / y));
} else {
tmp = 1.0 - (y + Math.log1p(-x));
}
return tmp;
}
def code(x, y): tmp = 0 if y <= -1.85: tmp = 1.0 - math.log((-1.0 / y)) else: tmp = 1.0 - (y + math.log1p(-x)) return tmp
function code(x, y) tmp = 0.0 if (y <= -1.85) tmp = Float64(1.0 - log(Float64(-1.0 / y))); else tmp = Float64(1.0 - Float64(y + log1p(Float64(-x)))); end return tmp end
code[x_, y_] := If[LessEqual[y, -1.85], N[(1.0 - N[Log[N[(-1.0 / y), $MachinePrecision]], $MachinePrecision]), $MachinePrecision], N[(1.0 - N[(y + N[Log[1 + (-x)], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]]
\begin{array}{l}
\\
\begin{array}{l}
\mathbf{if}\;y \leq -1.85:\\
\;\;\;\;1 - \log \left(\frac{-1}{y}\right)\\
\mathbf{else}:\\
\;\;\;\;1 - \left(y + \mathsf{log1p}\left(-x\right)\right)\\
\end{array}
\end{array}
if y < -1.8500000000000001Initial program 19.7%
sub-neg19.7%
log1p-def19.7%
distribute-neg-frac19.7%
sub-neg19.7%
distribute-neg-in19.7%
remove-double-neg19.7%
+-commutative19.7%
sub-neg19.7%
Simplified19.7%
Taylor expanded in y around -inf 98.3%
sub-neg98.3%
metadata-eval98.3%
distribute-lft-in98.3%
metadata-eval98.3%
+-commutative98.3%
log1p-def98.3%
mul-1-neg98.3%
Simplified98.3%
Taylor expanded in x around 0 71.0%
if -1.8500000000000001 < y Initial program 94.8%
sub-neg94.8%
log1p-def94.8%
distribute-neg-frac94.8%
sub-neg94.8%
distribute-neg-in94.8%
remove-double-neg94.8%
+-commutative94.8%
sub-neg94.8%
Simplified94.8%
Taylor expanded in y around 0 84.2%
+-commutative84.2%
div-sub84.2%
*-commutative84.2%
fma-def84.2%
sub-neg84.2%
mul-1-neg84.2%
*-inverses84.2%
metadata-eval84.2%
distribute-lft-in84.2%
+-commutative84.2%
metadata-eval84.2%
sub-neg84.2%
fma-def84.2%
*-lft-identity84.2%
sub-neg84.2%
metadata-eval84.2%
distribute-lft-in84.2%
metadata-eval84.2%
+-commutative84.2%
Simplified84.3%
Final simplification80.3%
(FPCore (x y) :precision binary64 (if (<= y -1.85e+21) (- 1.0 (log (/ -1.0 y))) (- 1.0 (log1p (- x)))))
double code(double x, double y) {
double tmp;
if (y <= -1.85e+21) {
tmp = 1.0 - log((-1.0 / y));
} else {
tmp = 1.0 - log1p(-x);
}
return tmp;
}
public static double code(double x, double y) {
double tmp;
if (y <= -1.85e+21) {
tmp = 1.0 - Math.log((-1.0 / y));
} else {
tmp = 1.0 - Math.log1p(-x);
}
return tmp;
}
def code(x, y): tmp = 0 if y <= -1.85e+21: tmp = 1.0 - math.log((-1.0 / y)) else: tmp = 1.0 - math.log1p(-x) return tmp
function code(x, y) tmp = 0.0 if (y <= -1.85e+21) tmp = Float64(1.0 - log(Float64(-1.0 / y))); else tmp = Float64(1.0 - log1p(Float64(-x))); end return tmp end
code[x_, y_] := If[LessEqual[y, -1.85e+21], N[(1.0 - N[Log[N[(-1.0 / y), $MachinePrecision]], $MachinePrecision]), $MachinePrecision], N[(1.0 - N[Log[1 + (-x)], $MachinePrecision]), $MachinePrecision]]
\begin{array}{l}
\\
\begin{array}{l}
\mathbf{if}\;y \leq -1.85 \cdot 10^{+21}:\\
\;\;\;\;1 - \log \left(\frac{-1}{y}\right)\\
\mathbf{else}:\\
\;\;\;\;1 - \mathsf{log1p}\left(-x\right)\\
\end{array}
\end{array}
if y < -1.85e21Initial program 16.4%
sub-neg16.4%
log1p-def16.4%
distribute-neg-frac16.4%
sub-neg16.4%
distribute-neg-in16.4%
remove-double-neg16.4%
+-commutative16.4%
sub-neg16.4%
Simplified16.4%
Taylor expanded in y around -inf 99.4%
sub-neg99.4%
metadata-eval99.4%
distribute-lft-in99.4%
metadata-eval99.4%
+-commutative99.4%
log1p-def99.4%
mul-1-neg99.4%
Simplified99.4%
Taylor expanded in x around 0 73.6%
if -1.85e21 < y Initial program 94.9%
sub-neg94.9%
log1p-def94.9%
distribute-neg-frac94.9%
sub-neg94.9%
distribute-neg-in94.9%
remove-double-neg94.9%
+-commutative94.9%
sub-neg94.9%
Simplified94.9%
Taylor expanded in y around 0 82.2%
log1p-def82.2%
mul-1-neg82.2%
Simplified82.2%
Final simplification79.7%
(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 72.5%
sub-neg72.5%
log1p-def72.5%
distribute-neg-frac72.5%
sub-neg72.5%
distribute-neg-in72.5%
remove-double-neg72.5%
+-commutative72.5%
sub-neg72.5%
Simplified72.5%
Taylor expanded in y around 0 62.3%
log1p-def62.3%
mul-1-neg62.3%
Simplified62.3%
Final simplification62.3%
(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 2024013
(FPCore (x y)
:name "Numeric.SpecFunctions:invIncompleteGamma from math-functions-0.1.5.2, B"
:precision binary64
:herbie-target
(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))))))