Quadratic roots, wide range

Average Error: 52.9 → 0.4

Time: 2.4min

Precision: binary64

Cost: 14276

\[\left(\left(4.930380657631324 \cdot 10^{-32} < a \land a < 2.028240960365167 \cdot 10^{+31}\right) \land \left(4.930380657631324 \cdot 10^{-32} < b \land b < 2.028240960365167 \cdot 10^{+31}\right)\right) \land \left(4.930380657631324 \cdot 10^{-32} < c \land c < 2.028240960365167 \cdot 10^{+31}\right)\]

Math

\[\frac{\left(-b\right) + \sqrt{b \cdot b - \left(4 \cdot a\right) \cdot c}}{2 \cdot a} \]

↓

\[\frac{\begin{array}{l} \mathbf{if}\;b \cdot 2 \ne 0:\\ \;\;\;\;\left(-4 \cdot c\right) \cdot \frac{a}{b + \sqrt{\mathsf{fma}\left(b, b, \left(c \cdot a\right) \cdot -4\right)}}\\ \mathbf{else}:\\ \;\;\;\;\sqrt{\mathsf{fma}\left(b, b, -4 \cdot \left(c \cdot a\right)\right)} - b\\ \end{array}}{2 \cdot a} \]

FPCore

(FPCore (a b c)
 :precision binary64
 (/ (+ (- b) (sqrt (- (* b b) (* (* 4.0 a) c)))) (* 2.0 a)))

↓

(FPCore (a b c)
 :precision binary64
 (/
  (if (!= (* b 2.0) 0.0)
    (* (* -4.0 c) (/ a (+ b (sqrt (fma b b (* (* c a) -4.0))))))
    (- (sqrt (fma b b (* -4.0 (* c a)))) b))
  (* 2.0 a)))

C

double code(double a, double b, double c) {
	return (-b + sqrt(((b * b) - ((4.0 * a) * c)))) / (2.0 * a);
}

↓

double code(double a, double b, double c) {
	double tmp;
	if ((b * 2.0) != 0.0) {
		tmp = (-4.0 * c) * (a / (b + sqrt(fma(b, b, ((c * a) * -4.0)))));
	} else {
		tmp = sqrt(fma(b, b, (-4.0 * (c * a)))) - b;
	}
	return tmp / (2.0 * a);
}

Julia

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 code(a, b, c)
	tmp = 0.0
	if (Float64(b * 2.0) != 0.0)
		tmp = Float64(Float64(-4.0 * c) * Float64(a / Float64(b + sqrt(fma(b, b, Float64(Float64(c * a) * -4.0))))));
	else
		tmp = Float64(sqrt(fma(b, b, Float64(-4.0 * Float64(c * a)))) - b);
	end
	return Float64(tmp / Float64(2.0 * a))
end

Wolfram

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]

↓

