
(FPCore (a b c) :precision binary64 (/ (+ (- b) (sqrt (- (* b b) (* (* 4.0 a) c)))) (* 2.0 a)))
double code(double a, double b, double c) {
return (-b + sqrt(((b * b) - ((4.0 * a) * c)))) / (2.0 * a);
}
real(8) function code(a, b, c)
real(8), intent (in) :: a
real(8), intent (in) :: b
real(8), intent (in) :: c
code = (-b + sqrt(((b * b) - ((4.0d0 * a) * c)))) / (2.0d0 * a)
end function
public static double code(double a, double b, double c) {
return (-b + Math.sqrt(((b * b) - ((4.0 * a) * c)))) / (2.0 * a);
}
def code(a, b, c): return (-b + math.sqrt(((b * b) - ((4.0 * a) * c)))) / (2.0 * a)
function code(a, b, c) return Float64(Float64(Float64(-b) + sqrt(Float64(Float64(b * b) - Float64(Float64(4.0 * a) * c)))) / Float64(2.0 * a)) end
function tmp = code(a, b, c) tmp = (-b + sqrt(((b * b) - ((4.0 * a) * c)))) / (2.0 * a); end
code[a_, b_, c_] := N[(N[((-b) + N[Sqrt[N[(N[(b * b), $MachinePrecision] - N[(N[(4.0 * a), $MachinePrecision] * c), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision] / N[(2.0 * a), $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
\frac{\left(-b\right) + \sqrt{b \cdot b - \left(4 \cdot a\right) \cdot c}}{2 \cdot a}
\end{array}
Sampling outcomes in binary64 precision:
Herbie found 5 alternatives:
| Alternative | Accuracy | Speedup |
|---|
(FPCore (a b c) :precision binary64 (/ (+ (- b) (sqrt (- (* b b) (* (* 4.0 a) c)))) (* 2.0 a)))
double code(double a, double b, double c) {
return (-b + sqrt(((b * b) - ((4.0 * a) * c)))) / (2.0 * a);
}
real(8) function code(a, b, c)
real(8), intent (in) :: a
real(8), intent (in) :: b
real(8), intent (in) :: c
code = (-b + sqrt(((b * b) - ((4.0d0 * a) * c)))) / (2.0d0 * a)
end function
public static double code(double a, double b, double c) {
return (-b + Math.sqrt(((b * b) - ((4.0 * a) * c)))) / (2.0 * a);
}
def code(a, b, c): return (-b + math.sqrt(((b * b) - ((4.0 * a) * c)))) / (2.0 * a)
function code(a, b, c) return Float64(Float64(Float64(-b) + sqrt(Float64(Float64(b * b) - Float64(Float64(4.0 * a) * c)))) / Float64(2.0 * a)) end
function tmp = code(a, b, c) tmp = (-b + sqrt(((b * b) - ((4.0 * a) * c)))) / (2.0 * a); end
code[a_, b_, c_] := N[(N[((-b) + N[Sqrt[N[(N[(b * b), $MachinePrecision] - N[(N[(4.0 * a), $MachinePrecision] * c), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision] / N[(2.0 * a), $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
\frac{\left(-b\right) + \sqrt{b \cdot b - \left(4 \cdot a\right) \cdot c}}{2 \cdot a}
\end{array}
(FPCore (a b c)
:precision binary64
(fma
-2.0
(/ (* (pow a 2.0) (pow c 3.0)) (pow b 5.0))
(-
(-
(* (/ -0.25 a) (/ (* (pow (* a c) 4.0) 20.0) (pow b 7.0)))
(/ (* a (pow c 2.0)) (pow b 3.0)))
(/ c b))))
double code(double a, double b, double c) {
return fma(-2.0, ((pow(a, 2.0) * pow(c, 3.0)) / pow(b, 5.0)), ((((-0.25 / a) * ((pow((a * c), 4.0) * 20.0) / pow(b, 7.0))) - ((a * pow(c, 2.0)) / pow(b, 3.0))) - (c / b)));
}
function code(a, b, c) return fma(-2.0, Float64(Float64((a ^ 2.0) * (c ^ 3.0)) / (b ^ 5.0)), Float64(Float64(Float64(Float64(-0.25 / a) * Float64(Float64((Float64(a * c) ^ 4.0) * 20.0) / (b ^ 7.0))) - Float64(Float64(a * (c ^ 2.0)) / (b ^ 3.0))) - Float64(c / b))) end
code[a_, b_, c_] := N[(-2.0 * N[(N[(N[Power[a, 2.0], $MachinePrecision] * N[Power[c, 3.0], $MachinePrecision]), $MachinePrecision] / N[Power[b, 5.0], $MachinePrecision]), $MachinePrecision] + N[(N[(N[(N[(-0.25 / a), $MachinePrecision] * N[(N[(N[Power[N[(a * c), $MachinePrecision], 4.0], $MachinePrecision] * 20.0), $MachinePrecision] / N[Power[b, 7.0], $MachinePrecision]), $MachinePrecision]), $MachinePrecision] - N[(N[(a * N[Power[c, 2.0], $MachinePrecision]), $MachinePrecision] / N[Power[b, 3.0], $MachinePrecision]), $MachinePrecision]), $MachinePrecision] - N[(c / b), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
\mathsf{fma}\left(-2, \frac{{a}^{2} \cdot {c}^{3}}{{b}^{5}}, \left(\frac{-0.25}{a} \cdot \frac{{\left(a \cdot c\right)}^{4} \cdot 20}{{b}^{7}} - \frac{a \cdot {c}^{2}}{{b}^{3}}\right) - \frac{c}{b}\right)
\end{array}
Initial program 36.4%
*-commutative36.4%
Simplified36.4%
fma-neg36.4%
*-commutative36.4%
distribute-rgt-neg-in36.4%
distribute-lft-neg-in36.4%
metadata-eval36.4%
*-commutative36.4%
add-cbrt-cube36.2%
pow336.2%
*-commutative36.2%
associate-*r*36.2%
Applied egg-rr36.2%
Taylor expanded in b around inf 94.3%
Simplified94.3%
Final simplification94.3%
(FPCore (a b c) :precision binary64 (- (- (/ -2.0 (/ (pow b 5.0) (* (pow a 2.0) (pow c 3.0)))) (/ c b)) (* (pow c 2.0) (/ a (pow b 3.0)))))
double code(double a, double b, double c) {
return ((-2.0 / (pow(b, 5.0) / (pow(a, 2.0) * pow(c, 3.0)))) - (c / b)) - (pow(c, 2.0) * (a / pow(b, 3.0)));
}
real(8) function code(a, b, c)
real(8), intent (in) :: a
real(8), intent (in) :: b
real(8), intent (in) :: c
code = (((-2.0d0) / ((b ** 5.0d0) / ((a ** 2.0d0) * (c ** 3.0d0)))) - (c / b)) - ((c ** 2.0d0) * (a / (b ** 3.0d0)))
end function
public static double code(double a, double b, double c) {
return ((-2.0 / (Math.pow(b, 5.0) / (Math.pow(a, 2.0) * Math.pow(c, 3.0)))) - (c / b)) - (Math.pow(c, 2.0) * (a / Math.pow(b, 3.0)));
}
def code(a, b, c): return ((-2.0 / (math.pow(b, 5.0) / (math.pow(a, 2.0) * math.pow(c, 3.0)))) - (c / b)) - (math.pow(c, 2.0) * (a / math.pow(b, 3.0)))
function code(a, b, c) return Float64(Float64(Float64(-2.0 / Float64((b ^ 5.0) / Float64((a ^ 2.0) * (c ^ 3.0)))) - Float64(c / b)) - Float64((c ^ 2.0) * Float64(a / (b ^ 3.0)))) end
function tmp = code(a, b, c) tmp = ((-2.0 / ((b ^ 5.0) / ((a ^ 2.0) * (c ^ 3.0)))) - (c / b)) - ((c ^ 2.0) * (a / (b ^ 3.0))); end
code[a_, b_, c_] := N[(N[(N[(-2.0 / N[(N[Power[b, 5.0], $MachinePrecision] / N[(N[Power[a, 2.0], $MachinePrecision] * N[Power[c, 3.0], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision] - N[(c / b), $MachinePrecision]), $MachinePrecision] - N[(N[Power[c, 2.0], $MachinePrecision] * N[(a / N[Power[b, 3.0], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
\left(\frac{-2}{\frac{{b}^{5}}{{a}^{2} \cdot {c}^{3}}} - \frac{c}{b}\right) - {c}^{2} \cdot \frac{a}{{b}^{3}}
\end{array}
Initial program 36.4%
*-commutative36.4%
Simplified36.4%
Taylor expanded in b around inf 92.4%
associate-+r+92.4%
mul-1-neg92.4%
unsub-neg92.4%
mul-1-neg92.4%
unsub-neg92.4%
associate-*r/92.4%
associate-/l*92.4%
*-commutative92.4%
associate-/l*92.4%
associate-/r/92.4%
Simplified92.4%
Final simplification92.4%
(FPCore (a b c) :precision binary64 (let* ((t_0 (/ (- (sqrt (- (* b b) (* c (* a 4.0)))) b) (* a 2.0)))) (if (<= t_0 -8e-6) t_0 (/ (- c) b))))
double code(double a, double b, double c) {
double t_0 = (sqrt(((b * b) - (c * (a * 4.0)))) - b) / (a * 2.0);
double tmp;
if (t_0 <= -8e-6) {
tmp = t_0;
} else {
tmp = -c / b;
}
return tmp;
}
real(8) function code(a, b, c)
real(8), intent (in) :: a
real(8), intent (in) :: b
real(8), intent (in) :: c
real(8) :: t_0
real(8) :: tmp
t_0 = (sqrt(((b * b) - (c * (a * 4.0d0)))) - b) / (a * 2.0d0)
if (t_0 <= (-8d-6)) then
tmp = t_0
else
tmp = -c / b
end if
code = tmp
end function
public static double code(double a, double b, double c) {
double t_0 = (Math.sqrt(((b * b) - (c * (a * 4.0)))) - b) / (a * 2.0);
double tmp;
if (t_0 <= -8e-6) {
tmp = t_0;
} else {
tmp = -c / b;
}
return tmp;
}
def code(a, b, c): t_0 = (math.sqrt(((b * b) - (c * (a * 4.0)))) - b) / (a * 2.0) tmp = 0 if t_0 <= -8e-6: tmp = t_0 else: tmp = -c / b return tmp
function code(a, b, c) t_0 = Float64(Float64(sqrt(Float64(Float64(b * b) - Float64(c * Float64(a * 4.0)))) - b) / Float64(a * 2.0)) tmp = 0.0 if (t_0 <= -8e-6) tmp = t_0; else tmp = Float64(Float64(-c) / b); end return tmp end
function tmp_2 = code(a, b, c) t_0 = (sqrt(((b * b) - (c * (a * 4.0)))) - b) / (a * 2.0); tmp = 0.0; if (t_0 <= -8e-6) tmp = t_0; else tmp = -c / b; end tmp_2 = tmp; end
code[a_, b_, c_] := Block[{t$95$0 = N[(N[(N[Sqrt[N[(N[(b * b), $MachinePrecision] - N[(c * N[(a * 4.0), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision] - b), $MachinePrecision] / N[(a * 2.0), $MachinePrecision]), $MachinePrecision]}, If[LessEqual[t$95$0, -8e-6], t$95$0, N[((-c) / b), $MachinePrecision]]]
\begin{array}{l}
\\
\begin{array}{l}
t_0 := \frac{\sqrt{b \cdot b - c \cdot \left(a \cdot 4\right)} - b}{a \cdot 2}\\
\mathbf{if}\;t_0 \leq -8 \cdot 10^{-6}:\\
\;\;\;\;t_0\\
\mathbf{else}:\\
\;\;\;\;\frac{-c}{b}\\
\end{array}
\end{array}
if (/.f64 (+.f64 (neg.f64 b) (sqrt.f64 (-.f64 (*.f64 b b) (*.f64 (*.f64 4 a) c)))) (*.f64 2 a)) < -7.99999999999999964e-6Initial program 69.2%
if -7.99999999999999964e-6 < (/.f64 (+.f64 (neg.f64 b) (sqrt.f64 (-.f64 (*.f64 b b) (*.f64 (*.f64 4 a) c)))) (*.f64 2 a)) Initial program 22.3%
*-commutative22.3%
Simplified22.3%
Taylor expanded in b around inf 88.0%
mul-1-neg88.0%
distribute-neg-frac88.0%
Simplified88.0%
Final simplification82.4%
(FPCore (a b c) :precision binary64 (- (* (pow c 2.0) (/ (- a) (pow b 3.0))) (/ c b)))
double code(double a, double b, double c) {
return (pow(c, 2.0) * (-a / pow(b, 3.0))) - (c / b);
}
real(8) function code(a, b, c)
real(8), intent (in) :: a
real(8), intent (in) :: b
real(8), intent (in) :: c
code = ((c ** 2.0d0) * (-a / (b ** 3.0d0))) - (c / b)
end function
public static double code(double a, double b, double c) {
return (Math.pow(c, 2.0) * (-a / Math.pow(b, 3.0))) - (c / b);
}
def code(a, b, c): return (math.pow(c, 2.0) * (-a / math.pow(b, 3.0))) - (c / b)
function code(a, b, c) return Float64(Float64((c ^ 2.0) * Float64(Float64(-a) / (b ^ 3.0))) - Float64(c / b)) end
function tmp = code(a, b, c) tmp = ((c ^ 2.0) * (-a / (b ^ 3.0))) - (c / b); end
code[a_, b_, c_] := N[(N[(N[Power[c, 2.0], $MachinePrecision] * N[((-a) / N[Power[b, 3.0], $MachinePrecision]), $MachinePrecision]), $MachinePrecision] - N[(c / b), $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
{c}^{2} \cdot \frac{-a}{{b}^{3}} - \frac{c}{b}
\end{array}
Initial program 36.4%
*-commutative36.4%
Simplified36.4%
Taylor expanded in b around inf 89.1%
mul-1-neg89.1%
unsub-neg89.1%
mul-1-neg89.1%
distribute-neg-frac89.1%
associate-/l*89.1%
associate-/r/89.1%
Simplified89.1%
Final simplification89.1%
(FPCore (a b c) :precision binary64 (/ (- c) b))
double code(double a, double b, double c) {
return -c / b;
}
real(8) function code(a, b, c)
real(8), intent (in) :: a
real(8), intent (in) :: b
real(8), intent (in) :: c
code = -c / b
end function
public static double code(double a, double b, double c) {
return -c / b;
}
def code(a, b, c): return -c / b
function code(a, b, c) return Float64(Float64(-c) / b) end
function tmp = code(a, b, c) tmp = -c / b; end
code[a_, b_, c_] := N[((-c) / b), $MachinePrecision]
\begin{array}{l}
\\
\frac{-c}{b}
\end{array}
Initial program 36.4%
*-commutative36.4%
Simplified36.4%
Taylor expanded in b around inf 77.7%
mul-1-neg77.7%
distribute-neg-frac77.7%
Simplified77.7%
Final simplification77.7%
herbie shell --seed 2024021
(FPCore (a b c)
:name "Quadratic roots, medium range"
:precision binary64
:pre (and (and (and (< 1.1102230246251565e-16 a) (< a 9007199254740992.0)) (and (< 1.1102230246251565e-16 b) (< b 9007199254740992.0))) (and (< 1.1102230246251565e-16 c) (< c 9007199254740992.0)))
(/ (+ (- b) (sqrt (- (* b b) (* (* 4.0 a) c)))) (* 2.0 a)))