
(FPCore (a b c) :precision binary64 (let* ((t_0 (sqrt (- (* b b) (* (* 4.0 a) c))))) (if (>= b 0.0) (/ (- (- b) t_0) (* 2.0 a)) (/ (* 2.0 c) (+ (- b) t_0)))))
double code(double a, double b, double c) {
double t_0 = sqrt(((b * b) - ((4.0 * a) * c)));
double tmp;
if (b >= 0.0) {
tmp = (-b - t_0) / (2.0 * a);
} else {
tmp = (2.0 * c) / (-b + t_0);
}
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) - ((4.0d0 * a) * c)))
if (b >= 0.0d0) then
tmp = (-b - t_0) / (2.0d0 * a)
else
tmp = (2.0d0 * c) / (-b + t_0)
end if
code = tmp
end function
public static double code(double a, double b, double c) {
double t_0 = Math.sqrt(((b * b) - ((4.0 * a) * c)));
double tmp;
if (b >= 0.0) {
tmp = (-b - t_0) / (2.0 * a);
} else {
tmp = (2.0 * c) / (-b + t_0);
}
return tmp;
}
def code(a, b, c): t_0 = math.sqrt(((b * b) - ((4.0 * a) * c))) tmp = 0 if b >= 0.0: tmp = (-b - t_0) / (2.0 * a) else: tmp = (2.0 * c) / (-b + t_0) return tmp
function code(a, b, c) t_0 = sqrt(Float64(Float64(b * b) - Float64(Float64(4.0 * a) * c))) tmp = 0.0 if (b >= 0.0) tmp = Float64(Float64(Float64(-b) - t_0) / Float64(2.0 * a)); else tmp = Float64(Float64(2.0 * c) / Float64(Float64(-b) + t_0)); end return tmp end
function tmp_2 = code(a, b, c) t_0 = sqrt(((b * b) - ((4.0 * a) * c))); tmp = 0.0; if (b >= 0.0) tmp = (-b - t_0) / (2.0 * a); else tmp = (2.0 * c) / (-b + t_0); end tmp_2 = tmp; end
code[a_, b_, c_] := Block[{t$95$0 = N[Sqrt[N[(N[(b * b), $MachinePrecision] - N[(N[(4.0 * a), $MachinePrecision] * c), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]}, If[GreaterEqual[b, 0.0], N[(N[((-b) - t$95$0), $MachinePrecision] / N[(2.0 * a), $MachinePrecision]), $MachinePrecision], N[(N[(2.0 * c), $MachinePrecision] / N[((-b) + t$95$0), $MachinePrecision]), $MachinePrecision]]]
\begin{array}{l}
\\
\begin{array}{l}
t_0 := \sqrt{b \cdot b - \left(4 \cdot a\right) \cdot c}\\
\mathbf{if}\;b \geq 0:\\
\;\;\;\;\frac{\left(-b\right) - t\_0}{2 \cdot a}\\
\mathbf{else}:\\
\;\;\;\;\frac{2 \cdot c}{\left(-b\right) + t\_0}\\
\end{array}
\end{array}
Sampling outcomes in binary64 precision:
Herbie found 12 alternatives:
| Alternative | Accuracy | Speedup |
|---|
(FPCore (a b c) :precision binary64 (let* ((t_0 (sqrt (- (* b b) (* (* 4.0 a) c))))) (if (>= b 0.0) (/ (- (- b) t_0) (* 2.0 a)) (/ (* 2.0 c) (+ (- b) t_0)))))
double code(double a, double b, double c) {
double t_0 = sqrt(((b * b) - ((4.0 * a) * c)));
double tmp;
if (b >= 0.0) {
tmp = (-b - t_0) / (2.0 * a);
} else {
tmp = (2.0 * c) / (-b + t_0);
}
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) - ((4.0d0 * a) * c)))
if (b >= 0.0d0) then
tmp = (-b - t_0) / (2.0d0 * a)
else
tmp = (2.0d0 * c) / (-b + t_0)
end if
code = tmp
end function
public static double code(double a, double b, double c) {
double t_0 = Math.sqrt(((b * b) - ((4.0 * a) * c)));
double tmp;
if (b >= 0.0) {
tmp = (-b - t_0) / (2.0 * a);
} else {
tmp = (2.0 * c) / (-b + t_0);
}
return tmp;
}
def code(a, b, c): t_0 = math.sqrt(((b * b) - ((4.0 * a) * c))) tmp = 0 if b >= 0.0: tmp = (-b - t_0) / (2.0 * a) else: tmp = (2.0 * c) / (-b + t_0) return tmp
function code(a, b, c) t_0 = sqrt(Float64(Float64(b * b) - Float64(Float64(4.0 * a) * c))) tmp = 0.0 if (b >= 0.0) tmp = Float64(Float64(Float64(-b) - t_0) / Float64(2.0 * a)); else tmp = Float64(Float64(2.0 * c) / Float64(Float64(-b) + t_0)); end return tmp end
function tmp_2 = code(a, b, c) t_0 = sqrt(((b * b) - ((4.0 * a) * c))); tmp = 0.0; if (b >= 0.0) tmp = (-b - t_0) / (2.0 * a); else tmp = (2.0 * c) / (-b + t_0); end tmp_2 = tmp; end
code[a_, b_, c_] := Block[{t$95$0 = N[Sqrt[N[(N[(b * b), $MachinePrecision] - N[(N[(4.0 * a), $MachinePrecision] * c), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]}, If[GreaterEqual[b, 0.0], N[(N[((-b) - t$95$0), $MachinePrecision] / N[(2.0 * a), $MachinePrecision]), $MachinePrecision], N[(N[(2.0 * c), $MachinePrecision] / N[((-b) + t$95$0), $MachinePrecision]), $MachinePrecision]]]
\begin{array}{l}
\\
\begin{array}{l}
t_0 := \sqrt{b \cdot b - \left(4 \cdot a\right) \cdot c}\\
\mathbf{if}\;b \geq 0:\\
\;\;\;\;\frac{\left(-b\right) - t\_0}{2 \cdot a}\\
\mathbf{else}:\\
\;\;\;\;\frac{2 \cdot c}{\left(-b\right) + t\_0}\\
\end{array}
\end{array}
(FPCore (a b c)
:precision binary64
(let* ((t_0 (- (/ b a))))
(if (<= b -4e+157)
(if (>= b 0.0) (/ c b) (/ c (- b)))
(if (<= b 1e+101)
(if (>= b 0.0)
(/ (- (- b) (sqrt (- (* b b) (* c (* 4.0 a))))) (* a 2.0))
(/ (* c 2.0) (- (sqrt (fma (- (* 4.0 c)) a (* b b))) b)))
(if (>= b 0.0) t_0 t_0)))))
double code(double a, double b, double c) {
double t_0 = -(b / a);
double tmp_1;
if (b <= -4e+157) {
double tmp_2;
if (b >= 0.0) {
tmp_2 = c / b;
} else {
tmp_2 = c / -b;
}
tmp_1 = tmp_2;
} else if (b <= 1e+101) {
double tmp_3;
if (b >= 0.0) {
tmp_3 = (-b - sqrt(((b * b) - (c * (4.0 * a))))) / (a * 2.0);
} else {
tmp_3 = (c * 2.0) / (sqrt(fma(-(4.0 * c), a, (b * b))) - b);
}
tmp_1 = tmp_3;
} else if (b >= 0.0) {
tmp_1 = t_0;
} else {
tmp_1 = t_0;
}
return tmp_1;
}
function code(a, b, c) t_0 = Float64(-Float64(b / a)) tmp_1 = 0.0 if (b <= -4e+157) tmp_2 = 0.0 if (b >= 0.0) tmp_2 = Float64(c / b); else tmp_2 = Float64(c / Float64(-b)); end tmp_1 = tmp_2; elseif (b <= 1e+101) tmp_3 = 0.0 if (b >= 0.0) tmp_3 = Float64(Float64(Float64(-b) - sqrt(Float64(Float64(b * b) - Float64(c * Float64(4.0 * a))))) / Float64(a * 2.0)); else tmp_3 = Float64(Float64(c * 2.0) / Float64(sqrt(fma(Float64(-Float64(4.0 * c)), a, Float64(b * b))) - b)); end tmp_1 = tmp_3; elseif (b >= 0.0) tmp_1 = t_0; else tmp_1 = t_0; end return tmp_1 end
code[a_, b_, c_] := Block[{t$95$0 = (-N[(b / a), $MachinePrecision])}, If[LessEqual[b, -4e+157], If[GreaterEqual[b, 0.0], N[(c / b), $MachinePrecision], N[(c / (-b)), $MachinePrecision]], If[LessEqual[b, 1e+101], If[GreaterEqual[b, 0.0], N[(N[((-b) - N[Sqrt[N[(N[(b * b), $MachinePrecision] - N[(c * N[(4.0 * a), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision] / N[(a * 2.0), $MachinePrecision]), $MachinePrecision], N[(N[(c * 2.0), $MachinePrecision] / N[(N[Sqrt[N[((-N[(4.0 * c), $MachinePrecision]) * a + N[(b * b), $MachinePrecision]), $MachinePrecision]], $MachinePrecision] - b), $MachinePrecision]), $MachinePrecision]], If[GreaterEqual[b, 0.0], t$95$0, t$95$0]]]]
\begin{array}{l}
\\
\begin{array}{l}
t_0 := -\frac{b}{a}\\
\mathbf{if}\;b \leq -4 \cdot 10^{+157}:\\
\;\;\;\;\begin{array}{l}
\mathbf{if}\;b \geq 0:\\
\;\;\;\;\frac{c}{b}\\
\mathbf{else}:\\
\;\;\;\;\frac{c}{-b}\\
\end{array}\\
\mathbf{elif}\;b \leq 10^{+101}:\\
\;\;\;\;\begin{array}{l}
\mathbf{if}\;b \geq 0:\\
\;\;\;\;\frac{\left(-b\right) - \sqrt{b \cdot b - c \cdot \left(4 \cdot a\right)}}{a \cdot 2}\\
\mathbf{else}:\\
\;\;\;\;\frac{c \cdot 2}{\sqrt{\mathsf{fma}\left(-4 \cdot c, a, b \cdot b\right)} - b}\\
\end{array}\\
\mathbf{elif}\;b \geq 0:\\
\;\;\;\;t\_0\\
\mathbf{else}:\\
\;\;\;\;t\_0\\
\end{array}
\end{array}
if b < -3.99999999999999993e157Initial program 50.4%
Taylor expanded in c around 0
+-commutativeN/A
mul-1-negN/A
unsub-negN/A
lower--.f64N/A
lower-/.f64N/A
lower-/.f6450.4
Applied rewrites50.4%
Taylor expanded in c around 0
mul-1-negN/A
distribute-neg-frac2N/A
mul-1-negN/A
lower-/.f64N/A
mul-1-negN/A
lower-neg.f642.1
Applied rewrites2.1%
Taylor expanded in c around inf
Applied rewrites2.1%
Taylor expanded in b around -inf
mul-1-negN/A
distribute-neg-frac2N/A
mul-1-negN/A
lower-/.f64N/A
mul-1-negN/A
lower-neg.f64100.0
Applied rewrites100.0%
if -3.99999999999999993e157 < b < 9.9999999999999998e100Initial program 89.9%
lift--.f64N/A
sub-negN/A
+-commutativeN/A
lift-*.f64N/A
*-commutativeN/A
distribute-lft-neg-inN/A
lift-*.f64N/A
associate-*r*N/A
lower-fma.f64N/A
lower-*.f64N/A
lower-neg.f6489.9
Applied rewrites89.9%
if 9.9999999999999998e100 < b Initial program 59.6%
Taylor expanded in c around 0
+-commutativeN/A
mul-1-negN/A
unsub-negN/A
lower--.f64N/A
lower-/.f64N/A
lower-/.f6495.3
Applied rewrites95.3%
Taylor expanded in c around 0
mul-1-negN/A
distribute-neg-frac2N/A
mul-1-negN/A
lower-/.f64N/A
mul-1-negN/A
lower-neg.f6495.3
Applied rewrites95.3%
Taylor expanded in b around inf
mul-1-negN/A
distribute-neg-frac2N/A
mul-1-negN/A
lower-/.f64N/A
mul-1-negN/A
lower-neg.f6495.3
Applied rewrites95.3%
Final simplification92.5%
(FPCore (a b c)
:precision binary64
(let* ((t_0 (- (/ b a))) (t_1 (sqrt (- (* b b) (* c (* 4.0 a))))))
(if (<= b -4e+157)
(if (>= b 0.0) (/ c b) (/ c (- b)))
(if (<= b 1e+101)
(if (>= b 0.0) (/ (- (- b) t_1) (* a 2.0)) (/ (* c 2.0) (- t_1 b)))
(if (>= b 0.0) t_0 t_0)))))
double code(double a, double b, double c) {
double t_0 = -(b / a);
double t_1 = sqrt(((b * b) - (c * (4.0 * a))));
double tmp_1;
if (b <= -4e+157) {
double tmp_2;
if (b >= 0.0) {
tmp_2 = c / b;
} else {
tmp_2 = c / -b;
}
tmp_1 = tmp_2;
} else if (b <= 1e+101) {
double tmp_3;
if (b >= 0.0) {
tmp_3 = (-b - t_1) / (a * 2.0);
} else {
tmp_3 = (c * 2.0) / (t_1 - b);
}
tmp_1 = tmp_3;
} else if (b >= 0.0) {
tmp_1 = t_0;
} else {
tmp_1 = t_0;
}
return tmp_1;
}
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) :: t_1
real(8) :: tmp
real(8) :: tmp_1
real(8) :: tmp_2
real(8) :: tmp_3
t_0 = -(b / a)
t_1 = sqrt(((b * b) - (c * (4.0d0 * a))))
if (b <= (-4d+157)) then
if (b >= 0.0d0) then
tmp_2 = c / b
else
tmp_2 = c / -b
end if
tmp_1 = tmp_2
else if (b <= 1d+101) then
if (b >= 0.0d0) then
tmp_3 = (-b - t_1) / (a * 2.0d0)
else
tmp_3 = (c * 2.0d0) / (t_1 - b)
end if
tmp_1 = tmp_3
else if (b >= 0.0d0) then
tmp_1 = t_0
else
tmp_1 = t_0
end if
code = tmp_1
end function
public static double code(double a, double b, double c) {
double t_0 = -(b / a);
double t_1 = Math.sqrt(((b * b) - (c * (4.0 * a))));
double tmp_1;
if (b <= -4e+157) {
double tmp_2;
if (b >= 0.0) {
tmp_2 = c / b;
} else {
tmp_2 = c / -b;
}
tmp_1 = tmp_2;
} else if (b <= 1e+101) {
double tmp_3;
if (b >= 0.0) {
tmp_3 = (-b - t_1) / (a * 2.0);
} else {
tmp_3 = (c * 2.0) / (t_1 - b);
}
tmp_1 = tmp_3;
} else if (b >= 0.0) {
tmp_1 = t_0;
} else {
tmp_1 = t_0;
}
return tmp_1;
}
def code(a, b, c): t_0 = -(b / a) t_1 = math.sqrt(((b * b) - (c * (4.0 * a)))) tmp_1 = 0 if b <= -4e+157: tmp_2 = 0 if b >= 0.0: tmp_2 = c / b else: tmp_2 = c / -b tmp_1 = tmp_2 elif b <= 1e+101: tmp_3 = 0 if b >= 0.0: tmp_3 = (-b - t_1) / (a * 2.0) else: tmp_3 = (c * 2.0) / (t_1 - b) tmp_1 = tmp_3 elif b >= 0.0: tmp_1 = t_0 else: tmp_1 = t_0 return tmp_1
function code(a, b, c) t_0 = Float64(-Float64(b / a)) t_1 = sqrt(Float64(Float64(b * b) - Float64(c * Float64(4.0 * a)))) tmp_1 = 0.0 if (b <= -4e+157) tmp_2 = 0.0 if (b >= 0.0) tmp_2 = Float64(c / b); else tmp_2 = Float64(c / Float64(-b)); end tmp_1 = tmp_2; elseif (b <= 1e+101) tmp_3 = 0.0 if (b >= 0.0) tmp_3 = Float64(Float64(Float64(-b) - t_1) / Float64(a * 2.0)); else tmp_3 = Float64(Float64(c * 2.0) / Float64(t_1 - b)); end tmp_1 = tmp_3; elseif (b >= 0.0) tmp_1 = t_0; else tmp_1 = t_0; end return tmp_1 end
function tmp_5 = code(a, b, c) t_0 = -(b / a); t_1 = sqrt(((b * b) - (c * (4.0 * a)))); tmp_2 = 0.0; if (b <= -4e+157) tmp_3 = 0.0; if (b >= 0.0) tmp_3 = c / b; else tmp_3 = c / -b; end tmp_2 = tmp_3; elseif (b <= 1e+101) tmp_4 = 0.0; if (b >= 0.0) tmp_4 = (-b - t_1) / (a * 2.0); else tmp_4 = (c * 2.0) / (t_1 - b); end tmp_2 = tmp_4; elseif (b >= 0.0) tmp_2 = t_0; else tmp_2 = t_0; end tmp_5 = tmp_2; end
code[a_, b_, c_] := Block[{t$95$0 = (-N[(b / a), $MachinePrecision])}, Block[{t$95$1 = N[Sqrt[N[(N[(b * b), $MachinePrecision] - N[(c * N[(4.0 * a), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]}, If[LessEqual[b, -4e+157], If[GreaterEqual[b, 0.0], N[(c / b), $MachinePrecision], N[(c / (-b)), $MachinePrecision]], If[LessEqual[b, 1e+101], If[GreaterEqual[b, 0.0], N[(N[((-b) - t$95$1), $MachinePrecision] / N[(a * 2.0), $MachinePrecision]), $MachinePrecision], N[(N[(c * 2.0), $MachinePrecision] / N[(t$95$1 - b), $MachinePrecision]), $MachinePrecision]], If[GreaterEqual[b, 0.0], t$95$0, t$95$0]]]]]
\begin{array}{l}
\\
\begin{array}{l}
t_0 := -\frac{b}{a}\\
t_1 := \sqrt{b \cdot b - c \cdot \left(4 \cdot a\right)}\\
\mathbf{if}\;b \leq -4 \cdot 10^{+157}:\\
\;\;\;\;\begin{array}{l}
\mathbf{if}\;b \geq 0:\\
\;\;\;\;\frac{c}{b}\\
\mathbf{else}:\\
\;\;\;\;\frac{c}{-b}\\
\end{array}\\
\mathbf{elif}\;b \leq 10^{+101}:\\
\;\;\;\;\begin{array}{l}
\mathbf{if}\;b \geq 0:\\
\;\;\;\;\frac{\left(-b\right) - t\_1}{a \cdot 2}\\
\mathbf{else}:\\
\;\;\;\;\frac{c \cdot 2}{t\_1 - b}\\
\end{array}\\
\mathbf{elif}\;b \geq 0:\\
\;\;\;\;t\_0\\
\mathbf{else}:\\
\;\;\;\;t\_0\\
\end{array}
\end{array}
if b < -3.99999999999999993e157Initial program 37.7%
Taylor expanded in c around 0
+-commutativeN/A
mul-1-negN/A
unsub-negN/A
lower--.f64N/A
lower-/.f64N/A
lower-/.f6437.7
Applied rewrites37.7%
Taylor expanded in c around 0
mul-1-negN/A
distribute-neg-frac2N/A
mul-1-negN/A
lower-/.f64N/A
mul-1-negN/A
lower-neg.f642.3
Applied rewrites2.3%
Taylor expanded in c around inf
Applied rewrites2.3%
Taylor expanded in b around -inf
mul-1-negN/A
distribute-neg-frac2N/A
mul-1-negN/A
lower-/.f64N/A
mul-1-negN/A
lower-neg.f6498.5
Applied rewrites98.5%
if -3.99999999999999993e157 < b < 9.9999999999999998e100Initial program 86.7%
if 9.9999999999999998e100 < b Initial program 55.4%
Taylor expanded in c around 0
+-commutativeN/A
mul-1-negN/A
unsub-negN/A
lower--.f64N/A
lower-/.f64N/A
lower-/.f6497.0
Applied rewrites97.0%
Taylor expanded in c around 0
mul-1-negN/A
distribute-neg-frac2N/A
mul-1-negN/A
lower-/.f64N/A
mul-1-negN/A
lower-neg.f6497.0
Applied rewrites97.0%
Taylor expanded in b around inf
mul-1-negN/A
distribute-neg-frac2N/A
mul-1-negN/A
lower-/.f64N/A
mul-1-negN/A
lower-neg.f6496.8
Applied rewrites96.8%
Final simplification91.0%
herbie shell --seed 2024227
(FPCore (a b c)
:name "jeff quadratic root 1"
:precision binary64
(if (>= b 0.0) (/ (- (- b) (sqrt (- (* b b) (* (* 4.0 a) c)))) (* 2.0 a)) (/ (* 2.0 c) (+ (- b) (sqrt (- (* b b) (* (* 4.0 a) c)))))))