code[a_, b_, c_] := N[(If[Unequal[N[(b * 2.0), $MachinePrecision], 0.0], N[(N[(-4.0 * c), $MachinePrecision] * N[(a / N[(b + N[Sqrt[N[(b * b + N[(N[(c * a), $MachinePrecision] * -4.0), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision], N[(N[Sqrt[N[(b * b + N[(-4.0 * N[(c * a), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision] - b), $MachinePrecision]] / N[(2.0 * a), $MachinePrecision]), $MachinePrecision]

Derivation

1. Initial program 52.9

   \[\frac{\left(-b\right) + \sqrt{b \cdot b - \left(4 \cdot a\right) \cdot c}}{2 \cdot a} \]

2. Applied egg-rr 52.9

   \[\leadsto \frac{\color{blue}{\begin{array}{l} \color{blue}{\mathbf{if}\;b + \sqrt{\mathsf{fma}\left(a \cdot -4, c, {b}^{2}\right)} \ne 0:\\ \;\;\;\;\frac{{\left(\sqrt{\mathsf{fma}\left(a \cdot -4, c, {b}^{2}\right)}\right)}^{2} - {b}^{2}}{b + \sqrt{\mathsf{fma}\left(a \cdot -4, c, {b}^{2}\right)}}\\ \mathbf{else}:\\ \;\;\;\;\sqrt{\mathsf{fma}\left(a \cdot -4, c, {b}^{2}\right)} - b\\ } \end{array}}}{2 \cdot a} \]

3. Simplified 52.9

   \[\leadsto \frac{\color{blue}{\begin{array}{l} \color{blue}{\mathbf{if}\;b + \sqrt{\mathsf{fma}\left(b, b, -4 \cdot \left(c \cdot a\right)\right)} \ne 0:\\ \;\;\;\;\frac{{\left(\sqrt{\mathsf{fma}\left(b, b, -4 \cdot \left(c \cdot a\right)\right)}\right)}^{2} - {b}^{2}}{b + \sqrt{\mathsf{fma}\left(b, b, -4 \cdot \left(c \cdot a\right)\right)}}\\ \mathbf{else}:\\ \;\;\;\;\sqrt{\mathsf{fma}\left(b, b, -4 \cdot \left(c \cdot a\right)\right)} - b\\ } \end{array}}}{2 \cdot a} \]
   Proof

4. Taylor expanded in b around inf 0.4

   \[\leadsto \frac{\begin{array}{l} \mathbf{if}\;b + \sqrt{\mathsf{fma}\left(b, b, -4 \cdot \left(c \cdot a\right)\right)} \ne 0:\\ \;\;\;\;\frac{\color{blue}{-4 \cdot \left(c \cdot a\right)}}{b + \sqrt{\mathsf{fma}\left(b, b, -4 \cdot \left(c \cdot a\right)\right)}}\\ \mathbf{else}:\\ \;\;\;\;\sqrt{\mathsf{fma}\left(b, b, -4 \cdot \left(c \cdot a\right)\right)} - b\\ \end{array}}{2 \cdot a} \]

5. Taylor expanded in b around inf 0.4

   \[\leadsto \frac{\begin{array}{l} \mathbf{if}\;\color{blue}{2 \cdot b} \ne 0:\\ \;\;\;\;\frac{-4 \cdot \left(c \cdot a\right)}{b + \sqrt{\mathsf{fma}\left(b, b, -4 \cdot \left(c \cdot a\right)\right)}}\\ \mathbf{else}:\\ \;\;\;\;\sqrt{\mathsf{fma}\left(b, b, -4 \cdot \left(c \cdot a\right)\right)} - b\\ \end{array}}{2 \cdot a} \]

6. Simplified 0.4

   \[\leadsto \frac{\begin{array}{l} \mathbf{if}\;\color{blue}{b \cdot 2} \ne 0:\\ \;\;\;\;\frac{-4 \cdot \left(c \cdot a\right)}{b + \sqrt{\mathsf{fma}\left(b, b, -4 \cdot \left(c \cdot a\right)\right)}}\\ \mathbf{else}:\\ \;\;\;\;\sqrt{\mathsf{fma}\left(b, b, -4 \cdot \left(c \cdot a\right)\right)} - b\\ \end{array}}{2 \cdot a} \]
   Proof

7. Applied egg-rr 0.4

   \[\leadsto \frac{\begin{array}{l} \mathbf{if}\;b \cdot 2 \ne 0:\\ \;\;\;\;\color{blue}{\left(-4 \cdot c\right) \cdot \frac{a}{b + \sqrt{\mathsf{fma}\left(b, b, \left(c \cdot a\right) \cdot -4\right)}}}\\ \mathbf{else}:\\ \;\;\;\;\sqrt{\mathsf{fma}\left(b, b, -4 \cdot \left(c \cdot a\right)\right)} - b\\ \end{array}}{2 \cdot a} \]

Alternatives

Alternative 1
Error	0.5
Cost	14276

\[\begin{array}{l} t_0 := \left(c \cdot a\right) \cdot -4\\ t_1 := \sqrt{\mathsf{fma}\left(b, b, t_0\right)}\\ \begin{array}{l} \mathbf{if}\;2 \cdot b \ne 0:\\ \;\;\;\;\frac{t_0}{b + t_1}\\ \mathbf{else}:\\ \;\;\;\;t_1 - b\\ \end{array} \cdot \frac{0.5}{a} \end{array} \]

Alternative 2
Error	2.8
Cost	7232

\[\frac{-\left(c \cdot c\right) \cdot a}{{b}^{3}} - \frac{c}{b} \]

Alternative 3
Error	5.9
Cost	256

\[\frac{-c}{b} \]

Reproduce

herbie shell --seed 2023033 
(FPCore (a b c)
  :name "Quadratic roots, wide range"
  :precision binary64
  :pre (and (and (and (< 4.930380657631324e-32 a) (< a 2.028240960365167e+31)) (and (< 4.930380657631324e-32 b) (< b 2.028240960365167e+31))) (and (< 4.930380657631324e-32 c) (< c 2.028240960365167e+31)))
  (/ (+ (- b) (sqrt (- (* b b) (* (* 4.0 a) c)))) (* 2.0 a)))