
(FPCore (c_p c_n t s)
:precision binary64
(let* ((t_1 (/ 1.0 (+ 1.0 (exp (- t))))) (t_2 (/ 1.0 (+ 1.0 (exp (- s))))))
(/
(* (pow t_2 c_p) (pow (- 1.0 t_2) c_n))
(* (pow t_1 c_p) (pow (- 1.0 t_1) c_n)))))
double code(double c_p, double c_n, double t, double s) {
double t_1 = 1.0 / (1.0 + exp(-t));
double t_2 = 1.0 / (1.0 + exp(-s));
return (pow(t_2, c_p) * pow((1.0 - t_2), c_n)) / (pow(t_1, c_p) * pow((1.0 - t_1), c_n));
}
real(8) function code(c_p, c_n, t, s)
real(8), intent (in) :: c_p
real(8), intent (in) :: c_n
real(8), intent (in) :: t
real(8), intent (in) :: s
real(8) :: t_1
real(8) :: t_2
t_1 = 1.0d0 / (1.0d0 + exp(-t))
t_2 = 1.0d0 / (1.0d0 + exp(-s))
code = ((t_2 ** c_p) * ((1.0d0 - t_2) ** c_n)) / ((t_1 ** c_p) * ((1.0d0 - t_1) ** c_n))
end function
public static double code(double c_p, double c_n, double t, double s) {
double t_1 = 1.0 / (1.0 + Math.exp(-t));
double t_2 = 1.0 / (1.0 + Math.exp(-s));
return (Math.pow(t_2, c_p) * Math.pow((1.0 - t_2), c_n)) / (Math.pow(t_1, c_p) * Math.pow((1.0 - t_1), c_n));
}
def code(c_p, c_n, t, s): t_1 = 1.0 / (1.0 + math.exp(-t)) t_2 = 1.0 / (1.0 + math.exp(-s)) return (math.pow(t_2, c_p) * math.pow((1.0 - t_2), c_n)) / (math.pow(t_1, c_p) * math.pow((1.0 - t_1), c_n))
function code(c_p, c_n, t, s) t_1 = Float64(1.0 / Float64(1.0 + exp(Float64(-t)))) t_2 = Float64(1.0 / Float64(1.0 + exp(Float64(-s)))) return Float64(Float64((t_2 ^ c_p) * (Float64(1.0 - t_2) ^ c_n)) / Float64((t_1 ^ c_p) * (Float64(1.0 - t_1) ^ c_n))) end
function tmp = code(c_p, c_n, t, s) t_1 = 1.0 / (1.0 + exp(-t)); t_2 = 1.0 / (1.0 + exp(-s)); tmp = ((t_2 ^ c_p) * ((1.0 - t_2) ^ c_n)) / ((t_1 ^ c_p) * ((1.0 - t_1) ^ c_n)); end
code[c$95$p_, c$95$n_, t_, s_] := Block[{t$95$1 = N[(1.0 / N[(1.0 + N[Exp[(-t)], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]}, Block[{t$95$2 = N[(1.0 / N[(1.0 + N[Exp[(-s)], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]}, N[(N[(N[Power[t$95$2, c$95$p], $MachinePrecision] * N[Power[N[(1.0 - t$95$2), $MachinePrecision], c$95$n], $MachinePrecision]), $MachinePrecision] / N[(N[Power[t$95$1, c$95$p], $MachinePrecision] * N[Power[N[(1.0 - t$95$1), $MachinePrecision], c$95$n], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]]]
\begin{array}{l}
\\
\begin{array}{l}
t_1 := \frac{1}{1 + e^{-t}}\\
t_2 := \frac{1}{1 + e^{-s}}\\
\frac{{t\_2}^{c\_p} \cdot {\left(1 - t\_2\right)}^{c\_n}}{{t\_1}^{c\_p} \cdot {\left(1 - t\_1\right)}^{c\_n}}
\end{array}
\end{array}
Sampling outcomes in binary64 precision:
Herbie found 5 alternatives:
| Alternative | Accuracy | Speedup |
|---|
(FPCore (c_p c_n t s)
:precision binary64
(let* ((t_1 (/ 1.0 (+ 1.0 (exp (- t))))) (t_2 (/ 1.0 (+ 1.0 (exp (- s))))))
(/
(* (pow t_2 c_p) (pow (- 1.0 t_2) c_n))
(* (pow t_1 c_p) (pow (- 1.0 t_1) c_n)))))
double code(double c_p, double c_n, double t, double s) {
double t_1 = 1.0 / (1.0 + exp(-t));
double t_2 = 1.0 / (1.0 + exp(-s));
return (pow(t_2, c_p) * pow((1.0 - t_2), c_n)) / (pow(t_1, c_p) * pow((1.0 - t_1), c_n));
}
real(8) function code(c_p, c_n, t, s)
real(8), intent (in) :: c_p
real(8), intent (in) :: c_n
real(8), intent (in) :: t
real(8), intent (in) :: s
real(8) :: t_1
real(8) :: t_2
t_1 = 1.0d0 / (1.0d0 + exp(-t))
t_2 = 1.0d0 / (1.0d0 + exp(-s))
code = ((t_2 ** c_p) * ((1.0d0 - t_2) ** c_n)) / ((t_1 ** c_p) * ((1.0d0 - t_1) ** c_n))
end function
public static double code(double c_p, double c_n, double t, double s) {
double t_1 = 1.0 / (1.0 + Math.exp(-t));
double t_2 = 1.0 / (1.0 + Math.exp(-s));
return (Math.pow(t_2, c_p) * Math.pow((1.0 - t_2), c_n)) / (Math.pow(t_1, c_p) * Math.pow((1.0 - t_1), c_n));
}
def code(c_p, c_n, t, s): t_1 = 1.0 / (1.0 + math.exp(-t)) t_2 = 1.0 / (1.0 + math.exp(-s)) return (math.pow(t_2, c_p) * math.pow((1.0 - t_2), c_n)) / (math.pow(t_1, c_p) * math.pow((1.0 - t_1), c_n))
function code(c_p, c_n, t, s) t_1 = Float64(1.0 / Float64(1.0 + exp(Float64(-t)))) t_2 = Float64(1.0 / Float64(1.0 + exp(Float64(-s)))) return Float64(Float64((t_2 ^ c_p) * (Float64(1.0 - t_2) ^ c_n)) / Float64((t_1 ^ c_p) * (Float64(1.0 - t_1) ^ c_n))) end
function tmp = code(c_p, c_n, t, s) t_1 = 1.0 / (1.0 + exp(-t)); t_2 = 1.0 / (1.0 + exp(-s)); tmp = ((t_2 ^ c_p) * ((1.0 - t_2) ^ c_n)) / ((t_1 ^ c_p) * ((1.0 - t_1) ^ c_n)); end
code[c$95$p_, c$95$n_, t_, s_] := Block[{t$95$1 = N[(1.0 / N[(1.0 + N[Exp[(-t)], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]}, Block[{t$95$2 = N[(1.0 / N[(1.0 + N[Exp[(-s)], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]}, N[(N[(N[Power[t$95$2, c$95$p], $MachinePrecision] * N[Power[N[(1.0 - t$95$2), $MachinePrecision], c$95$n], $MachinePrecision]), $MachinePrecision] / N[(N[Power[t$95$1, c$95$p], $MachinePrecision] * N[Power[N[(1.0 - t$95$1), $MachinePrecision], c$95$n], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]]]
\begin{array}{l}
\\
\begin{array}{l}
t_1 := \frac{1}{1 + e^{-t}}\\
t_2 := \frac{1}{1 + e^{-s}}\\
\frac{{t\_2}^{c\_p} \cdot {\left(1 - t\_2\right)}^{c\_n}}{{t\_1}^{c\_p} \cdot {\left(1 - t\_1\right)}^{c\_n}}
\end{array}
\end{array}
(FPCore (c_p c_n t s)
:precision binary64
(let* ((t_1 (exp (- s))) (t_2 (exp (- t))) (t_3 (log1p (/ -1.0 (+ 1.0 t_2)))))
(if (<= c_p 5e-21)
(exp (* c_n (- (log1p (/ -1.0 (+ 1.0 t_1))) t_3)))
(exp
(fma
c_p
(- (log1p t_2) (log1p t_1))
(* c_n (- (log1p (/ 1.0 (- -1.0 t_1))) t_3)))))))
double code(double c_p, double c_n, double t, double s) {
double t_1 = exp(-s);
double t_2 = exp(-t);
double t_3 = log1p((-1.0 / (1.0 + t_2)));
double tmp;
if (c_p <= 5e-21) {
tmp = exp((c_n * (log1p((-1.0 / (1.0 + t_1))) - t_3)));
} else {
tmp = exp(fma(c_p, (log1p(t_2) - log1p(t_1)), (c_n * (log1p((1.0 / (-1.0 - t_1))) - t_3))));
}
return tmp;
}
function code(c_p, c_n, t, s) t_1 = exp(Float64(-s)) t_2 = exp(Float64(-t)) t_3 = log1p(Float64(-1.0 / Float64(1.0 + t_2))) tmp = 0.0 if (c_p <= 5e-21) tmp = exp(Float64(c_n * Float64(log1p(Float64(-1.0 / Float64(1.0 + t_1))) - t_3))); else tmp = exp(fma(c_p, Float64(log1p(t_2) - log1p(t_1)), Float64(c_n * Float64(log1p(Float64(1.0 / Float64(-1.0 - t_1))) - t_3)))); end return tmp end
code[c$95$p_, c$95$n_, t_, s_] := Block[{t$95$1 = N[Exp[(-s)], $MachinePrecision]}, Block[{t$95$2 = N[Exp[(-t)], $MachinePrecision]}, Block[{t$95$3 = N[Log[1 + N[(-1.0 / N[(1.0 + t$95$2), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]}, If[LessEqual[c$95$p, 5e-21], N[Exp[N[(c$95$n * N[(N[Log[1 + N[(-1.0 / N[(1.0 + t$95$1), $MachinePrecision]), $MachinePrecision]], $MachinePrecision] - t$95$3), $MachinePrecision]), $MachinePrecision]], $MachinePrecision], N[Exp[N[(c$95$p * N[(N[Log[1 + t$95$2], $MachinePrecision] - N[Log[1 + t$95$1], $MachinePrecision]), $MachinePrecision] + N[(c$95$n * N[(N[Log[1 + N[(1.0 / N[(-1.0 - t$95$1), $MachinePrecision]), $MachinePrecision]], $MachinePrecision] - t$95$3), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]]]]]
\begin{array}{l}
\\
\begin{array}{l}
t_1 := e^{-s}\\
t_2 := e^{-t}\\
t_3 := \mathsf{log1p}\left(\frac{-1}{1 + t\_2}\right)\\
\mathbf{if}\;c\_p \leq 5 \cdot 10^{-21}:\\
\;\;\;\;e^{c\_n \cdot \left(\mathsf{log1p}\left(\frac{-1}{1 + t\_1}\right) - t\_3\right)}\\
\mathbf{else}:\\
\;\;\;\;e^{\mathsf{fma}\left(c\_p, \mathsf{log1p}\left(t\_2\right) - \mathsf{log1p}\left(t\_1\right), c\_n \cdot \left(\mathsf{log1p}\left(\frac{1}{-1 - t\_1}\right) - t\_3\right)\right)}\\
\end{array}
\end{array}
if c_p < 4.99999999999999973e-21Initial program 92.3%
Applied rewrites95.3%
Taylor expanded in c_p around 0
lower-*.f64N/A
lower--.f64N/A
sub-negN/A
lower-log1p.f64N/A
distribute-neg-fracN/A
metadata-evalN/A
lower-/.f64N/A
lower-+.f64N/A
lower-exp.f64N/A
lower-neg.f64N/A
sub-negN/A
lower-log1p.f64N/A
distribute-neg-fracN/A
Applied rewrites98.7%
if 4.99999999999999973e-21 < c_p Initial program 59.2%
Applied rewrites100.0%
Final simplification98.9%
(FPCore (c_p c_n t s)
:precision binary64
(let* ((t_1 (+ 1.0 (exp (- t)))) (t_2 (+ 1.0 (exp (- s)))))
(if (<= c_p 2e-27)
(exp (* c_n (- (log1p (/ -1.0 t_2)) (log1p (/ -1.0 t_1)))))
(pow (* t_1 (/ 1.0 t_2)) c_p))))
double code(double c_p, double c_n, double t, double s) {
double t_1 = 1.0 + exp(-t);
double t_2 = 1.0 + exp(-s);
double tmp;
if (c_p <= 2e-27) {
tmp = exp((c_n * (log1p((-1.0 / t_2)) - log1p((-1.0 / t_1)))));
} else {
tmp = pow((t_1 * (1.0 / t_2)), c_p);
}
return tmp;
}
public static double code(double c_p, double c_n, double t, double s) {
double t_1 = 1.0 + Math.exp(-t);
double t_2 = 1.0 + Math.exp(-s);
double tmp;
if (c_p <= 2e-27) {
tmp = Math.exp((c_n * (Math.log1p((-1.0 / t_2)) - Math.log1p((-1.0 / t_1)))));
} else {
tmp = Math.pow((t_1 * (1.0 / t_2)), c_p);
}
return tmp;
}
def code(c_p, c_n, t, s): t_1 = 1.0 + math.exp(-t) t_2 = 1.0 + math.exp(-s) tmp = 0 if c_p <= 2e-27: tmp = math.exp((c_n * (math.log1p((-1.0 / t_2)) - math.log1p((-1.0 / t_1))))) else: tmp = math.pow((t_1 * (1.0 / t_2)), c_p) return tmp
function code(c_p, c_n, t, s) t_1 = Float64(1.0 + exp(Float64(-t))) t_2 = Float64(1.0 + exp(Float64(-s))) tmp = 0.0 if (c_p <= 2e-27) tmp = exp(Float64(c_n * Float64(log1p(Float64(-1.0 / t_2)) - log1p(Float64(-1.0 / t_1))))); else tmp = Float64(t_1 * Float64(1.0 / t_2)) ^ c_p; end return tmp end
code[c$95$p_, c$95$n_, t_, s_] := Block[{t$95$1 = N[(1.0 + N[Exp[(-t)], $MachinePrecision]), $MachinePrecision]}, Block[{t$95$2 = N[(1.0 + N[Exp[(-s)], $MachinePrecision]), $MachinePrecision]}, If[LessEqual[c$95$p, 2e-27], N[Exp[N[(c$95$n * N[(N[Log[1 + N[(-1.0 / t$95$2), $MachinePrecision]], $MachinePrecision] - N[Log[1 + N[(-1.0 / t$95$1), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision], N[Power[N[(t$95$1 * N[(1.0 / t$95$2), $MachinePrecision]), $MachinePrecision], c$95$p], $MachinePrecision]]]]
\begin{array}{l}
\\
\begin{array}{l}
t_1 := 1 + e^{-t}\\
t_2 := 1 + e^{-s}\\
\mathbf{if}\;c\_p \leq 2 \cdot 10^{-27}:\\
\;\;\;\;e^{c\_n \cdot \left(\mathsf{log1p}\left(\frac{-1}{t\_2}\right) - \mathsf{log1p}\left(\frac{-1}{t\_1}\right)\right)}\\
\mathbf{else}:\\
\;\;\;\;{\left(t\_1 \cdot \frac{1}{t\_2}\right)}^{c\_p}\\
\end{array}
\end{array}
if c_p < 2.0000000000000001e-27Initial program 94.0%
Applied rewrites96.5%
Taylor expanded in c_p around 0
lower-*.f64N/A
lower--.f64N/A
sub-negN/A
lower-log1p.f64N/A
distribute-neg-fracN/A
metadata-evalN/A
lower-/.f64N/A
lower-+.f64N/A
lower-exp.f64N/A
lower-neg.f64N/A
sub-negN/A
lower-log1p.f64N/A
distribute-neg-fracN/A
Applied rewrites99.4%
if 2.0000000000000001e-27 < c_p Initial program 76.1%
Taylor expanded in c_n around 0
lower-/.f64N/A
lower-pow.f64N/A
lower-/.f64N/A
lower-+.f64N/A
lower-exp.f64N/A
lower-neg.f64N/A
lower-pow.f64N/A
lower-/.f64N/A
lower-+.f64N/A
lower-exp.f64N/A
lower-neg.f6478.7
Applied rewrites78.7%
Applied rewrites96.3%
Final simplification98.9%
(FPCore (c_p c_n t s) :precision binary64 (* (pow (/ (+ 1.0 (exp (- t))) (+ 1.0 (exp (- s)))) c_p) (pow (/ (+ 1.0 (exp t)) (+ 1.0 (exp s))) c_n)))
double code(double c_p, double c_n, double t, double s) {
return pow(((1.0 + exp(-t)) / (1.0 + exp(-s))), c_p) * pow(((1.0 + exp(t)) / (1.0 + exp(s))), c_n);
}
real(8) function code(c_p, c_n, t, s)
real(8), intent (in) :: c_p
real(8), intent (in) :: c_n
real(8), intent (in) :: t
real(8), intent (in) :: s
code = (((1.0d0 + exp(-t)) / (1.0d0 + exp(-s))) ** c_p) * (((1.0d0 + exp(t)) / (1.0d0 + exp(s))) ** c_n)
end function
public static double code(double c_p, double c_n, double t, double s) {
return Math.pow(((1.0 + Math.exp(-t)) / (1.0 + Math.exp(-s))), c_p) * Math.pow(((1.0 + Math.exp(t)) / (1.0 + Math.exp(s))), c_n);
}
def code(c_p, c_n, t, s): return math.pow(((1.0 + math.exp(-t)) / (1.0 + math.exp(-s))), c_p) * math.pow(((1.0 + math.exp(t)) / (1.0 + math.exp(s))), c_n)
function code(c_p, c_n, t, s) return Float64((Float64(Float64(1.0 + exp(Float64(-t))) / Float64(1.0 + exp(Float64(-s)))) ^ c_p) * (Float64(Float64(1.0 + exp(t)) / Float64(1.0 + exp(s))) ^ c_n)) end
function tmp = code(c_p, c_n, t, s) tmp = (((1.0 + exp(-t)) / (1.0 + exp(-s))) ^ c_p) * (((1.0 + exp(t)) / (1.0 + exp(s))) ^ c_n); end
code[c$95$p_, c$95$n_, t_, s_] := N[(N[Power[N[(N[(1.0 + N[Exp[(-t)], $MachinePrecision]), $MachinePrecision] / N[(1.0 + N[Exp[(-s)], $MachinePrecision]), $MachinePrecision]), $MachinePrecision], c$95$p], $MachinePrecision] * N[Power[N[(N[(1.0 + N[Exp[t], $MachinePrecision]), $MachinePrecision] / N[(1.0 + N[Exp[s], $MachinePrecision]), $MachinePrecision]), $MachinePrecision], c$95$n], $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
{\left(\frac{1 + e^{-t}}{1 + e^{-s}}\right)}^{c\_p} \cdot {\left(\frac{1 + e^{t}}{1 + e^{s}}\right)}^{c\_n}
\end{array}
herbie shell --seed 2024223
(FPCore (c_p c_n t s)
:name "Harley's example"
:precision binary64
:pre (and (< 0.0 c_p) (< 0.0 c_n))
:alt
(! :herbie-platform default (* (pow (/ (+ 1 (exp (- t))) (+ 1 (exp (- s)))) c_p) (pow (/ (+ 1 (exp t)) (+ 1 (exp s))) c_n)))
(/ (* (pow (/ 1.0 (+ 1.0 (exp (- s)))) c_p) (pow (- 1.0 (/ 1.0 (+ 1.0 (exp (- s))))) c_n)) (* (pow (/ 1.0 (+ 1.0 (exp (- t)))) c_p) (pow (- 1.0 (/ 1.0 (+ 1.0 (exp (- t))))) c_n))))