quadp (p42, positive)

Average Error: 34.4 → 8.0

Time: 1.7min

Precision: binary64

Cost: 63248

Math

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

↓

\[\begin{array}{l} t_0 := \frac{\left(-b\right) + \sqrt{b \cdot b - 4 \cdot \left(a \cdot c\right)}}{2 \cdot a}\\ t_1 := \mathsf{fma}\left(0.5, \frac{\left|b\right| - b}{a}, \frac{\left(-0.25 \cdot a\right) \cdot {\left(\frac{2 \cdot c}{\left|b\right|}\right)}^{2} - c}{\left|b\right|}\right)\\ \mathbf{if}\;t_0 \leq -\infty:\\ \;\;\;\;t_1\\ \mathbf{elif}\;t_0 \leq -4 \cdot 10^{-240}:\\ \;\;\;\;t_0\\ \mathbf{elif}\;t_0 \leq 0:\\ \;\;\;\;\frac{-c}{b}\\ \mathbf{elif}\;t_0 \leq 4 \cdot 10^{+216}:\\ \;\;\;\;\frac{\sqrt{\mathsf{fma}\left(b, b, -4 \cdot \left(a \cdot c\right)\right)}}{a + a} - \frac{b}{a + a}\\ \mathbf{else}:\\ \;\;\;\;t_1\\ \end{array} \]

FPCore

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

↓

(FPCore (a b c)
 :precision binary64
 (let* ((t_0 (/ (+ (- b) (sqrt (- (* b b) (* 4.0 (* a c))))) (* 2.0 a)))
        (t_1
         (fma
          0.5
          (/ (- (fabs b) b) a)
          (/
           (- (* (* -0.25 a) (pow (/ (* 2.0 c) (fabs b)) 2.0)) c)
           (fabs b)))))
   (if (<= t_0 (- INFINITY))
     t_1
     (if (<= t_0 -4e-240)
       t_0
       (if (<= t_0 0.0)
         (/ (- c) b)
         (if (<= t_0 4e+216)
           (- (/ (sqrt (fma b b (* -4.0 (* a c)))) (+ a a)) (/ b (+ a a)))
           t_1))))))

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 t_0 = (-b + sqrt(((b * b) - (4.0 * (a * c))))) / (2.0 * a);
	double t_1 = fma(0.5, ((fabs(b) - b) / a), ((((-0.25 * a) * pow(((2.0 * c) / fabs(b)), 2.0)) - c) / fabs(b)));
	double tmp;
	if (t_0 <= -((double) INFINITY)) {
		tmp = t_1;
	} else if (t_0 <= -4e-240) {
		tmp = t_0;
	} else if (t_0 <= 0.0) {
		tmp = -c / b;
	} else if (t_0 <= 4e+216) {
		tmp = (sqrt(fma(b, b, (-4.0 * (a * c)))) / (a + a)) - (b / (a + a));
	} else {
		tmp = t_1;
	}
	return tmp;
}

Julia

function code(a, b, c)
	return Float64(Float64(Float64(-b) + sqrt(Float64(Float64(b * b) - Float64(4.0 * Float64(a * c))))) / Float64(2.0 * a))
end

↓

function code(a, b, c)
	t_0 = Float64(Float64(Float64(-b) + sqrt(Float64(Float64(b * b) - Float64(4.0 * Float64(a * c))))) / Float64(2.0 * a))
	t_1 = fma(0.5, Float64(Float64(abs(b) - b) / a), Float64(Float64(Float64(Float64(-0.25 * a) * (Float64(Float64(2.0 * c) / abs(b)) ^ 2.0)) - c) / abs(b)))
	tmp = 0.0
	if (t_0 <= Float64(-Inf))
		tmp = t_1;
	elseif (t_0 <= -4e-240)
		tmp = t_0;
	elseif (t_0 <= 0.0)
		tmp = Float64(Float64(-c) / b);
	elseif (t_0 <= 4e+216)
		tmp = Float64(Float64(sqrt(fma(b, b, Float64(-4.0 * Float64(a * c)))) / Float64(a + a)) - Float64(b / Float64(a + a)));
	else
		tmp = t_1;
	end
	return tmp
