
(FPCore (g h a) :precision binary64 (let* ((t_0 (/ 1.0 (* 2.0 a))) (t_1 (sqrt (- (* g g) (* h h))))) (+ (cbrt (* t_0 (+ (- g) t_1))) (cbrt (* t_0 (- (- g) t_1))))))
double code(double g, double h, double a) {
double t_0 = 1.0 / (2.0 * a);
double t_1 = sqrt(((g * g) - (h * h)));
return cbrt((t_0 * (-g + t_1))) + cbrt((t_0 * (-g - t_1)));
}
public static double code(double g, double h, double a) {
double t_0 = 1.0 / (2.0 * a);
double t_1 = Math.sqrt(((g * g) - (h * h)));
return Math.cbrt((t_0 * (-g + t_1))) + Math.cbrt((t_0 * (-g - t_1)));
}
function code(g, h, a) t_0 = Float64(1.0 / Float64(2.0 * a)) t_1 = sqrt(Float64(Float64(g * g) - Float64(h * h))) return Float64(cbrt(Float64(t_0 * Float64(Float64(-g) + t_1))) + cbrt(Float64(t_0 * Float64(Float64(-g) - t_1)))) end
code[g_, h_, a_] := Block[{t$95$0 = N[(1.0 / N[(2.0 * a), $MachinePrecision]), $MachinePrecision]}, Block[{t$95$1 = N[Sqrt[N[(N[(g * g), $MachinePrecision] - N[(h * h), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]}, N[(N[Power[N[(t$95$0 * N[((-g) + t$95$1), $MachinePrecision]), $MachinePrecision], 1/3], $MachinePrecision] + N[Power[N[(t$95$0 * N[((-g) - t$95$1), $MachinePrecision]), $MachinePrecision], 1/3], $MachinePrecision]), $MachinePrecision]]]
\begin{array}{l}
\\
\begin{array}{l}
t_0 := \frac{1}{2 \cdot a}\\
t_1 := \sqrt{g \cdot g - h \cdot h}\\
\sqrt[3]{t\_0 \cdot \left(\left(-g\right) + t\_1\right)} + \sqrt[3]{t\_0 \cdot \left(\left(-g\right) - t\_1\right)}
\end{array}
\end{array}
Sampling outcomes in binary64 precision:
Herbie found 6 alternatives:
| Alternative | Accuracy | Speedup |
|---|
(FPCore (g h a) :precision binary64 (let* ((t_0 (/ 1.0 (* 2.0 a))) (t_1 (sqrt (- (* g g) (* h h))))) (+ (cbrt (* t_0 (+ (- g) t_1))) (cbrt (* t_0 (- (- g) t_1))))))
double code(double g, double h, double a) {
double t_0 = 1.0 / (2.0 * a);
double t_1 = sqrt(((g * g) - (h * h)));
return cbrt((t_0 * (-g + t_1))) + cbrt((t_0 * (-g - t_1)));
}
public static double code(double g, double h, double a) {
double t_0 = 1.0 / (2.0 * a);
double t_1 = Math.sqrt(((g * g) - (h * h)));
return Math.cbrt((t_0 * (-g + t_1))) + Math.cbrt((t_0 * (-g - t_1)));
}
function code(g, h, a) t_0 = Float64(1.0 / Float64(2.0 * a)) t_1 = sqrt(Float64(Float64(g * g) - Float64(h * h))) return Float64(cbrt(Float64(t_0 * Float64(Float64(-g) + t_1))) + cbrt(Float64(t_0 * Float64(Float64(-g) - t_1)))) end
code[g_, h_, a_] := Block[{t$95$0 = N[(1.0 / N[(2.0 * a), $MachinePrecision]), $MachinePrecision]}, Block[{t$95$1 = N[Sqrt[N[(N[(g * g), $MachinePrecision] - N[(h * h), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]}, N[(N[Power[N[(t$95$0 * N[((-g) + t$95$1), $MachinePrecision]), $MachinePrecision], 1/3], $MachinePrecision] + N[Power[N[(t$95$0 * N[((-g) - t$95$1), $MachinePrecision]), $MachinePrecision], 1/3], $MachinePrecision]), $MachinePrecision]]]
\begin{array}{l}
\\
\begin{array}{l}
t_0 := \frac{1}{2 \cdot a}\\
t_1 := \sqrt{g \cdot g - h \cdot h}\\
\sqrt[3]{t\_0 \cdot \left(\left(-g\right) + t\_1\right)} + \sqrt[3]{t\_0 \cdot \left(\left(-g\right) - t\_1\right)}
\end{array}
\end{array}
(FPCore (g h a)
:precision binary64
(if (<= g 8.2e-122)
(+ (cbrt (/ (- g) a)) (* (cbrt h) (cbrt (* (/ h g) (/ -0.25 a)))))
(*
(cbrt (/ 0.5 a))
(+
(cbrt (- (fma (sqrt (- g h)) (sqrt (+ h g)) g)))
(cbrt (* (* (/ h g) h) -0.5))))))
double code(double g, double h, double a) {
double tmp;
if (g <= 8.2e-122) {
tmp = cbrt((-g / a)) + (cbrt(h) * cbrt(((h / g) * (-0.25 / a))));
} else {
tmp = cbrt((0.5 / a)) * (cbrt(-fma(sqrt((g - h)), sqrt((h + g)), g)) + cbrt((((h / g) * h) * -0.5)));
}
return tmp;
}
function code(g, h, a) tmp = 0.0 if (g <= 8.2e-122) tmp = Float64(cbrt(Float64(Float64(-g) / a)) + Float64(cbrt(h) * cbrt(Float64(Float64(h / g) * Float64(-0.25 / a))))); else tmp = Float64(cbrt(Float64(0.5 / a)) * Float64(cbrt(Float64(-fma(sqrt(Float64(g - h)), sqrt(Float64(h + g)), g))) + cbrt(Float64(Float64(Float64(h / g) * h) * -0.5)))); end return tmp end
code[g_, h_, a_] := If[LessEqual[g, 8.2e-122], N[(N[Power[N[((-g) / a), $MachinePrecision], 1/3], $MachinePrecision] + N[(N[Power[h, 1/3], $MachinePrecision] * N[Power[N[(N[(h / g), $MachinePrecision] * N[(-0.25 / a), $MachinePrecision]), $MachinePrecision], 1/3], $MachinePrecision]), $MachinePrecision]), $MachinePrecision], N[(N[Power[N[(0.5 / a), $MachinePrecision], 1/3], $MachinePrecision] * N[(N[Power[(-N[(N[Sqrt[N[(g - h), $MachinePrecision]], $MachinePrecision] * N[Sqrt[N[(h + g), $MachinePrecision]], $MachinePrecision] + g), $MachinePrecision]), 1/3], $MachinePrecision] + N[Power[N[(N[(N[(h / g), $MachinePrecision] * h), $MachinePrecision] * -0.5), $MachinePrecision], 1/3], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]]
\begin{array}{l}
\\
\begin{array}{l}
\mathbf{if}\;g \leq 8.2 \cdot 10^{-122}:\\
\;\;\;\;\sqrt[3]{\frac{-g}{a}} + \sqrt[3]{h} \cdot \sqrt[3]{\frac{h}{g} \cdot \frac{-0.25}{a}}\\
\mathbf{else}:\\
\;\;\;\;\sqrt[3]{\frac{0.5}{a}} \cdot \left(\sqrt[3]{-\mathsf{fma}\left(\sqrt{g - h}, \sqrt{h + g}, g\right)} + \sqrt[3]{\left(\frac{h}{g} \cdot h\right) \cdot -0.5}\right)\\
\end{array}
\end{array}
if g < 8.2000000000000001e-122Initial program 48.3%
Taylor expanded in g around inf
associate-*r/N/A
mul-1-negN/A
lower-/.f64N/A
lower-neg.f6414.5
Applied rewrites14.5%
Taylor expanded in g around inf
lower-*.f64N/A
lower-cbrt.f64N/A
unpow2N/A
times-fracN/A
lower-*.f64N/A
lower-/.f64N/A
lower-/.f64N/A
lower-*.f64N/A
lower-cbrt.f64N/A
lower-cbrt.f6482.4
Applied rewrites82.4%
lift-+.f64N/A
+-commutativeN/A
lower-+.f6482.4
Applied rewrites82.4%
Applied rewrites82.5%
if 8.2000000000000001e-122 < g Initial program 42.8%
Taylor expanded in g around inf
*-commutativeN/A
lower-*.f64N/A
lower-/.f64N/A
unpow2N/A
lower-*.f6444.9
Applied rewrites44.9%
Applied rewrites98.1%
(FPCore (g h a) :precision binary64 (+ (cbrt (/ (- g) a)) (* (cbrt (/ -0.25 a)) (cbrt (* (/ h g) h)))))
double code(double g, double h, double a) {
return cbrt((-g / a)) + (cbrt((-0.25 / a)) * cbrt(((h / g) * h)));
}
public static double code(double g, double h, double a) {
return Math.cbrt((-g / a)) + (Math.cbrt((-0.25 / a)) * Math.cbrt(((h / g) * h)));
}
function code(g, h, a) return Float64(cbrt(Float64(Float64(-g) / a)) + Float64(cbrt(Float64(-0.25 / a)) * cbrt(Float64(Float64(h / g) * h)))) end
code[g_, h_, a_] := N[(N[Power[N[((-g) / a), $MachinePrecision], 1/3], $MachinePrecision] + N[(N[Power[N[(-0.25 / a), $MachinePrecision], 1/3], $MachinePrecision] * N[Power[N[(N[(h / g), $MachinePrecision] * h), $MachinePrecision], 1/3], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
\sqrt[3]{\frac{-g}{a}} + \sqrt[3]{\frac{-0.25}{a}} \cdot \sqrt[3]{\frac{h}{g} \cdot h}
\end{array}
Initial program 45.9%
Taylor expanded in g around inf
associate-*r/N/A
mul-1-negN/A
lower-/.f64N/A
lower-neg.f6426.6
Applied rewrites26.6%
Taylor expanded in g around inf
lower-*.f64N/A
lower-cbrt.f64N/A
unpow2N/A
times-fracN/A
lower-*.f64N/A
lower-/.f64N/A
lower-/.f64N/A
lower-*.f64N/A
lower-cbrt.f64N/A
lower-cbrt.f6479.0
Applied rewrites79.0%
lift-+.f64N/A
+-commutativeN/A
lower-+.f6479.0
Applied rewrites79.0%
Applied rewrites79.1%
(FPCore (g h a) :precision binary64 (+ (cbrt (/ (- g) a)) (* (cbrt h) (cbrt (* (/ h g) (/ -0.25 a))))))
double code(double g, double h, double a) {
return cbrt((-g / a)) + (cbrt(h) * cbrt(((h / g) * (-0.25 / a))));
}
public static double code(double g, double h, double a) {
return Math.cbrt((-g / a)) + (Math.cbrt(h) * Math.cbrt(((h / g) * (-0.25 / a))));
}
function code(g, h, a) return Float64(cbrt(Float64(Float64(-g) / a)) + Float64(cbrt(h) * cbrt(Float64(Float64(h / g) * Float64(-0.25 / a))))) end
code[g_, h_, a_] := N[(N[Power[N[((-g) / a), $MachinePrecision], 1/3], $MachinePrecision] + N[(N[Power[h, 1/3], $MachinePrecision] * N[Power[N[(N[(h / g), $MachinePrecision] * N[(-0.25 / a), $MachinePrecision]), $MachinePrecision], 1/3], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
\sqrt[3]{\frac{-g}{a}} + \sqrt[3]{h} \cdot \sqrt[3]{\frac{h}{g} \cdot \frac{-0.25}{a}}
\end{array}
Initial program 45.9%
Taylor expanded in g around inf
associate-*r/N/A
mul-1-negN/A
lower-/.f64N/A
lower-neg.f6426.6
Applied rewrites26.6%
Taylor expanded in g around inf
lower-*.f64N/A
lower-cbrt.f64N/A
unpow2N/A
times-fracN/A
lower-*.f64N/A
lower-/.f64N/A
lower-/.f64N/A
lower-*.f64N/A
lower-cbrt.f64N/A
lower-cbrt.f6479.0
Applied rewrites79.0%
lift-+.f64N/A
+-commutativeN/A
lower-+.f6479.0
Applied rewrites79.0%
Applied rewrites79.1%
(FPCore (g h a) :precision binary64 (+ (cbrt (* (pow (* 2.0 a) -1.0) (* -2.0 g))) (cbrt (/ (- g) a))))
double code(double g, double h, double a) {
return cbrt((pow((2.0 * a), -1.0) * (-2.0 * g))) + cbrt((-g / a));
}
public static double code(double g, double h, double a) {
return Math.cbrt((Math.pow((2.0 * a), -1.0) * (-2.0 * g))) + Math.cbrt((-g / a));
}
function code(g, h, a) return Float64(cbrt(Float64((Float64(2.0 * a) ^ -1.0) * Float64(-2.0 * g))) + cbrt(Float64(Float64(-g) / a))) end
code[g_, h_, a_] := N[(N[Power[N[(N[Power[N[(2.0 * a), $MachinePrecision], -1.0], $MachinePrecision] * N[(-2.0 * g), $MachinePrecision]), $MachinePrecision], 1/3], $MachinePrecision] + N[Power[N[((-g) / a), $MachinePrecision], 1/3], $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
\sqrt[3]{{\left(2 \cdot a\right)}^{-1} \cdot \left(-2 \cdot g\right)} + \sqrt[3]{\frac{-g}{a}}
\end{array}
Initial program 45.9%
Taylor expanded in g around inf
associate-*r/N/A
mul-1-negN/A
lower-/.f64N/A
lower-neg.f6426.6
Applied rewrites26.6%
Taylor expanded in g around -inf
lower-*.f6415.8
Applied rewrites15.8%
Final simplification15.8%
(FPCore (g h a) :precision binary64 (+ (cbrt (/ (- g) a)) (cbrt (* -0.25 (* (/ h g) (/ h a))))))
double code(double g, double h, double a) {
return cbrt((-g / a)) + cbrt((-0.25 * ((h / g) * (h / a))));
}
public static double code(double g, double h, double a) {
return Math.cbrt((-g / a)) + Math.cbrt((-0.25 * ((h / g) * (h / a))));
}
function code(g, h, a) return Float64(cbrt(Float64(Float64(-g) / a)) + cbrt(Float64(-0.25 * Float64(Float64(h / g) * Float64(h / a))))) end
code[g_, h_, a_] := N[(N[Power[N[((-g) / a), $MachinePrecision], 1/3], $MachinePrecision] + N[Power[N[(-0.25 * N[(N[(h / g), $MachinePrecision] * N[(h / a), $MachinePrecision]), $MachinePrecision]), $MachinePrecision], 1/3], $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
\sqrt[3]{\frac{-g}{a}} + \sqrt[3]{-0.25 \cdot \left(\frac{h}{g} \cdot \frac{h}{a}\right)}
\end{array}
Initial program 45.9%
Taylor expanded in g around inf
associate-*r/N/A
mul-1-negN/A
lower-/.f64N/A
lower-neg.f6426.6
Applied rewrites26.6%
Taylor expanded in g around inf
lower-*.f64N/A
lower-cbrt.f64N/A
unpow2N/A
times-fracN/A
lower-*.f64N/A
lower-/.f64N/A
lower-/.f64N/A
lower-*.f64N/A
lower-cbrt.f64N/A
lower-cbrt.f6479.0
Applied rewrites79.0%
lift-+.f64N/A
+-commutativeN/A
lower-+.f6479.0
Applied rewrites79.0%
(FPCore (g h a) :precision binary64 0.0)
double code(double g, double h, double a) {
return 0.0;
}
real(8) function code(g, h, a)
real(8), intent (in) :: g
real(8), intent (in) :: h
real(8), intent (in) :: a
code = 0.0d0
end function
public static double code(double g, double h, double a) {
return 0.0;
}
def code(g, h, a): return 0.0
function code(g, h, a) return 0.0 end
function tmp = code(g, h, a) tmp = 0.0; end
code[g_, h_, a_] := 0.0
\begin{array}{l}
\\
0
\end{array}
Initial program 45.9%
lift-cbrt.f64N/A
lift-*.f64N/A
lift-/.f64N/A
associate-*l/N/A
lift-*.f64N/A
times-fracN/A
metadata-evalN/A
cbrt-prodN/A
pow1/3N/A
lower-*.f64N/A
pow1/3N/A
lower-cbrt.f64N/A
lower-cbrt.f64N/A
lower-/.f6445.9
Applied rewrites45.9%
Taylor expanded in g around -inf
mul-1-negN/A
lower-neg.f64N/A
lower-*.f64N/A
lower-cbrt.f64N/A
lower-/.f64N/A
+-commutativeN/A
unpow2N/A
rem-square-sqrtN/A
metadata-evalN/A
lower-*.f64N/A
lower-cbrt.f643.0
Applied rewrites3.0%
Applied rewrites3.0%
herbie shell --seed 2024307
(FPCore (g h a)
:name "2-ancestry mixing, positive discriminant"
:precision binary64
(+ (cbrt (* (/ 1.0 (* 2.0 a)) (+ (- g) (sqrt (- (* g g) (* h h)))))) (cbrt (* (/ 1.0 (* 2.0 a)) (- (- g) (sqrt (- (* g g) (* h h))))))))