
(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 (* (/ 1.0 (+ b (sqrt (fma c (* a -4.0) (* b b))))) (* c -2.0)))
double code(double a, double b, double c) {
return (1.0 / (b + sqrt(fma(c, (a * -4.0), (b * b))))) * (c * -2.0);
}
function code(a, b, c) return Float64(Float64(1.0 / Float64(b + sqrt(fma(c, Float64(a * -4.0), Float64(b * b))))) * Float64(c * -2.0)) end
code[a_, b_, c_] := N[(N[(1.0 / N[(b + N[Sqrt[N[(c * N[(a * -4.0), $MachinePrecision] + N[(b * b), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]), $MachinePrecision] * N[(c * -2.0), $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
\frac{1}{b + \sqrt{\mathsf{fma}\left(c, a \cdot -4, b \cdot b\right)}} \cdot \left(c \cdot -2\right)
\end{array}
Initial program 31.7%
lift-+.f64N/A
+-commutativeN/A
lift-neg.f64N/A
unsub-negN/A
lower--.f6431.7
lift--.f64N/A
sub-negN/A
+-commutativeN/A
lift-*.f64N/A
distribute-lft-neg-inN/A
lift-*.f64N/A
*-commutativeN/A
distribute-rgt-neg-inN/A
associate-*l*N/A
lower-fma.f64N/A
lower-*.f64N/A
metadata-eval31.7
Applied rewrites31.7%
lift--.f64N/A
flip--N/A
div-invN/A
lower-*.f64N/A
Applied rewrites32.5%
lift-/.f64N/A
lift-*.f64N/A
*-commutativeN/A
associate-/l*N/A
lower-*.f64N/A
Applied rewrites99.5%
Taylor expanded in c around 0
lower-*.f6499.6
Applied rewrites99.6%
Final simplification99.6%
(FPCore (a b c) :precision binary64 (if (<= (/ (- (sqrt (- (* b b) (* c (* a 4.0)))) b) (* a 2.0)) -500000.0) (/ (- (sqrt (fma b b (* -4.0 (* c a)))) b) (* a 2.0)) (- (fma a (/ (* c c) (* b (* b b))) (/ c b)))))
double code(double a, double b, double c) {
double tmp;
if (((sqrt(((b * b) - (c * (a * 4.0)))) - b) / (a * 2.0)) <= -500000.0) {
tmp = (sqrt(fma(b, b, (-4.0 * (c * a)))) - b) / (a * 2.0);
} else {
tmp = -fma(a, ((c * c) / (b * (b * b))), (c / b));
}
return tmp;
}
function code(a, b, c) tmp = 0.0 if (Float64(Float64(sqrt(Float64(Float64(b * b) - Float64(c * Float64(a * 4.0)))) - b) / Float64(a * 2.0)) <= -500000.0) tmp = Float64(Float64(sqrt(fma(b, b, Float64(-4.0 * Float64(c * a)))) - b) / Float64(a * 2.0)); else tmp = Float64(-fma(a, Float64(Float64(c * c) / Float64(b * Float64(b * b))), Float64(c / b))); end return tmp end
code[a_, b_, c_] := If[LessEqual[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], -500000.0], N[(N[(N[Sqrt[N[(b * b + N[(-4.0 * N[(c * a), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision] - b), $MachinePrecision] / N[(a * 2.0), $MachinePrecision]), $MachinePrecision], (-N[(a * N[(N[(c * c), $MachinePrecision] / N[(b * N[(b * b), $MachinePrecision]), $MachinePrecision]), $MachinePrecision] + N[(c / b), $MachinePrecision]), $MachinePrecision])]
\begin{array}{l}
\\
\begin{array}{l}
\mathbf{if}\;\frac{\sqrt{b \cdot b - c \cdot \left(a \cdot 4\right)} - b}{a \cdot 2} \leq -500000:\\
\;\;\;\;\frac{\sqrt{\mathsf{fma}\left(b, b, -4 \cdot \left(c \cdot a\right)\right)} - b}{a \cdot 2}\\
\mathbf{else}:\\
\;\;\;\;-\mathsf{fma}\left(a, \frac{c \cdot c}{b \cdot \left(b \cdot b\right)}, \frac{c}{b}\right)\\
\end{array}
\end{array}
if (/.f64 (+.f64 (neg.f64 b) (sqrt.f64 (-.f64 (*.f64 b b) (*.f64 (*.f64 #s(literal 4 binary64) a) c)))) (*.f64 #s(literal 2 binary64) a)) < -5e5Initial program 76.6%
lift-+.f64N/A
+-commutativeN/A
lift-neg.f64N/A
unsub-negN/A
lower--.f6476.6
lift--.f64N/A
sub-negN/A
+-commutativeN/A
lift-*.f64N/A
distribute-lft-neg-inN/A
lift-*.f64N/A
*-commutativeN/A
distribute-rgt-neg-inN/A
associate-*l*N/A
lower-fma.f64N/A
lower-*.f64N/A
metadata-eval76.6
Applied rewrites76.6%
lift-fma.f64N/A
+-commutativeN/A
lift-*.f64N/A
lower-fma.f64N/A
*-commutativeN/A
lift-*.f64N/A
associate-*l*N/A
*-commutativeN/A
lift-*.f64N/A
lower-*.f6476.7
lift-*.f64N/A
*-commutativeN/A
lower-*.f6476.7
Applied rewrites76.7%
if -5e5 < (/.f64 (+.f64 (neg.f64 b) (sqrt.f64 (-.f64 (*.f64 b b) (*.f64 (*.f64 #s(literal 4 binary64) a) c)))) (*.f64 #s(literal 2 binary64) a)) Initial program 29.9%
lift-+.f64N/A
+-commutativeN/A
lift-neg.f64N/A
unsub-negN/A
lower--.f6429.9
lift--.f64N/A
sub-negN/A
+-commutativeN/A
lift-*.f64N/A
distribute-lft-neg-inN/A
lift-*.f64N/A
*-commutativeN/A
distribute-rgt-neg-inN/A
associate-*l*N/A
lower-fma.f64N/A
lower-*.f64N/A
metadata-eval29.9
Applied rewrites29.9%
Taylor expanded in a around 0
+-commutativeN/A
mul-1-negN/A
mul-1-negN/A
distribute-neg-outN/A
lower-neg.f64N/A
associate-/l*N/A
lower-fma.f64N/A
lower-/.f64N/A
unpow2N/A
lower-*.f64N/A
cube-multN/A
unpow2N/A
lower-*.f64N/A
unpow2N/A
lower-*.f64N/A
lower-/.f6491.8
Applied rewrites91.8%
Final simplification91.2%
herbie shell --seed 2024226
(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)))