Complex division, imag part

Average Error: 26.5 → 12.4

Time: 9.3s

Precision: binary64

Math

\[\frac{b \cdot c - a \cdot d}{c \cdot c + d \cdot d} \]

↓

\[\begin{array}{l} t_0 := {\left(\mathsf{hypot}\left(d, c\right)\right)}^{2}\\ t_1 := \mathsf{fma}\left(\frac{b}{d}, \frac{c}{d}, \frac{-a}{d}\right)\\ t_2 := \frac{b}{c} - \frac{a}{c} \cdot \frac{d}{c}\\ \mathbf{if}\;c \leq -5.174874392946877 \cdot 10^{+44}:\\ \;\;\;\;t_2\\ \mathbf{elif}\;c \leq -1.6974732352827345 \cdot 10^{+24}:\\ \;\;\;\;t_1\\ \mathbf{elif}\;c \leq -1.7029803518844395 \cdot 10^{-105}:\\ \;\;\;\;\left(a \cdot d - c \cdot b\right) \cdot \frac{1}{-t_0}\\ \mathbf{elif}\;c \leq 1.700724934216277 \cdot 10^{-110}:\\ \;\;\;\;\frac{b \cdot \frac{c}{d} - a}{d}\\ \mathbf{elif}\;c \leq 1.948418831344656 \cdot 10^{+28}:\\ \;\;\;\;\frac{c \cdot b - a \cdot d}{t_0}\\ \mathbf{elif}\;c \leq 6.489345394500517 \cdot 10^{+62}:\\ \;\;\;\;t_1\\ \mathbf{else}:\\ \;\;\;\;t_2\\ \end{array} \]

FPCore

(FPCore (a b c d)
 :precision binary64
 (/ (- (* b c) (* a d)) (+ (* c c) (* d d))))

↓

(FPCore (a b c d)
 :precision binary64
 (let* ((t_0 (pow (hypot d c) 2.0))
        (t_1 (fma (/ b d) (/ c d) (/ (- a) d)))
        (t_2 (- (/ b c) (* (/ a c) (/ d c)))))
   (if (<= c -5.174874392946877e+44)
     t_2
     (if (<= c -1.6974732352827345e+24)
       t_1
       (if (<= c -1.7029803518844395e-105)
         (* (- (* a d) (* c b)) (/ 1.0 (- t_0)))
         (if (<= c 1.700724934216277e-110)
           (/ (- (* b (/ c d)) a) d)
           (if (<= c 1.948418831344656e+28)
             (/ (- (* c b) (* a d)) t_0)
             (if (<= c 6.489345394500517e+62) t_1 t_2))))))))

C

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

↓

double code(double a, double b, double c, double d) {
	double t_0 = pow(hypot(d, c), 2.0);
	double t_1 = fma((b / d), (c / d), (-a / d));
	double t_2 = (b / c) - ((a / c) * (d / c));
	double tmp;
	if (c <= -5.174874392946877e+44) {
		tmp = t_2;
	} else if (c <= -1.6974732352827345e+24) {
		tmp = t_1;
	} else if (c <= -1.7029803518844395e-105) {
		tmp = ((a * d) - (c * b)) * (1.0 / -t_0);
	} else if (c <= 1.700724934216277e-110) {
		tmp = ((b * (c / d)) - a) / d;
	} else if (c <= 1.948418831344656e+28) {
		tmp = ((c * b) - (a * d)) / t_0;
	} else if (c <= 6.489345394500517e+62) {
		tmp = t_1;
	} else {
		tmp = t_2;
	}
	return tmp;
}

Julia

function code(a, b, c, d)
	return Float64(Float64(Float64(b * c) - Float64(a * d)) / Float64(Float64(c * c) + Float64(d * d)))
end

↓

function code(a, b, c, d)
	t_0 = hypot(d, c) ^ 2.0
	t_1 = fma(Float64(b / d), Float64(c / d), Float64(Float64(-a) / d))
	t_2 = Float64(Float64(b / c) - Float64(Float64(a / c) * Float64(d / c)))
	tmp = 0.0
	if (c <= -5.174874392946877e+44)
		tmp = t_2;
	elseif (c <= -1.6974732352827345e+24)
		tmp = t_1;
	elseif (c <= -1.7029803518844395e-105)
		tmp = Float64(Float64(Float64(a * d) - Float64(c * b)) * Float64(1.0 / Float64(-t_0)));
	elseif (c <= 1.700724934216277e-110)
		tmp = Float64(Float64(Float64(b * Float64(c / d)) - a) / d);
	elseif (c <= 1.948418831344656e+28)
		tmp = Float64(Float64(Float64(c * b) - Float64(a * d)) / t_0);
	elseif (c <= 6.489345394500517e+62)
		tmp = t_1;
	else
		tmp = t_2;
	end
	return tmp