end

Wolfram

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

↓

code[a_, b_, c_] := Block[{t$95$0 = N[(N[((-b) + N[Sqrt[N[(N[(b * b), $MachinePrecision] - N[(4.0 * N[(a * c), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision] / N[(2.0 * a), $MachinePrecision]), $MachinePrecision]}, Block[{t$95$1 = N[(0.5 * N[(N[(N[Abs[b], $MachinePrecision] - b), $MachinePrecision] / a), $MachinePrecision] + N[(N[(N[(N[(-0.25 * a), $MachinePrecision] * N[Power[N[(N[(2.0 * c), $MachinePrecision] / N[Abs[b], $MachinePrecision]), $MachinePrecision], 2.0], $MachinePrecision]), $MachinePrecision] - c), $MachinePrecision] / N[Abs[b], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]}, If[LessEqual[t$95$0, (-Infinity)], t$95$1, If[LessEqual[t$95$0, -4e-240], t$95$0, If[LessEqual[t$95$0, 0.0], N[((-c) / b), $MachinePrecision], If[LessEqual[t$95$0, 4e+216], N[(N[(N[Sqrt[N[(b * b + N[(-4.0 * N[(a * c), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision] / N[(a + a), $MachinePrecision]), $MachinePrecision] - N[(b / N[(a + a), $MachinePrecision]), $MachinePrecision]), $MachinePrecision], t$95$1]]]]]]

Target

Original	34.4
Target	21.3
Herbie	8.0

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

Derivation

1. Split input into 4 regimes

2. if (/.f64 (+.f64 (neg.f64 b) (sqrt.f64 (-.f64 (*.f64 b b) (*.f64 4 (*.f64 a c))))) (*.f64 2 a)) < -inf.0 or 4.0000000000000001e216 < (/.f64 (+.f64 (neg.f64 b) (sqrt.f64 (-.f64 (*.f64 b b) (*.f64 4 (*.f64 a c))))) (*.f64 2 a))

   1. Initial program 59.8

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

   2. Taylor expanded in a around 0 61.7

      \[\leadsto \color{blue}{0.5 \cdot \frac{\sqrt{{b}^{2}} - b}{a} + \left(-1 \cdot \frac{c}{\sqrt{{b}^{2}}} + -0.25 \cdot \frac{a \cdot {\left(-2 \cdot \frac{c}{\sqrt{{b}^{2}}}\right)}^{2}}{\sqrt{{b}^{2}}}\right)} \]

   3. Simplified 15.0

      \[\leadsto \color{blue}{\mathsf{fma}\left(0.5, \frac{\left|b\right| - b}{a}, \frac{\left(-c\right) + -0.25 \cdot \left(a \cdot {\left(\frac{-2 \cdot c}{\left|b\right|}\right)}^{2}\right)}{\left|b\right|}\right)} \]
      Proof

   4. Applied egg-rr 22.6

      \[\leadsto \mathsf{fma}\left(0.5, \frac{\left|b\right| - b}{a}, \color{blue}{\sqrt[3]{\frac{\left(-0.25 \cdot a\right) \cdot {\left(\frac{2 \cdot c}{\left|b\right|}\right)}^{2} - c}{\left|b\right|}} \cdot \sqrt[3]{{\left(\frac{c - \left(-0.25 \cdot a\right) \cdot {\left(\frac{2 \cdot c}{\left|b\right|}\right)}^{2}}{\left|b\right|}\right)}^{2}}}\right) \]

   5. Simplified 15.0

      \[\leadsto \mathsf{fma}\left(0.5, \frac{\left|b\right| - b}{a}, \color{blue}{\frac{\left(-0.25 \cdot a\right) \cdot {\left(\frac{2 \cdot c}{\left|b\right|}\right)}^{2} - c}{\left|b\right|}}\right) \]
      Proof

   if -inf.0 < (/.f64 (+.f64 (neg.f64 b) (sqrt.f64 (-.f64 (*.f64 b b) (*.f64 4 (*.f64 a c))))) (*.f64 2 a)) < -3.9999999999999999e-240

   1. Initial program 4.1

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

   if -3.9999999999999999e-240 < (/.f64 (+.f64 (neg.f64 b) (sqrt.f64 (-.f64 (*.f64 b b) (*.f64 4 (*.f64 a c))))) (*.f64 2 a)) < -0.0

   1. Initial program 53.0

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

   2. Taylor expanded in b around inf 1.1

      \[\leadsto \color{blue}{-1 \cdot \frac{c}{b}} \]

   3. Simplified 1.1

      \[\leadsto \color{blue}{\frac{-c}{b}} \]
      Proof

   if -0.0 < (/.f64 (+.f64 (neg.f64 b) (sqrt.f64 (-.f64 (*.f64 b b) (*.f64 4 (*.f64 a c))))) (*.f64 2 a)) < 4.0000000000000001e216

   1. Initial program 5.0

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

   2. Applied egg-rr 5.0

      \[\leadsto \color{blue}{\frac{\sqrt{\mathsf{fma}\left(b, b, -4 \cdot \left(a \cdot c\right)\right)}}{a + a} - \frac{b}{a + a}} \]
3. Recombined 4 regimes into one program.

Alternatives

Alternative 1
Error	10.5
Cost	7688

\[\begin{array}{l} \mathbf{if}\;b \leq -2 \cdot 10^{+145}:\\ \;\;\;\;\frac{-b}{a}\\ \mathbf{elif}\;b \leq 2.8 \cdot 10^{-27}:\\ \;\;\;\;\frac{\left(-b\right) + \sqrt{b \cdot b - 4 \cdot \left(a \cdot c\right)}}{2 \cdot a}\\ \mathbf{else}:\\ \;\;\;\;\frac{-c}{b}\\ \end{array} \]

Alternative 2
Error	14.0
Cost	7368

\[\begin{array}{l} \mathbf{if}\;b \leq -1.05 \cdot 10^{-95}:\\ \;\;\;\;\frac{c}{b} - \frac{b}{a}\\ \mathbf{elif}\;b \leq 4.5 \cdot 10^{-26}:\\ \;\;\;\;\frac{\sqrt{\left(-4 \cdot a\right) \cdot c} - b}{2 \cdot a}\\ \mathbf{else}:\\ \;\;\;\;\frac{-c}{b}\\ \end{array} \]

Alternative 3
Error	39.3
Cost	388

\[\begin{array}{l} \mathbf{if}\;b \leq 1.16 \cdot 10^{-23}:\\ \;\;\;\;\frac{-b}{a}\\ \mathbf{else}:\\ \;\;\;\;\frac{c}{b}\\ \end{array} \]

Alternative 4
Error	22.1
Cost	388

\[\begin{array}{l} \mathbf{if}\;b \leq 4.2 \cdot 10^{-256}:\\ \;\;\;\;\frac{-b}{a}\\ \mathbf{else}:\\ \;\;\;\;\frac{-c}{b}\\ \end{array} \]

Alternative 5
Error	56.8
Cost	192

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

Reproduce

herbie shell --seed 2023033 
(FPCore (a b c)
  :name "quadp (p42, positive)"
  :precision binary64

  :herbie-target
  (if (< b 0.0) (/ (+ (- b) (sqrt (- (* b b) (* 4.0 (* a c))))) (* 2.0 a)) (/ c (* a (/ (- (- b) (sqrt (- (* b b) (* 4.0 (* a c))))) (* 2.0 a)))))

  (/ (+ (- b) (sqrt (- (* b b) (* 4.0 (* a c))))) (* 2.0 a)))