end

Wolfram

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

↓

code[a_, b_, c_, d_] := Block[{t$95$0 = N[Power[N[Sqrt[d ^ 2 + c ^ 2], $MachinePrecision], 2.0], $MachinePrecision]}, Block[{t$95$1 = N[(N[(b / d), $MachinePrecision] * N[(c / d), $MachinePrecision] + N[((-a) / d), $MachinePrecision]), $MachinePrecision]}, Block[{t$95$2 = N[(N[(b / c), $MachinePrecision] - N[(N[(a / c), $MachinePrecision] * N[(d / c), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]}, If[LessEqual[c, -5.174874392946877e+44], t$95$2, If[LessEqual[c, -1.6974732352827345e+24], t$95$1, If[LessEqual[c, -1.7029803518844395e-105], N[(N[(N[(a * d), $MachinePrecision] - N[(c * b), $MachinePrecision]), $MachinePrecision] * N[(1.0 / (-t$95$0)), $MachinePrecision]), $MachinePrecision], If[LessEqual[c, 1.700724934216277e-110], N[(N[(N[(b * N[(c / d), $MachinePrecision]), $MachinePrecision] - a), $MachinePrecision] / d), $MachinePrecision], If[LessEqual[c, 1.948418831344656e+28], N[(N[(N[(c * b), $MachinePrecision] - N[(a * d), $MachinePrecision]), $MachinePrecision] / t$95$0), $MachinePrecision], If[LessEqual[c, 6.489345394500517e+62], t$95$1, t$95$2]]]]]]]]]

Error

Bits error versus a

Bits error versus b

Bits error versus c

Bits error versus d

Target

Original	26.5
Target	0.4
Herbie	12.4

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

Derivation

1. Split input into 5 regimes

2. if c < -5.1748743929468772e44 or 6.48934539450051737e62 < c

   1. Initial program 36.6

      \[\frac{b \cdot c - a \cdot d}{c \cdot c + d \cdot d} \]

   2. Simplified 36.6

      \[\leadsto \color{blue}{\frac{\mathsf{fma}\left(d, -a, b \cdot c\right)}{\mathsf{fma}\left(c, c, d \cdot d\right)}} \]

   3. Applied egg-rr 36.6

      \[\leadsto \color{blue}{\left(d \cdot a - b \cdot c\right) \cdot \frac{1}{-{\left(\mathsf{hypot}\left(d, c\right)\right)}^{2}}} \]

   4. Applied egg-rr 36.8

      \[\leadsto \color{blue}{{\left(\sqrt[3]{d \cdot a - b \cdot c}\right)}^{3}} \cdot \frac{1}{-{\left(\mathsf{hypot}\left(d, c\right)\right)}^{2}} \]

   5. Taylor expanded in d around 0 17.7

      \[\leadsto \color{blue}{\frac{b}{c} - \frac{a \cdot d}{{c}^{2}}} \]

   6. Simplified 12.6

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

   if -5.1748743929468772e44 < c < -1.69747323528273453e24 or 1.9484188313446559e28 < c < 6.48934539450051737e62

   1. Initial program 20.7

      \[\frac{b \cdot c - a \cdot d}{c \cdot c + d \cdot d} \]

   2. Simplified 20.7

      \[\leadsto \color{blue}{\frac{\mathsf{fma}\left(d, -a, b \cdot c\right)}{\mathsf{fma}\left(c, c, d \cdot d\right)}} \]

   3. Applied egg-rr 20.7

      \[\leadsto \color{blue}{\left(d \cdot a - b \cdot c\right) \cdot \frac{1}{-{\left(\mathsf{hypot}\left(d, c\right)\right)}^{2}}} \]

   4. Applied egg-rr 21.3

      \[\leadsto \color{blue}{{\left(\sqrt[3]{d \cdot a - b \cdot c}\right)}^{3}} \cdot \frac{1}{-{\left(\mathsf{hypot}\left(d, c\right)\right)}^{2}} \]

   5. Taylor expanded in d around inf 37.3

      \[\leadsto \color{blue}{\frac{c \cdot b}{{d}^{2}} - \frac{a}{d}} \]

   6. Simplified 31.7

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

   7. Applied egg-rr 31.7

      \[\leadsto \color{blue}{\mathsf{fma}\left(\frac{\frac{b}{d}}{1}, \frac{c}{d}, \frac{-a}{d} \cdot 1\right)} \]

   if -1.69747323528273453e24 < c < -1.70298035188443947e-105

   1. Initial program 14.1

      \[\frac{b \cdot c - a \cdot d}{c \cdot c + d \cdot d} \]

   2. Simplified 14.1

      \[\leadsto \color{blue}{\frac{\mathsf{fma}\left(d, -a, b \cdot c\right)}{\mathsf{fma}\left(c, c, d \cdot d\right)}} \]

   3. Applied egg-rr 14.1

      \[\leadsto \color{blue}{\left(d \cdot a - b \cdot c\right) \cdot \frac{1}{-{\left(\mathsf{hypot}\left(d, c\right)\right)}^{2}}} \]

   if -1.70298035188443947e-105 < c < 1.70072493421627706e-110

   1. Initial program 22.2

      \[\frac{b \cdot c - a \cdot d}{c \cdot c + d \cdot d} \]

   2. Simplified 22.2

      \[\leadsto \color{blue}{\frac{\mathsf{fma}\left(d, -a, b \cdot c\right)}{\mathsf{fma}\left(c, c, d \cdot d\right)}} \]

   3. Applied egg-rr 22.7

      \[\leadsto \color{blue}{\left(d \cdot a - b \cdot c\right) \cdot \frac{1}{-{\left(\mathsf{hypot}\left(d, c\right)\right)}^{2}}} \]

   4. Applied egg-rr 23.2

      \[\leadsto \color{blue}{{\left(\sqrt[3]{d \cdot a - b \cdot c}\right)}^{3}} \cdot \frac{1}{-{\left(\mathsf{hypot}\left(d, c\right)\right)}^{2}} \]

   5. Taylor expanded in d around inf 11.4

      \[\leadsto \color{blue}{\frac{c \cdot b}{{d}^{2}} - \frac{a}{d}} \]

   6. Simplified 10.1

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

   7. Applied egg-rr 8.0

      \[\leadsto \color{blue}{{\left(\frac{b \cdot \frac{c}{d} - a}{d}\right)}^{1}} \]

   if 1.70072493421627706e-110 < c < 1.9484188313446559e28

   1. Initial program 13.8

      \[\frac{b \cdot c - a \cdot d}{c \cdot c + d \cdot d} \]

   2. Simplified 13.8

      \[\leadsto \color{blue}{\frac{\mathsf{fma}\left(d, -a, b \cdot c\right)}{\mathsf{fma}\left(c, c, d \cdot d\right)}} \]

   3. Applied egg-rr 13.8

      \[\leadsto \color{blue}{\left(d \cdot a - b \cdot c\right) \cdot \frac{1}{-{\left(\mathsf{hypot}\left(d, c\right)\right)}^{2}}} \]

   4. Applied egg-rr 13.8

      \[\leadsto \color{blue}{\frac{\left(d \cdot a - b \cdot c\right) \cdot -1}{{\left(\mathsf{hypot}\left(d, c\right)\right)}^{2}}} \]
3. Recombined 5 regimes into one program.

4. Final simplification 12.4

   \[\leadsto \begin{array}{l} \mathbf{if}\;c \leq -5.174874392946877 \cdot 10^{+44}:\\ \;\;\;\;\frac{b}{c} - \frac{a}{c} \cdot \frac{d}{c}\\ \mathbf{elif}\;c \leq -1.6974732352827345 \cdot 10^{+24}:\\ \;\;\;\;\mathsf{fma}\left(\frac{b}{d}, \frac{c}{d}, \frac{-a}{d}\right)\\ \mathbf{elif}\;c \leq -1.7029803518844395 \cdot 10^{-105}:\\ \;\;\;\;\left(a \cdot d - c \cdot b\right) \cdot \frac{1}{-{\left(\mathsf{hypot}\left(d, c\right)\right)}^{2}}\\ \mathbf{elif}\;c \leq 1.700724934216277 \cdot 10^{-110}:\\ \;\;\;\;\frac{b \cdot \frac{c}{d} - a}{d}\\ \mathbf{elif}\;c \leq 1.948418831344656 \cdot 10^{+28}:\\ \;\;\;\;\frac{c \cdot b - a \cdot d}{{\left(\mathsf{hypot}\left(d, c\right)\right)}^{2}}\\ \mathbf{elif}\;c \leq 6.489345394500517 \cdot 10^{+62}:\\ \;\;\;\;\mathsf{fma}\left(\frac{b}{d}, \frac{c}{d}, \frac{-a}{d}\right)\\ \mathbf{else}:\\ \;\;\;\;\frac{b}{c} - \frac{a}{c} \cdot \frac{d}{c}\\ \end{array} \]

Reproduce

herbie shell --seed 2022150 
(FPCore (a b c d)
  :name "Complex division, imag part"
  :precision binary64

  :herbie-target
  (if (< (fabs d) (fabs c)) (/ (- b (* a (/ d c))) (+ c (* d (/ d c)))) (/ (+ (- a) (* b (/ c d))) (+ d (* c (/ c d)))))

  (/ (- (* b c) (* a d)) (+ (* c c) (* d d))))