Diagrams.TwoD.Apollonian:initialConfig from diagrams-contrib-1.3.0.5, A

Percentage Accurate: 68.5% → 97.1%
Time: 6.9s
Alternatives: 9
Speedup: 1.3×

Specification

?
\[\begin{array}{l} \\ \frac{\left(x \cdot x + y \cdot y\right) - z \cdot z}{y \cdot 2} \end{array} \]
(FPCore (x y z)
 :precision binary64
 (/ (- (+ (* x x) (* y y)) (* z z)) (* y 2.0)))
double code(double x, double y, double z) {
	return (((x * x) + (y * y)) - (z * z)) / (y * 2.0);
}
real(8) function code(x, y, z)
    real(8), intent (in) :: x
    real(8), intent (in) :: y
    real(8), intent (in) :: z
    code = (((x * x) + (y * y)) - (z * z)) / (y * 2.0d0)
end function
public static double code(double x, double y, double z) {
	return (((x * x) + (y * y)) - (z * z)) / (y * 2.0);
}
def code(x, y, z):
	return (((x * x) + (y * y)) - (z * z)) / (y * 2.0)
function code(x, y, z)
	return Float64(Float64(Float64(Float64(x * x) + Float64(y * y)) - Float64(z * z)) / Float64(y * 2.0))
end
function tmp = code(x, y, z)
	tmp = (((x * x) + (y * y)) - (z * z)) / (y * 2.0);
end
code[x_, y_, z_] := N[(N[(N[(N[(x * x), $MachinePrecision] + N[(y * y), $MachinePrecision]), $MachinePrecision] - N[(z * z), $MachinePrecision]), $MachinePrecision] / N[(y * 2.0), $MachinePrecision]), $MachinePrecision]
\begin{array}{l}

\\
\frac{\left(x \cdot x + y \cdot y\right) - z \cdot z}{y \cdot 2}
\end{array}

Sampling outcomes in binary64 precision:

Local Percentage Accuracy vs ?

The average percentage accuracy by input value. Horizontal axis shows value of an input variable; the variable is choosen in the title. Vertical axis is accuracy; higher is better. Red represent the original program, while blue represents Herbie's suggestion. These can be toggled with buttons below the plot. The line is an average while dots represent individual samples.

Accuracy vs Speed?

Herbie found 9 alternatives:

AlternativeAccuracySpeedup
The accuracy (vertical axis) and speed (horizontal axis) of each alternatives. Up and to the right is better. The red square shows the initial program, and each blue circle shows an alternative.The line shows the best available speed-accuracy tradeoffs.

Initial Program: 68.5% accurate, 1.0× speedup?

\[\begin{array}{l} \\ \frac{\left(x \cdot x + y \cdot y\right) - z \cdot z}{y \cdot 2} \end{array} \]
(FPCore (x y z)
 :precision binary64
 (/ (- (+ (* x x) (* y y)) (* z z)) (* y 2.0)))
double code(double x, double y, double z) {
	return (((x * x) + (y * y)) - (z * z)) / (y * 2.0);
}
real(8) function code(x, y, z)
    real(8), intent (in) :: x
    real(8), intent (in) :: y
    real(8), intent (in) :: z
    code = (((x * x) + (y * y)) - (z * z)) / (y * 2.0d0)
end function
public static double code(double x, double y, double z) {
	return (((x * x) + (y * y)) - (z * z)) / (y * 2.0);
}
def code(x, y, z):
	return (((x * x) + (y * y)) - (z * z)) / (y * 2.0)
function code(x, y, z)
	return Float64(Float64(Float64(Float64(x * x) + Float64(y * y)) - Float64(z * z)) / Float64(y * 2.0))
end
function tmp = code(x, y, z)
	tmp = (((x * x) + (y * y)) - (z * z)) / (y * 2.0);
end
code[x_, y_, z_] := N[(N[(N[(N[(x * x), $MachinePrecision] + N[(y * y), $MachinePrecision]), $MachinePrecision] - N[(z * z), $MachinePrecision]), $MachinePrecision] / N[(y * 2.0), $MachinePrecision]), $MachinePrecision]
\begin{array}{l}

\\
\frac{\left(x \cdot x + y \cdot y\right) - z \cdot z}{y \cdot 2}
\end{array}

Alternative 1: 97.1% accurate, 0.3× speedup?

\[\begin{array}{l} x_m = \left|x\right| \\ y\_m = \left|y\right| \\ y\_s = \mathsf{copysign}\left(1, y\right) \\ \begin{array}{l} t_0 := \frac{\left(x\_m \cdot x\_m + y\_m \cdot y\_m\right) - z \cdot z}{y\_m \cdot 2}\\ y\_s \cdot \begin{array}{l} \mathbf{if}\;t\_0 \leq 0:\\ \;\;\;\;\left(-0.5 \cdot \frac{z}{y\_m}\right) \cdot z\\ \mathbf{elif}\;t\_0 \leq \infty:\\ \;\;\;\;\mathsf{fma}\left(\frac{x\_m}{y\_m}, x\_m, y\_m\right) \cdot 0.5\\ \mathbf{else}:\\ \;\;\;\;\mathsf{fma}\left(x\_m - z, \frac{z}{y\_m}, y\_m\right) \cdot 0.5\\ \end{array} \end{array} \end{array} \]
x_m = (fabs.f64 x)
y\_m = (fabs.f64 y)
y\_s = (copysign.f64 #s(literal 1 binary64) y)
(FPCore (y_s x_m y_m z)
 :precision binary64
 (let* ((t_0 (/ (- (+ (* x_m x_m) (* y_m y_m)) (* z z)) (* y_m 2.0))))
   (*
    y_s
    (if (<= t_0 0.0)
      (* (* -0.5 (/ z y_m)) z)
      (if (<= t_0 INFINITY)
        (* (fma (/ x_m y_m) x_m y_m) 0.5)
        (* (fma (- x_m z) (/ z y_m) y_m) 0.5))))))
x_m = fabs(x);
y\_m = fabs(y);
y\_s = copysign(1.0, y);
double code(double y_s, double x_m, double y_m, double z) {
	double t_0 = (((x_m * x_m) + (y_m * y_m)) - (z * z)) / (y_m * 2.0);
	double tmp;
	if (t_0 <= 0.0) {
		tmp = (-0.5 * (z / y_m)) * z;
	} else if (t_0 <= ((double) INFINITY)) {
		tmp = fma((x_m / y_m), x_m, y_m) * 0.5;
	} else {
		tmp = fma((x_m - z), (z / y_m), y_m) * 0.5;
	}
	return y_s * tmp;
}
x_m = abs(x)
y\_m = abs(y)
y\_s = copysign(1.0, y)
function code(y_s, x_m, y_m, z)
	t_0 = Float64(Float64(Float64(Float64(x_m * x_m) + Float64(y_m * y_m)) - Float64(z * z)) / Float64(y_m * 2.0))
	tmp = 0.0
	if (t_0 <= 0.0)
		tmp = Float64(Float64(-0.5 * Float64(z / y_m)) * z);
	elseif (t_0 <= Inf)
		tmp = Float64(fma(Float64(x_m / y_m), x_m, y_m) * 0.5);
	else
		tmp = Float64(fma(Float64(x_m - z), Float64(z / y_m), y_m) * 0.5);
	end
	return Float64(y_s * tmp)
end
x_m = N[Abs[x], $MachinePrecision]
y\_m = N[Abs[y], $MachinePrecision]
y\_s = N[With[{TMP1 = Abs[1.0], TMP2 = Sign[y]}, TMP1 * If[TMP2 == 0, 1, TMP2]], $MachinePrecision]
code[y$95$s_, x$95$m_, y$95$m_, z_] := Block[{t$95$0 = N[(N[(N[(N[(x$95$m * x$95$m), $MachinePrecision] + N[(y$95$m * y$95$m), $MachinePrecision]), $MachinePrecision] - N[(z * z), $MachinePrecision]), $MachinePrecision] / N[(y$95$m * 2.0), $MachinePrecision]), $MachinePrecision]}, N[(y$95$s * If[LessEqual[t$95$0, 0.0], N[(N[(-0.5 * N[(z / y$95$m), $MachinePrecision]), $MachinePrecision] * z), $MachinePrecision], If[LessEqual[t$95$0, Infinity], N[(N[(N[(x$95$m / y$95$m), $MachinePrecision] * x$95$m + y$95$m), $MachinePrecision] * 0.5), $MachinePrecision], N[(N[(N[(x$95$m - z), $MachinePrecision] * N[(z / y$95$m), $MachinePrecision] + y$95$m), $MachinePrecision] * 0.5), $MachinePrecision]]]), $MachinePrecision]]
\begin{array}{l}
x_m = \left|x\right|
\\
y\_m = \left|y\right|
\\
y\_s = \mathsf{copysign}\left(1, y\right)

\\
\begin{array}{l}
t_0 := \frac{\left(x\_m \cdot x\_m + y\_m \cdot y\_m\right) - z \cdot z}{y\_m \cdot 2}\\
y\_s \cdot \begin{array}{l}
\mathbf{if}\;t\_0 \leq 0:\\
\;\;\;\;\left(-0.5 \cdot \frac{z}{y\_m}\right) \cdot z\\

\mathbf{elif}\;t\_0 \leq \infty:\\
\;\;\;\;\mathsf{fma}\left(\frac{x\_m}{y\_m}, x\_m, y\_m\right) \cdot 0.5\\

\mathbf{else}:\\
\;\;\;\;\mathsf{fma}\left(x\_m - z, \frac{z}{y\_m}, y\_m\right) \cdot 0.5\\


\end{array}
\end{array}
\end{array}
Derivation
  1. Split input into 3 regimes
  2. if (/.f64 (-.f64 (+.f64 (*.f64 x x) (*.f64 y y)) (*.f64 z z)) (*.f64 y #s(literal 2 binary64))) < 0.0

    1. Initial program 83.8%

      \[\frac{\left(x \cdot x + y \cdot y\right) - z \cdot z}{y \cdot 2} \]
    2. Add Preprocessing
    3. Taylor expanded in x around 0

      \[\leadsto \color{blue}{\frac{1}{2} \cdot \frac{{y}^{2} - {z}^{2}}{y} + \frac{1}{2} \cdot \frac{{x}^{2}}{y}} \]
    4. Step-by-step derivation
      1. distribute-lft-outN/A

        \[\leadsto \color{blue}{\frac{1}{2} \cdot \left(\frac{{y}^{2} - {z}^{2}}{y} + \frac{{x}^{2}}{y}\right)} \]
      2. *-commutativeN/A

        \[\leadsto \color{blue}{\left(\frac{{y}^{2} - {z}^{2}}{y} + \frac{{x}^{2}}{y}\right) \cdot \frac{1}{2}} \]
      3. lower-*.f64N/A

        \[\leadsto \color{blue}{\left(\frac{{y}^{2} - {z}^{2}}{y} + \frac{{x}^{2}}{y}\right) \cdot \frac{1}{2}} \]
    5. Applied rewrites99.3%

      \[\leadsto \color{blue}{\mathsf{fma}\left(z + x, \frac{x - z}{y}, y\right) \cdot 0.5} \]
    6. Step-by-step derivation
      1. Applied rewrites99.3%

        \[\leadsto \left(\frac{x - z}{y} \cdot \left(z + x\right) + y\right) \cdot 0.5 \]
      2. Taylor expanded in z around inf

        \[\leadsto {z}^{2} \cdot \color{blue}{\left(\frac{1}{2} \cdot \frac{-1 \cdot \frac{x}{y} + \frac{x}{y}}{z} - \frac{1}{2} \cdot \frac{1}{y}\right)} \]
      3. Step-by-step derivation
        1. Applied rewrites32.7%

          \[\leadsto \left(\left(0 - \frac{0.5}{y}\right) \cdot z\right) \cdot \color{blue}{z} \]
        2. Taylor expanded in y around 0

          \[\leadsto \left(\frac{-1}{2} \cdot \frac{z}{y}\right) \cdot z \]
        3. Step-by-step derivation
          1. Applied rewrites32.7%

            \[\leadsto \left(-0.5 \cdot \frac{z}{y}\right) \cdot z \]

          if 0.0 < (/.f64 (-.f64 (+.f64 (*.f64 x x) (*.f64 y y)) (*.f64 z z)) (*.f64 y #s(literal 2 binary64))) < +inf.0

          1. Initial program 76.9%

            \[\frac{\left(x \cdot x + y \cdot y\right) - z \cdot z}{y \cdot 2} \]
          2. Add Preprocessing
          3. Taylor expanded in z around 0

            \[\leadsto \color{blue}{\frac{1}{2} \cdot \frac{{x}^{2} + {y}^{2}}{y}} \]
          4. Step-by-step derivation
            1. div-addN/A

              \[\leadsto \frac{1}{2} \cdot \color{blue}{\left(\frac{{x}^{2}}{y} + \frac{{y}^{2}}{y}\right)} \]
            2. distribute-lft-inN/A

              \[\leadsto \color{blue}{\frac{1}{2} \cdot \frac{{x}^{2}}{y} + \frac{1}{2} \cdot \frac{{y}^{2}}{y}} \]
            3. unpow2N/A

              \[\leadsto \frac{1}{2} \cdot \frac{{x}^{2}}{y} + \frac{1}{2} \cdot \frac{\color{blue}{y \cdot y}}{y} \]
            4. associate-/l*N/A

              \[\leadsto \frac{1}{2} \cdot \frac{{x}^{2}}{y} + \frac{1}{2} \cdot \color{blue}{\left(y \cdot \frac{y}{y}\right)} \]
            5. *-inversesN/A

              \[\leadsto \frac{1}{2} \cdot \frac{{x}^{2}}{y} + \frac{1}{2} \cdot \left(y \cdot \color{blue}{1}\right) \]
            6. *-rgt-identityN/A

              \[\leadsto \frac{1}{2} \cdot \frac{{x}^{2}}{y} + \frac{1}{2} \cdot \color{blue}{y} \]
            7. distribute-lft-inN/A

              \[\leadsto \color{blue}{\frac{1}{2} \cdot \left(\frac{{x}^{2}}{y} + y\right)} \]
            8. +-commutativeN/A

              \[\leadsto \frac{1}{2} \cdot \color{blue}{\left(y + \frac{{x}^{2}}{y}\right)} \]
            9. *-commutativeN/A

              \[\leadsto \color{blue}{\left(y + \frac{{x}^{2}}{y}\right) \cdot \frac{1}{2}} \]
            10. lower-*.f64N/A

              \[\leadsto \color{blue}{\left(y + \frac{{x}^{2}}{y}\right) \cdot \frac{1}{2}} \]
            11. +-commutativeN/A

              \[\leadsto \color{blue}{\left(\frac{{x}^{2}}{y} + y\right)} \cdot \frac{1}{2} \]
            12. unpow2N/A

              \[\leadsto \left(\frac{\color{blue}{x \cdot x}}{y} + y\right) \cdot \frac{1}{2} \]
            13. associate-/l*N/A

              \[\leadsto \left(\color{blue}{x \cdot \frac{x}{y}} + y\right) \cdot \frac{1}{2} \]
            14. *-commutativeN/A

              \[\leadsto \left(\color{blue}{\frac{x}{y} \cdot x} + y\right) \cdot \frac{1}{2} \]
            15. lower-fma.f64N/A

              \[\leadsto \color{blue}{\mathsf{fma}\left(\frac{x}{y}, x, y\right)} \cdot \frac{1}{2} \]
            16. lower-/.f6470.1

              \[\leadsto \mathsf{fma}\left(\color{blue}{\frac{x}{y}}, x, y\right) \cdot 0.5 \]
          5. Applied rewrites70.1%

            \[\leadsto \color{blue}{\mathsf{fma}\left(\frac{x}{y}, x, y\right) \cdot 0.5} \]

          if +inf.0 < (/.f64 (-.f64 (+.f64 (*.f64 x x) (*.f64 y y)) (*.f64 z z)) (*.f64 y #s(literal 2 binary64)))

          1. Initial program 0.0%

            \[\frac{\left(x \cdot x + y \cdot y\right) - z \cdot z}{y \cdot 2} \]
          2. Add Preprocessing
          3. Taylor expanded in x around 0

            \[\leadsto \color{blue}{\frac{1}{2} \cdot \frac{{y}^{2} - {z}^{2}}{y} + \frac{1}{2} \cdot \frac{{x}^{2}}{y}} \]
          4. Step-by-step derivation
            1. distribute-lft-outN/A

              \[\leadsto \color{blue}{\frac{1}{2} \cdot \left(\frac{{y}^{2} - {z}^{2}}{y} + \frac{{x}^{2}}{y}\right)} \]
            2. *-commutativeN/A

              \[\leadsto \color{blue}{\left(\frac{{y}^{2} - {z}^{2}}{y} + \frac{{x}^{2}}{y}\right) \cdot \frac{1}{2}} \]
            3. lower-*.f64N/A

              \[\leadsto \color{blue}{\left(\frac{{y}^{2} - {z}^{2}}{y} + \frac{{x}^{2}}{y}\right) \cdot \frac{1}{2}} \]
          5. Applied rewrites99.9%

            \[\leadsto \color{blue}{\mathsf{fma}\left(z + x, \frac{x - z}{y}, y\right) \cdot 0.5} \]
          6. Taylor expanded in x around 0

            \[\leadsto \mathsf{fma}\left(z + x, -1 \cdot \frac{z}{y}, y\right) \cdot \frac{1}{2} \]
          7. Step-by-step derivation
            1. Applied rewrites92.7%

              \[\leadsto \mathsf{fma}\left(z + x, \frac{-z}{y}, y\right) \cdot 0.5 \]
            2. Step-by-step derivation
              1. Applied rewrites96.3%

                \[\leadsto \mathsf{fma}\left(x - z, \frac{z}{y}, y\right) \cdot 0.5 \]
            3. Recombined 3 regimes into one program.
            4. Add Preprocessing

            Alternative 2: 72.0% accurate, 0.3× speedup?

            \[\begin{array}{l} x_m = \left|x\right| \\ y\_m = \left|y\right| \\ y\_s = \mathsf{copysign}\left(1, y\right) \\ \begin{array}{l} t_0 := \left(-0.5 \cdot \frac{z}{y\_m}\right) \cdot z\\ t_1 := \frac{\left(x\_m \cdot x\_m + y\_m \cdot y\_m\right) - z \cdot z}{y\_m \cdot 2}\\ y\_s \cdot \begin{array}{l} \mathbf{if}\;t\_1 \leq 0:\\ \;\;\;\;t\_0\\ \mathbf{elif}\;t\_1 \leq 2 \cdot 10^{+151}:\\ \;\;\;\;0.5 \cdot y\_m\\ \mathbf{elif}\;t\_1 \leq \infty:\\ \;\;\;\;\left(\frac{x\_m}{y\_m} \cdot x\_m\right) \cdot 0.5\\ \mathbf{else}:\\ \;\;\;\;t\_0\\ \end{array} \end{array} \end{array} \]
            x_m = (fabs.f64 x)
            y\_m = (fabs.f64 y)
            y\_s = (copysign.f64 #s(literal 1 binary64) y)
            (FPCore (y_s x_m y_m z)
             :precision binary64
             (let* ((t_0 (* (* -0.5 (/ z y_m)) z))
                    (t_1 (/ (- (+ (* x_m x_m) (* y_m y_m)) (* z z)) (* y_m 2.0))))
               (*
                y_s
                (if (<= t_1 0.0)
                  t_0
                  (if (<= t_1 2e+151)
                    (* 0.5 y_m)
                    (if (<= t_1 INFINITY) (* (* (/ x_m y_m) x_m) 0.5) t_0))))))
            x_m = fabs(x);
            y\_m = fabs(y);
            y\_s = copysign(1.0, y);
            double code(double y_s, double x_m, double y_m, double z) {
            	double t_0 = (-0.5 * (z / y_m)) * z;
            	double t_1 = (((x_m * x_m) + (y_m * y_m)) - (z * z)) / (y_m * 2.0);
            	double tmp;
            	if (t_1 <= 0.0) {
            		tmp = t_0;
            	} else if (t_1 <= 2e+151) {
            		tmp = 0.5 * y_m;
            	} else if (t_1 <= ((double) INFINITY)) {
            		tmp = ((x_m / y_m) * x_m) * 0.5;
            	} else {
            		tmp = t_0;
            	}
            	return y_s * tmp;
            }
            
            x_m = Math.abs(x);
            y\_m = Math.abs(y);
            y\_s = Math.copySign(1.0, y);
            public static double code(double y_s, double x_m, double y_m, double z) {
            	double t_0 = (-0.5 * (z / y_m)) * z;
            	double t_1 = (((x_m * x_m) + (y_m * y_m)) - (z * z)) / (y_m * 2.0);
            	double tmp;
            	if (t_1 <= 0.0) {
            		tmp = t_0;
            	} else if (t_1 <= 2e+151) {
            		tmp = 0.5 * y_m;
            	} else if (t_1 <= Double.POSITIVE_INFINITY) {
            		tmp = ((x_m / y_m) * x_m) * 0.5;
            	} else {
            		tmp = t_0;
            	}
            	return y_s * tmp;
            }
            
            x_m = math.fabs(x)
            y\_m = math.fabs(y)
            y\_s = math.copysign(1.0, y)
            def code(y_s, x_m, y_m, z):
            	t_0 = (-0.5 * (z / y_m)) * z
            	t_1 = (((x_m * x_m) + (y_m * y_m)) - (z * z)) / (y_m * 2.0)
            	tmp = 0
            	if t_1 <= 0.0:
            		tmp = t_0
            	elif t_1 <= 2e+151:
            		tmp = 0.5 * y_m
            	elif t_1 <= math.inf:
            		tmp = ((x_m / y_m) * x_m) * 0.5
            	else:
            		tmp = t_0
            	return y_s * tmp
            
            x_m = abs(x)
            y\_m = abs(y)
            y\_s = copysign(1.0, y)
            function code(y_s, x_m, y_m, z)
            	t_0 = Float64(Float64(-0.5 * Float64(z / y_m)) * z)
            	t_1 = Float64(Float64(Float64(Float64(x_m * x_m) + Float64(y_m * y_m)) - Float64(z * z)) / Float64(y_m * 2.0))
            	tmp = 0.0
            	if (t_1 <= 0.0)
            		tmp = t_0;
            	elseif (t_1 <= 2e+151)
            		tmp = Float64(0.5 * y_m);
            	elseif (t_1 <= Inf)
            		tmp = Float64(Float64(Float64(x_m / y_m) * x_m) * 0.5);
            	else
            		tmp = t_0;
            	end
            	return Float64(y_s * tmp)
            end
            
            x_m = abs(x);
            y\_m = abs(y);
            y\_s = sign(y) * abs(1.0);
            function tmp_2 = code(y_s, x_m, y_m, z)
            	t_0 = (-0.5 * (z / y_m)) * z;
            	t_1 = (((x_m * x_m) + (y_m * y_m)) - (z * z)) / (y_m * 2.0);
            	tmp = 0.0;
            	if (t_1 <= 0.0)
            		tmp = t_0;
            	elseif (t_1 <= 2e+151)
            		tmp = 0.5 * y_m;
            	elseif (t_1 <= Inf)
            		tmp = ((x_m / y_m) * x_m) * 0.5;
            	else
            		tmp = t_0;
            	end
            	tmp_2 = y_s * tmp;
            end
            
            x_m = N[Abs[x], $MachinePrecision]
            y\_m = N[Abs[y], $MachinePrecision]
            y\_s = N[With[{TMP1 = Abs[1.0], TMP2 = Sign[y]}, TMP1 * If[TMP2 == 0, 1, TMP2]], $MachinePrecision]
            code[y$95$s_, x$95$m_, y$95$m_, z_] := Block[{t$95$0 = N[(N[(-0.5 * N[(z / y$95$m), $MachinePrecision]), $MachinePrecision] * z), $MachinePrecision]}, Block[{t$95$1 = N[(N[(N[(N[(x$95$m * x$95$m), $MachinePrecision] + N[(y$95$m * y$95$m), $MachinePrecision]), $MachinePrecision] - N[(z * z), $MachinePrecision]), $MachinePrecision] / N[(y$95$m * 2.0), $MachinePrecision]), $MachinePrecision]}, N[(y$95$s * If[LessEqual[t$95$1, 0.0], t$95$0, If[LessEqual[t$95$1, 2e+151], N[(0.5 * y$95$m), $MachinePrecision], If[LessEqual[t$95$1, Infinity], N[(N[(N[(x$95$m / y$95$m), $MachinePrecision] * x$95$m), $MachinePrecision] * 0.5), $MachinePrecision], t$95$0]]]), $MachinePrecision]]]
            
            \begin{array}{l}
            x_m = \left|x\right|
            \\
            y\_m = \left|y\right|
            \\
            y\_s = \mathsf{copysign}\left(1, y\right)
            
            \\
            \begin{array}{l}
            t_0 := \left(-0.5 \cdot \frac{z}{y\_m}\right) \cdot z\\
            t_1 := \frac{\left(x\_m \cdot x\_m + y\_m \cdot y\_m\right) - z \cdot z}{y\_m \cdot 2}\\
            y\_s \cdot \begin{array}{l}
            \mathbf{if}\;t\_1 \leq 0:\\
            \;\;\;\;t\_0\\
            
            \mathbf{elif}\;t\_1 \leq 2 \cdot 10^{+151}:\\
            \;\;\;\;0.5 \cdot y\_m\\
            
            \mathbf{elif}\;t\_1 \leq \infty:\\
            \;\;\;\;\left(\frac{x\_m}{y\_m} \cdot x\_m\right) \cdot 0.5\\
            
            \mathbf{else}:\\
            \;\;\;\;t\_0\\
            
            
            \end{array}
            \end{array}
            \end{array}
            
            Derivation
            1. Split input into 3 regimes
            2. if (/.f64 (-.f64 (+.f64 (*.f64 x x) (*.f64 y y)) (*.f64 z z)) (*.f64 y #s(literal 2 binary64))) < 0.0 or +inf.0 < (/.f64 (-.f64 (+.f64 (*.f64 x x) (*.f64 y y)) (*.f64 z z)) (*.f64 y #s(literal 2 binary64)))

              1. Initial program 68.5%

                \[\frac{\left(x \cdot x + y \cdot y\right) - z \cdot z}{y \cdot 2} \]
              2. Add Preprocessing
              3. Taylor expanded in x around 0

                \[\leadsto \color{blue}{\frac{1}{2} \cdot \frac{{y}^{2} - {z}^{2}}{y} + \frac{1}{2} \cdot \frac{{x}^{2}}{y}} \]
              4. Step-by-step derivation
                1. distribute-lft-outN/A

                  \[\leadsto \color{blue}{\frac{1}{2} \cdot \left(\frac{{y}^{2} - {z}^{2}}{y} + \frac{{x}^{2}}{y}\right)} \]
                2. *-commutativeN/A

                  \[\leadsto \color{blue}{\left(\frac{{y}^{2} - {z}^{2}}{y} + \frac{{x}^{2}}{y}\right) \cdot \frac{1}{2}} \]
                3. lower-*.f64N/A

                  \[\leadsto \color{blue}{\left(\frac{{y}^{2} - {z}^{2}}{y} + \frac{{x}^{2}}{y}\right) \cdot \frac{1}{2}} \]
              5. Applied rewrites99.4%

                \[\leadsto \color{blue}{\mathsf{fma}\left(z + x, \frac{x - z}{y}, y\right) \cdot 0.5} \]
              6. Step-by-step derivation
                1. Applied rewrites99.4%

                  \[\leadsto \left(\frac{x - z}{y} \cdot \left(z + x\right) + y\right) \cdot 0.5 \]
                2. Taylor expanded in z around inf

                  \[\leadsto {z}^{2} \cdot \color{blue}{\left(\frac{1}{2} \cdot \frac{-1 \cdot \frac{x}{y} + \frac{x}{y}}{z} - \frac{1}{2} \cdot \frac{1}{y}\right)} \]
                3. Step-by-step derivation
                  1. Applied rewrites38.0%

                    \[\leadsto \left(\left(0 - \frac{0.5}{y}\right) \cdot z\right) \cdot \color{blue}{z} \]
                  2. Taylor expanded in y around 0

                    \[\leadsto \left(\frac{-1}{2} \cdot \frac{z}{y}\right) \cdot z \]
                  3. Step-by-step derivation
                    1. Applied rewrites38.0%

                      \[\leadsto \left(-0.5 \cdot \frac{z}{y}\right) \cdot z \]

                    if 0.0 < (/.f64 (-.f64 (+.f64 (*.f64 x x) (*.f64 y y)) (*.f64 z z)) (*.f64 y #s(literal 2 binary64))) < 2.00000000000000003e151

                    1. Initial program 98.3%

                      \[\frac{\left(x \cdot x + y \cdot y\right) - z \cdot z}{y \cdot 2} \]
                    2. Add Preprocessing
                    3. Taylor expanded in y around inf

                      \[\leadsto \color{blue}{\frac{1}{2} \cdot y} \]
                    4. Step-by-step derivation
                      1. lower-*.f6464.9

                        \[\leadsto \color{blue}{0.5 \cdot y} \]
                    5. Applied rewrites64.9%

                      \[\leadsto \color{blue}{0.5 \cdot y} \]

                    if 2.00000000000000003e151 < (/.f64 (-.f64 (+.f64 (*.f64 x x) (*.f64 y y)) (*.f64 z z)) (*.f64 y #s(literal 2 binary64))) < +inf.0

                    1. Initial program 71.6%

                      \[\frac{\left(x \cdot x + y \cdot y\right) - z \cdot z}{y \cdot 2} \]
                    2. Add Preprocessing
                    3. Taylor expanded in x around inf

                      \[\leadsto \color{blue}{\frac{1}{2} \cdot \frac{{x}^{2}}{y}} \]
                    4. Step-by-step derivation
                      1. *-commutativeN/A

                        \[\leadsto \color{blue}{\frac{{x}^{2}}{y} \cdot \frac{1}{2}} \]
                      2. lower-*.f64N/A

                        \[\leadsto \color{blue}{\frac{{x}^{2}}{y} \cdot \frac{1}{2}} \]
                      3. lower-/.f64N/A

                        \[\leadsto \color{blue}{\frac{{x}^{2}}{y}} \cdot \frac{1}{2} \]
                      4. unpow2N/A

                        \[\leadsto \frac{\color{blue}{x \cdot x}}{y} \cdot \frac{1}{2} \]
                      5. lower-*.f6440.2

                        \[\leadsto \frac{\color{blue}{x \cdot x}}{y} \cdot 0.5 \]
                    5. Applied rewrites40.2%

                      \[\leadsto \color{blue}{\frac{x \cdot x}{y} \cdot 0.5} \]
                    6. Step-by-step derivation
                      1. Applied rewrites42.3%

                        \[\leadsto \left(\frac{x}{y} \cdot x\right) \cdot 0.5 \]
                    7. Recombined 3 regimes into one program.
                    8. Add Preprocessing

                    Alternative 3: 70.2% accurate, 0.3× speedup?

                    \[\begin{array}{l} x_m = \left|x\right| \\ y\_m = \left|y\right| \\ y\_s = \mathsf{copysign}\left(1, y\right) \\ \begin{array}{l} t_0 := \left(-0.5 \cdot \frac{z}{y\_m}\right) \cdot z\\ t_1 := \frac{\left(x\_m \cdot x\_m + y\_m \cdot y\_m\right) - z \cdot z}{y\_m \cdot 2}\\ y\_s \cdot \begin{array}{l} \mathbf{if}\;t\_1 \leq 0:\\ \;\;\;\;t\_0\\ \mathbf{elif}\;t\_1 \leq 2 \cdot 10^{+151}:\\ \;\;\;\;0.5 \cdot y\_m\\ \mathbf{elif}\;t\_1 \leq \infty:\\ \;\;\;\;\frac{x\_m \cdot x\_m}{y\_m + y\_m}\\ \mathbf{else}:\\ \;\;\;\;t\_0\\ \end{array} \end{array} \end{array} \]
                    x_m = (fabs.f64 x)
                    y\_m = (fabs.f64 y)
                    y\_s = (copysign.f64 #s(literal 1 binary64) y)
                    (FPCore (y_s x_m y_m z)
                     :precision binary64
                     (let* ((t_0 (* (* -0.5 (/ z y_m)) z))
                            (t_1 (/ (- (+ (* x_m x_m) (* y_m y_m)) (* z z)) (* y_m 2.0))))
                       (*
                        y_s
                        (if (<= t_1 0.0)
                          t_0
                          (if (<= t_1 2e+151)
                            (* 0.5 y_m)
                            (if (<= t_1 INFINITY) (/ (* x_m x_m) (+ y_m y_m)) t_0))))))
                    x_m = fabs(x);
                    y\_m = fabs(y);
                    y\_s = copysign(1.0, y);
                    double code(double y_s, double x_m, double y_m, double z) {
                    	double t_0 = (-0.5 * (z / y_m)) * z;
                    	double t_1 = (((x_m * x_m) + (y_m * y_m)) - (z * z)) / (y_m * 2.0);
                    	double tmp;
                    	if (t_1 <= 0.0) {
                    		tmp = t_0;
                    	} else if (t_1 <= 2e+151) {
                    		tmp = 0.5 * y_m;
                    	} else if (t_1 <= ((double) INFINITY)) {
                    		tmp = (x_m * x_m) / (y_m + y_m);
                    	} else {
                    		tmp = t_0;
                    	}
                    	return y_s * tmp;
                    }
                    
                    x_m = Math.abs(x);
                    y\_m = Math.abs(y);
                    y\_s = Math.copySign(1.0, y);
                    public static double code(double y_s, double x_m, double y_m, double z) {
                    	double t_0 = (-0.5 * (z / y_m)) * z;
                    	double t_1 = (((x_m * x_m) + (y_m * y_m)) - (z * z)) / (y_m * 2.0);
                    	double tmp;
                    	if (t_1 <= 0.0) {
                    		tmp = t_0;
                    	} else if (t_1 <= 2e+151) {
                    		tmp = 0.5 * y_m;
                    	} else if (t_1 <= Double.POSITIVE_INFINITY) {
                    		tmp = (x_m * x_m) / (y_m + y_m);
                    	} else {
                    		tmp = t_0;
                    	}
                    	return y_s * tmp;
                    }
                    
                    x_m = math.fabs(x)
                    y\_m = math.fabs(y)
                    y\_s = math.copysign(1.0, y)
                    def code(y_s, x_m, y_m, z):
                    	t_0 = (-0.5 * (z / y_m)) * z
                    	t_1 = (((x_m * x_m) + (y_m * y_m)) - (z * z)) / (y_m * 2.0)
                    	tmp = 0
                    	if t_1 <= 0.0:
                    		tmp = t_0
                    	elif t_1 <= 2e+151:
                    		tmp = 0.5 * y_m
                    	elif t_1 <= math.inf:
                    		tmp = (x_m * x_m) / (y_m + y_m)
                    	else:
                    		tmp = t_0
                    	return y_s * tmp
                    
                    x_m = abs(x)
                    y\_m = abs(y)
                    y\_s = copysign(1.0, y)
                    function code(y_s, x_m, y_m, z)
                    	t_0 = Float64(Float64(-0.5 * Float64(z / y_m)) * z)
                    	t_1 = Float64(Float64(Float64(Float64(x_m * x_m) + Float64(y_m * y_m)) - Float64(z * z)) / Float64(y_m * 2.0))
                    	tmp = 0.0
                    	if (t_1 <= 0.0)
                    		tmp = t_0;
                    	elseif (t_1 <= 2e+151)
                    		tmp = Float64(0.5 * y_m);
                    	elseif (t_1 <= Inf)
                    		tmp = Float64(Float64(x_m * x_m) / Float64(y_m + y_m));
                    	else
                    		tmp = t_0;
                    	end
                    	return Float64(y_s * tmp)
                    end
                    
                    x_m = abs(x);
                    y\_m = abs(y);
                    y\_s = sign(y) * abs(1.0);
                    function tmp_2 = code(y_s, x_m, y_m, z)
                    	t_0 = (-0.5 * (z / y_m)) * z;
                    	t_1 = (((x_m * x_m) + (y_m * y_m)) - (z * z)) / (y_m * 2.0);
                    	tmp = 0.0;
                    	if (t_1 <= 0.0)
                    		tmp = t_0;
                    	elseif (t_1 <= 2e+151)
                    		tmp = 0.5 * y_m;
                    	elseif (t_1 <= Inf)
                    		tmp = (x_m * x_m) / (y_m + y_m);
                    	else
                    		tmp = t_0;
                    	end
                    	tmp_2 = y_s * tmp;
                    end
                    
                    x_m = N[Abs[x], $MachinePrecision]
                    y\_m = N[Abs[y], $MachinePrecision]
                    y\_s = N[With[{TMP1 = Abs[1.0], TMP2 = Sign[y]}, TMP1 * If[TMP2 == 0, 1, TMP2]], $MachinePrecision]
                    code[y$95$s_, x$95$m_, y$95$m_, z_] := Block[{t$95$0 = N[(N[(-0.5 * N[(z / y$95$m), $MachinePrecision]), $MachinePrecision] * z), $MachinePrecision]}, Block[{t$95$1 = N[(N[(N[(N[(x$95$m * x$95$m), $MachinePrecision] + N[(y$95$m * y$95$m), $MachinePrecision]), $MachinePrecision] - N[(z * z), $MachinePrecision]), $MachinePrecision] / N[(y$95$m * 2.0), $MachinePrecision]), $MachinePrecision]}, N[(y$95$s * If[LessEqual[t$95$1, 0.0], t$95$0, If[LessEqual[t$95$1, 2e+151], N[(0.5 * y$95$m), $MachinePrecision], If[LessEqual[t$95$1, Infinity], N[(N[(x$95$m * x$95$m), $MachinePrecision] / N[(y$95$m + y$95$m), $MachinePrecision]), $MachinePrecision], t$95$0]]]), $MachinePrecision]]]
                    
                    \begin{array}{l}
                    x_m = \left|x\right|
                    \\
                    y\_m = \left|y\right|
                    \\
                    y\_s = \mathsf{copysign}\left(1, y\right)
                    
                    \\
                    \begin{array}{l}
                    t_0 := \left(-0.5 \cdot \frac{z}{y\_m}\right) \cdot z\\
                    t_1 := \frac{\left(x\_m \cdot x\_m + y\_m \cdot y\_m\right) - z \cdot z}{y\_m \cdot 2}\\
                    y\_s \cdot \begin{array}{l}
                    \mathbf{if}\;t\_1 \leq 0:\\
                    \;\;\;\;t\_0\\
                    
                    \mathbf{elif}\;t\_1 \leq 2 \cdot 10^{+151}:\\
                    \;\;\;\;0.5 \cdot y\_m\\
                    
                    \mathbf{elif}\;t\_1 \leq \infty:\\
                    \;\;\;\;\frac{x\_m \cdot x\_m}{y\_m + y\_m}\\
                    
                    \mathbf{else}:\\
                    \;\;\;\;t\_0\\
                    
                    
                    \end{array}
                    \end{array}
                    \end{array}
                    
                    Derivation
                    1. Split input into 3 regimes
                    2. if (/.f64 (-.f64 (+.f64 (*.f64 x x) (*.f64 y y)) (*.f64 z z)) (*.f64 y #s(literal 2 binary64))) < 0.0 or +inf.0 < (/.f64 (-.f64 (+.f64 (*.f64 x x) (*.f64 y y)) (*.f64 z z)) (*.f64 y #s(literal 2 binary64)))

                      1. Initial program 68.5%

                        \[\frac{\left(x \cdot x + y \cdot y\right) - z \cdot z}{y \cdot 2} \]
                      2. Add Preprocessing
                      3. Taylor expanded in x around 0

                        \[\leadsto \color{blue}{\frac{1}{2} \cdot \frac{{y}^{2} - {z}^{2}}{y} + \frac{1}{2} \cdot \frac{{x}^{2}}{y}} \]
                      4. Step-by-step derivation
                        1. distribute-lft-outN/A

                          \[\leadsto \color{blue}{\frac{1}{2} \cdot \left(\frac{{y}^{2} - {z}^{2}}{y} + \frac{{x}^{2}}{y}\right)} \]
                        2. *-commutativeN/A

                          \[\leadsto \color{blue}{\left(\frac{{y}^{2} - {z}^{2}}{y} + \frac{{x}^{2}}{y}\right) \cdot \frac{1}{2}} \]
                        3. lower-*.f64N/A

                          \[\leadsto \color{blue}{\left(\frac{{y}^{2} - {z}^{2}}{y} + \frac{{x}^{2}}{y}\right) \cdot \frac{1}{2}} \]
                      5. Applied rewrites99.4%

                        \[\leadsto \color{blue}{\mathsf{fma}\left(z + x, \frac{x - z}{y}, y\right) \cdot 0.5} \]
                      6. Step-by-step derivation
                        1. Applied rewrites99.4%

                          \[\leadsto \left(\frac{x - z}{y} \cdot \left(z + x\right) + y\right) \cdot 0.5 \]
                        2. Taylor expanded in z around inf

                          \[\leadsto {z}^{2} \cdot \color{blue}{\left(\frac{1}{2} \cdot \frac{-1 \cdot \frac{x}{y} + \frac{x}{y}}{z} - \frac{1}{2} \cdot \frac{1}{y}\right)} \]
                        3. Step-by-step derivation
                          1. Applied rewrites38.0%

                            \[\leadsto \left(\left(0 - \frac{0.5}{y}\right) \cdot z\right) \cdot \color{blue}{z} \]
                          2. Taylor expanded in y around 0

                            \[\leadsto \left(\frac{-1}{2} \cdot \frac{z}{y}\right) \cdot z \]
                          3. Step-by-step derivation
                            1. Applied rewrites38.0%

                              \[\leadsto \left(-0.5 \cdot \frac{z}{y}\right) \cdot z \]

                            if 0.0 < (/.f64 (-.f64 (+.f64 (*.f64 x x) (*.f64 y y)) (*.f64 z z)) (*.f64 y #s(literal 2 binary64))) < 2.00000000000000003e151

                            1. Initial program 98.3%

                              \[\frac{\left(x \cdot x + y \cdot y\right) - z \cdot z}{y \cdot 2} \]
                            2. Add Preprocessing
                            3. Taylor expanded in y around inf

                              \[\leadsto \color{blue}{\frac{1}{2} \cdot y} \]
                            4. Step-by-step derivation
                              1. lower-*.f6464.9

                                \[\leadsto \color{blue}{0.5 \cdot y} \]
                            5. Applied rewrites64.9%

                              \[\leadsto \color{blue}{0.5 \cdot y} \]

                            if 2.00000000000000003e151 < (/.f64 (-.f64 (+.f64 (*.f64 x x) (*.f64 y y)) (*.f64 z z)) (*.f64 y #s(literal 2 binary64))) < +inf.0

                            1. Initial program 71.6%

                              \[\frac{\left(x \cdot x + y \cdot y\right) - z \cdot z}{y \cdot 2} \]
                            2. Add Preprocessing
                            3. Taylor expanded in z around inf

                              \[\leadsto \frac{\color{blue}{-1 \cdot {z}^{2}}}{y \cdot 2} \]
                            4. Step-by-step derivation
                              1. mul-1-negN/A

                                \[\leadsto \frac{\color{blue}{\mathsf{neg}\left({z}^{2}\right)}}{y \cdot 2} \]
                              2. unpow2N/A

                                \[\leadsto \frac{\mathsf{neg}\left(\color{blue}{z \cdot z}\right)}{y \cdot 2} \]
                              3. distribute-lft-neg-inN/A

                                \[\leadsto \frac{\color{blue}{\left(\mathsf{neg}\left(z\right)\right) \cdot z}}{y \cdot 2} \]
                              4. lower-*.f64N/A

                                \[\leadsto \frac{\color{blue}{\left(\mathsf{neg}\left(z\right)\right) \cdot z}}{y \cdot 2} \]
                              5. lower-neg.f6431.3

                                \[\leadsto \frac{\color{blue}{\left(-z\right)} \cdot z}{y \cdot 2} \]
                            5. Applied rewrites31.3%

                              \[\leadsto \frac{\color{blue}{\left(-z\right) \cdot z}}{y \cdot 2} \]
                            6. Step-by-step derivation
                              1. lift-*.f64N/A

                                \[\leadsto \frac{\left(-z\right) \cdot z}{\color{blue}{y \cdot 2}} \]
                              2. *-commutativeN/A

                                \[\leadsto \frac{\left(-z\right) \cdot z}{\color{blue}{2 \cdot y}} \]
                              3. count-2-revN/A

                                \[\leadsto \frac{\left(-z\right) \cdot z}{\color{blue}{y + y}} \]
                              4. lower-+.f6431.3

                                \[\leadsto \frac{\left(-z\right) \cdot z}{\color{blue}{y + y}} \]
                            7. Applied rewrites31.3%

                              \[\leadsto \frac{\left(-z\right) \cdot z}{\color{blue}{y + y}} \]
                            8. Taylor expanded in x around inf

                              \[\leadsto \frac{\color{blue}{{x}^{2}}}{y + y} \]
                            9. Step-by-step derivation
                              1. unpow2N/A

                                \[\leadsto \frac{\color{blue}{x \cdot x}}{y + y} \]
                              2. lower-*.f6440.2

                                \[\leadsto \frac{\color{blue}{x \cdot x}}{y + y} \]
                            10. Applied rewrites40.2%

                              \[\leadsto \frac{\color{blue}{x \cdot x}}{y + y} \]
                          4. Recombined 3 regimes into one program.
                          5. Add Preprocessing

                          Alternative 4: 95.9% accurate, 0.3× speedup?

                          \[\begin{array}{l} x_m = \left|x\right| \\ y\_m = \left|y\right| \\ y\_s = \mathsf{copysign}\left(1, y\right) \\ \begin{array}{l} t_0 := \frac{\left(x\_m \cdot x\_m + y\_m \cdot y\_m\right) - z \cdot z}{y\_m \cdot 2}\\ y\_s \cdot \begin{array}{l} \mathbf{if}\;t\_0 \leq 0:\\ \;\;\;\;\left(-0.5 \cdot \frac{z}{y\_m}\right) \cdot z\\ \mathbf{elif}\;t\_0 \leq \infty:\\ \;\;\;\;\mathsf{fma}\left(\frac{x\_m}{y\_m}, x\_m, y\_m\right) \cdot 0.5\\ \mathbf{else}:\\ \;\;\;\;\mathsf{fma}\left(-z, \frac{z}{y\_m}, y\_m\right) \cdot 0.5\\ \end{array} \end{array} \end{array} \]
                          x_m = (fabs.f64 x)
                          y\_m = (fabs.f64 y)
                          y\_s = (copysign.f64 #s(literal 1 binary64) y)
                          (FPCore (y_s x_m y_m z)
                           :precision binary64
                           (let* ((t_0 (/ (- (+ (* x_m x_m) (* y_m y_m)) (* z z)) (* y_m 2.0))))
                             (*
                              y_s
                              (if (<= t_0 0.0)
                                (* (* -0.5 (/ z y_m)) z)
                                (if (<= t_0 INFINITY)
                                  (* (fma (/ x_m y_m) x_m y_m) 0.5)
                                  (* (fma (- z) (/ z y_m) y_m) 0.5))))))
                          x_m = fabs(x);
                          y\_m = fabs(y);
                          y\_s = copysign(1.0, y);
                          double code(double y_s, double x_m, double y_m, double z) {
                          	double t_0 = (((x_m * x_m) + (y_m * y_m)) - (z * z)) / (y_m * 2.0);
                          	double tmp;
                          	if (t_0 <= 0.0) {
                          		tmp = (-0.5 * (z / y_m)) * z;
                          	} else if (t_0 <= ((double) INFINITY)) {
                          		tmp = fma((x_m / y_m), x_m, y_m) * 0.5;
                          	} else {
                          		tmp = fma(-z, (z / y_m), y_m) * 0.5;
                          	}
                          	return y_s * tmp;
                          }
                          
                          x_m = abs(x)
                          y\_m = abs(y)
                          y\_s = copysign(1.0, y)
                          function code(y_s, x_m, y_m, z)
                          	t_0 = Float64(Float64(Float64(Float64(x_m * x_m) + Float64(y_m * y_m)) - Float64(z * z)) / Float64(y_m * 2.0))
                          	tmp = 0.0
                          	if (t_0 <= 0.0)
                          		tmp = Float64(Float64(-0.5 * Float64(z / y_m)) * z);
                          	elseif (t_0 <= Inf)
                          		tmp = Float64(fma(Float64(x_m / y_m), x_m, y_m) * 0.5);
                          	else
                          		tmp = Float64(fma(Float64(-z), Float64(z / y_m), y_m) * 0.5);
                          	end
                          	return Float64(y_s * tmp)
                          end
                          
                          x_m = N[Abs[x], $MachinePrecision]
                          y\_m = N[Abs[y], $MachinePrecision]
                          y\_s = N[With[{TMP1 = Abs[1.0], TMP2 = Sign[y]}, TMP1 * If[TMP2 == 0, 1, TMP2]], $MachinePrecision]
                          code[y$95$s_, x$95$m_, y$95$m_, z_] := Block[{t$95$0 = N[(N[(N[(N[(x$95$m * x$95$m), $MachinePrecision] + N[(y$95$m * y$95$m), $MachinePrecision]), $MachinePrecision] - N[(z * z), $MachinePrecision]), $MachinePrecision] / N[(y$95$m * 2.0), $MachinePrecision]), $MachinePrecision]}, N[(y$95$s * If[LessEqual[t$95$0, 0.0], N[(N[(-0.5 * N[(z / y$95$m), $MachinePrecision]), $MachinePrecision] * z), $MachinePrecision], If[LessEqual[t$95$0, Infinity], N[(N[(N[(x$95$m / y$95$m), $MachinePrecision] * x$95$m + y$95$m), $MachinePrecision] * 0.5), $MachinePrecision], N[(N[((-z) * N[(z / y$95$m), $MachinePrecision] + y$95$m), $MachinePrecision] * 0.5), $MachinePrecision]]]), $MachinePrecision]]
                          
                          \begin{array}{l}
                          x_m = \left|x\right|
                          \\
                          y\_m = \left|y\right|
                          \\
                          y\_s = \mathsf{copysign}\left(1, y\right)
                          
                          \\
                          \begin{array}{l}
                          t_0 := \frac{\left(x\_m \cdot x\_m + y\_m \cdot y\_m\right) - z \cdot z}{y\_m \cdot 2}\\
                          y\_s \cdot \begin{array}{l}
                          \mathbf{if}\;t\_0 \leq 0:\\
                          \;\;\;\;\left(-0.5 \cdot \frac{z}{y\_m}\right) \cdot z\\
                          
                          \mathbf{elif}\;t\_0 \leq \infty:\\
                          \;\;\;\;\mathsf{fma}\left(\frac{x\_m}{y\_m}, x\_m, y\_m\right) \cdot 0.5\\
                          
                          \mathbf{else}:\\
                          \;\;\;\;\mathsf{fma}\left(-z, \frac{z}{y\_m}, y\_m\right) \cdot 0.5\\
                          
                          
                          \end{array}
                          \end{array}
                          \end{array}
                          
                          Derivation
                          1. Split input into 3 regimes
                          2. if (/.f64 (-.f64 (+.f64 (*.f64 x x) (*.f64 y y)) (*.f64 z z)) (*.f64 y #s(literal 2 binary64))) < 0.0

                            1. Initial program 83.8%

                              \[\frac{\left(x \cdot x + y \cdot y\right) - z \cdot z}{y \cdot 2} \]
                            2. Add Preprocessing
                            3. Taylor expanded in x around 0

                              \[\leadsto \color{blue}{\frac{1}{2} \cdot \frac{{y}^{2} - {z}^{2}}{y} + \frac{1}{2} \cdot \frac{{x}^{2}}{y}} \]
                            4. Step-by-step derivation
                              1. distribute-lft-outN/A

                                \[\leadsto \color{blue}{\frac{1}{2} \cdot \left(\frac{{y}^{2} - {z}^{2}}{y} + \frac{{x}^{2}}{y}\right)} \]
                              2. *-commutativeN/A

                                \[\leadsto \color{blue}{\left(\frac{{y}^{2} - {z}^{2}}{y} + \frac{{x}^{2}}{y}\right) \cdot \frac{1}{2}} \]
                              3. lower-*.f64N/A

                                \[\leadsto \color{blue}{\left(\frac{{y}^{2} - {z}^{2}}{y} + \frac{{x}^{2}}{y}\right) \cdot \frac{1}{2}} \]
                            5. Applied rewrites99.3%

                              \[\leadsto \color{blue}{\mathsf{fma}\left(z + x, \frac{x - z}{y}, y\right) \cdot 0.5} \]
                            6. Step-by-step derivation
                              1. Applied rewrites99.3%

                                \[\leadsto \left(\frac{x - z}{y} \cdot \left(z + x\right) + y\right) \cdot 0.5 \]
                              2. Taylor expanded in z around inf

                                \[\leadsto {z}^{2} \cdot \color{blue}{\left(\frac{1}{2} \cdot \frac{-1 \cdot \frac{x}{y} + \frac{x}{y}}{z} - \frac{1}{2} \cdot \frac{1}{y}\right)} \]
                              3. Step-by-step derivation
                                1. Applied rewrites32.7%

                                  \[\leadsto \left(\left(0 - \frac{0.5}{y}\right) \cdot z\right) \cdot \color{blue}{z} \]
                                2. Taylor expanded in y around 0

                                  \[\leadsto \left(\frac{-1}{2} \cdot \frac{z}{y}\right) \cdot z \]
                                3. Step-by-step derivation
                                  1. Applied rewrites32.7%

                                    \[\leadsto \left(-0.5 \cdot \frac{z}{y}\right) \cdot z \]

                                  if 0.0 < (/.f64 (-.f64 (+.f64 (*.f64 x x) (*.f64 y y)) (*.f64 z z)) (*.f64 y #s(literal 2 binary64))) < +inf.0

                                  1. Initial program 76.9%

                                    \[\frac{\left(x \cdot x + y \cdot y\right) - z \cdot z}{y \cdot 2} \]
                                  2. Add Preprocessing
                                  3. Taylor expanded in z around 0

                                    \[\leadsto \color{blue}{\frac{1}{2} \cdot \frac{{x}^{2} + {y}^{2}}{y}} \]
                                  4. Step-by-step derivation
                                    1. div-addN/A

                                      \[\leadsto \frac{1}{2} \cdot \color{blue}{\left(\frac{{x}^{2}}{y} + \frac{{y}^{2}}{y}\right)} \]
                                    2. distribute-lft-inN/A

                                      \[\leadsto \color{blue}{\frac{1}{2} \cdot \frac{{x}^{2}}{y} + \frac{1}{2} \cdot \frac{{y}^{2}}{y}} \]
                                    3. unpow2N/A

                                      \[\leadsto \frac{1}{2} \cdot \frac{{x}^{2}}{y} + \frac{1}{2} \cdot \frac{\color{blue}{y \cdot y}}{y} \]
                                    4. associate-/l*N/A

                                      \[\leadsto \frac{1}{2} \cdot \frac{{x}^{2}}{y} + \frac{1}{2} \cdot \color{blue}{\left(y \cdot \frac{y}{y}\right)} \]
                                    5. *-inversesN/A

                                      \[\leadsto \frac{1}{2} \cdot \frac{{x}^{2}}{y} + \frac{1}{2} \cdot \left(y \cdot \color{blue}{1}\right) \]
                                    6. *-rgt-identityN/A

                                      \[\leadsto \frac{1}{2} \cdot \frac{{x}^{2}}{y} + \frac{1}{2} \cdot \color{blue}{y} \]
                                    7. distribute-lft-inN/A

                                      \[\leadsto \color{blue}{\frac{1}{2} \cdot \left(\frac{{x}^{2}}{y} + y\right)} \]
                                    8. +-commutativeN/A

                                      \[\leadsto \frac{1}{2} \cdot \color{blue}{\left(y + \frac{{x}^{2}}{y}\right)} \]
                                    9. *-commutativeN/A

                                      \[\leadsto \color{blue}{\left(y + \frac{{x}^{2}}{y}\right) \cdot \frac{1}{2}} \]
                                    10. lower-*.f64N/A

                                      \[\leadsto \color{blue}{\left(y + \frac{{x}^{2}}{y}\right) \cdot \frac{1}{2}} \]
                                    11. +-commutativeN/A

                                      \[\leadsto \color{blue}{\left(\frac{{x}^{2}}{y} + y\right)} \cdot \frac{1}{2} \]
                                    12. unpow2N/A

                                      \[\leadsto \left(\frac{\color{blue}{x \cdot x}}{y} + y\right) \cdot \frac{1}{2} \]
                                    13. associate-/l*N/A

                                      \[\leadsto \left(\color{blue}{x \cdot \frac{x}{y}} + y\right) \cdot \frac{1}{2} \]
                                    14. *-commutativeN/A

                                      \[\leadsto \left(\color{blue}{\frac{x}{y} \cdot x} + y\right) \cdot \frac{1}{2} \]
                                    15. lower-fma.f64N/A

                                      \[\leadsto \color{blue}{\mathsf{fma}\left(\frac{x}{y}, x, y\right)} \cdot \frac{1}{2} \]
                                    16. lower-/.f6470.1

                                      \[\leadsto \mathsf{fma}\left(\color{blue}{\frac{x}{y}}, x, y\right) \cdot 0.5 \]
                                  5. Applied rewrites70.1%

                                    \[\leadsto \color{blue}{\mathsf{fma}\left(\frac{x}{y}, x, y\right) \cdot 0.5} \]

                                  if +inf.0 < (/.f64 (-.f64 (+.f64 (*.f64 x x) (*.f64 y y)) (*.f64 z z)) (*.f64 y #s(literal 2 binary64)))

                                  1. Initial program 0.0%

                                    \[\frac{\left(x \cdot x + y \cdot y\right) - z \cdot z}{y \cdot 2} \]
                                  2. Add Preprocessing
                                  3. Taylor expanded in x around 0

                                    \[\leadsto \color{blue}{\frac{1}{2} \cdot \frac{{y}^{2} - {z}^{2}}{y} + \frac{1}{2} \cdot \frac{{x}^{2}}{y}} \]
                                  4. Step-by-step derivation
                                    1. distribute-lft-outN/A

                                      \[\leadsto \color{blue}{\frac{1}{2} \cdot \left(\frac{{y}^{2} - {z}^{2}}{y} + \frac{{x}^{2}}{y}\right)} \]
                                    2. *-commutativeN/A

                                      \[\leadsto \color{blue}{\left(\frac{{y}^{2} - {z}^{2}}{y} + \frac{{x}^{2}}{y}\right) \cdot \frac{1}{2}} \]
                                    3. lower-*.f64N/A

                                      \[\leadsto \color{blue}{\left(\frac{{y}^{2} - {z}^{2}}{y} + \frac{{x}^{2}}{y}\right) \cdot \frac{1}{2}} \]
                                  5. Applied rewrites99.9%

                                    \[\leadsto \color{blue}{\mathsf{fma}\left(z + x, \frac{x - z}{y}, y\right) \cdot 0.5} \]
                                  6. Taylor expanded in x around 0

                                    \[\leadsto \mathsf{fma}\left(z + x, -1 \cdot \frac{z}{y}, y\right) \cdot \frac{1}{2} \]
                                  7. Step-by-step derivation
                                    1. Applied rewrites92.7%

                                      \[\leadsto \mathsf{fma}\left(z + x, \frac{-z}{y}, y\right) \cdot 0.5 \]
                                    2. Step-by-step derivation
                                      1. Applied rewrites96.3%

                                        \[\leadsto \mathsf{fma}\left(x - z, \frac{z}{y}, y\right) \cdot 0.5 \]
                                      2. Taylor expanded in x around 0

                                        \[\leadsto \mathsf{fma}\left(-1 \cdot z, \frac{z}{y}, y\right) \cdot \frac{1}{2} \]
                                      3. Step-by-step derivation
                                        1. Applied rewrites89.2%

                                          \[\leadsto \mathsf{fma}\left(-z, \frac{z}{y}, y\right) \cdot 0.5 \]
                                      4. Recombined 3 regimes into one program.
                                      5. Add Preprocessing

                                      Alternative 5: 93.0% accurate, 0.3× speedup?

                                      \[\begin{array}{l} x_m = \left|x\right| \\ y\_m = \left|y\right| \\ y\_s = \mathsf{copysign}\left(1, y\right) \\ \begin{array}{l} t_0 := \frac{\left(x\_m \cdot x\_m + y\_m \cdot y\_m\right) - z \cdot z}{y\_m \cdot 2}\\ y\_s \cdot \begin{array}{l} \mathbf{if}\;t\_0 \leq 0 \lor \neg \left(t\_0 \leq \infty\right):\\ \;\;\;\;\left(-0.5 \cdot \frac{z}{y\_m}\right) \cdot z\\ \mathbf{else}:\\ \;\;\;\;\mathsf{fma}\left(\frac{x\_m}{y\_m}, x\_m, y\_m\right) \cdot 0.5\\ \end{array} \end{array} \end{array} \]
                                      x_m = (fabs.f64 x)
                                      y\_m = (fabs.f64 y)
                                      y\_s = (copysign.f64 #s(literal 1 binary64) y)
                                      (FPCore (y_s x_m y_m z)
                                       :precision binary64
                                       (let* ((t_0 (/ (- (+ (* x_m x_m) (* y_m y_m)) (* z z)) (* y_m 2.0))))
                                         (*
                                          y_s
                                          (if (or (<= t_0 0.0) (not (<= t_0 INFINITY)))
                                            (* (* -0.5 (/ z y_m)) z)
                                            (* (fma (/ x_m y_m) x_m y_m) 0.5)))))
                                      x_m = fabs(x);
                                      y\_m = fabs(y);
                                      y\_s = copysign(1.0, y);
                                      double code(double y_s, double x_m, double y_m, double z) {
                                      	double t_0 = (((x_m * x_m) + (y_m * y_m)) - (z * z)) / (y_m * 2.0);
                                      	double tmp;
                                      	if ((t_0 <= 0.0) || !(t_0 <= ((double) INFINITY))) {
                                      		tmp = (-0.5 * (z / y_m)) * z;
                                      	} else {
                                      		tmp = fma((x_m / y_m), x_m, y_m) * 0.5;
                                      	}
                                      	return y_s * tmp;
                                      }
                                      
                                      x_m = abs(x)
                                      y\_m = abs(y)
                                      y\_s = copysign(1.0, y)
                                      function code(y_s, x_m, y_m, z)
                                      	t_0 = Float64(Float64(Float64(Float64(x_m * x_m) + Float64(y_m * y_m)) - Float64(z * z)) / Float64(y_m * 2.0))
                                      	tmp = 0.0
                                      	if ((t_0 <= 0.0) || !(t_0 <= Inf))
                                      		tmp = Float64(Float64(-0.5 * Float64(z / y_m)) * z);
                                      	else
                                      		tmp = Float64(fma(Float64(x_m / y_m), x_m, y_m) * 0.5);
                                      	end
                                      	return Float64(y_s * tmp)
                                      end
                                      
                                      x_m = N[Abs[x], $MachinePrecision]
                                      y\_m = N[Abs[y], $MachinePrecision]
                                      y\_s = N[With[{TMP1 = Abs[1.0], TMP2 = Sign[y]}, TMP1 * If[TMP2 == 0, 1, TMP2]], $MachinePrecision]
                                      code[y$95$s_, x$95$m_, y$95$m_, z_] := Block[{t$95$0 = N[(N[(N[(N[(x$95$m * x$95$m), $MachinePrecision] + N[(y$95$m * y$95$m), $MachinePrecision]), $MachinePrecision] - N[(z * z), $MachinePrecision]), $MachinePrecision] / N[(y$95$m * 2.0), $MachinePrecision]), $MachinePrecision]}, N[(y$95$s * If[Or[LessEqual[t$95$0, 0.0], N[Not[LessEqual[t$95$0, Infinity]], $MachinePrecision]], N[(N[(-0.5 * N[(z / y$95$m), $MachinePrecision]), $MachinePrecision] * z), $MachinePrecision], N[(N[(N[(x$95$m / y$95$m), $MachinePrecision] * x$95$m + y$95$m), $MachinePrecision] * 0.5), $MachinePrecision]]), $MachinePrecision]]
                                      
                                      \begin{array}{l}
                                      x_m = \left|x\right|
                                      \\
                                      y\_m = \left|y\right|
                                      \\
                                      y\_s = \mathsf{copysign}\left(1, y\right)
                                      
                                      \\
                                      \begin{array}{l}
                                      t_0 := \frac{\left(x\_m \cdot x\_m + y\_m \cdot y\_m\right) - z \cdot z}{y\_m \cdot 2}\\
                                      y\_s \cdot \begin{array}{l}
                                      \mathbf{if}\;t\_0 \leq 0 \lor \neg \left(t\_0 \leq \infty\right):\\
                                      \;\;\;\;\left(-0.5 \cdot \frac{z}{y\_m}\right) \cdot z\\
                                      
                                      \mathbf{else}:\\
                                      \;\;\;\;\mathsf{fma}\left(\frac{x\_m}{y\_m}, x\_m, y\_m\right) \cdot 0.5\\
                                      
                                      
                                      \end{array}
                                      \end{array}
                                      \end{array}
                                      
                                      Derivation
                                      1. Split input into 2 regimes
                                      2. if (/.f64 (-.f64 (+.f64 (*.f64 x x) (*.f64 y y)) (*.f64 z z)) (*.f64 y #s(literal 2 binary64))) < 0.0 or +inf.0 < (/.f64 (-.f64 (+.f64 (*.f64 x x) (*.f64 y y)) (*.f64 z z)) (*.f64 y #s(literal 2 binary64)))

                                        1. Initial program 68.5%

                                          \[\frac{\left(x \cdot x + y \cdot y\right) - z \cdot z}{y \cdot 2} \]
                                        2. Add Preprocessing
                                        3. Taylor expanded in x around 0

                                          \[\leadsto \color{blue}{\frac{1}{2} \cdot \frac{{y}^{2} - {z}^{2}}{y} + \frac{1}{2} \cdot \frac{{x}^{2}}{y}} \]
                                        4. Step-by-step derivation
                                          1. distribute-lft-outN/A

                                            \[\leadsto \color{blue}{\frac{1}{2} \cdot \left(\frac{{y}^{2} - {z}^{2}}{y} + \frac{{x}^{2}}{y}\right)} \]
                                          2. *-commutativeN/A

                                            \[\leadsto \color{blue}{\left(\frac{{y}^{2} - {z}^{2}}{y} + \frac{{x}^{2}}{y}\right) \cdot \frac{1}{2}} \]
                                          3. lower-*.f64N/A

                                            \[\leadsto \color{blue}{\left(\frac{{y}^{2} - {z}^{2}}{y} + \frac{{x}^{2}}{y}\right) \cdot \frac{1}{2}} \]
                                        5. Applied rewrites99.4%

                                          \[\leadsto \color{blue}{\mathsf{fma}\left(z + x, \frac{x - z}{y}, y\right) \cdot 0.5} \]
                                        6. Step-by-step derivation
                                          1. Applied rewrites99.4%

                                            \[\leadsto \left(\frac{x - z}{y} \cdot \left(z + x\right) + y\right) \cdot 0.5 \]
                                          2. Taylor expanded in z around inf

                                            \[\leadsto {z}^{2} \cdot \color{blue}{\left(\frac{1}{2} \cdot \frac{-1 \cdot \frac{x}{y} + \frac{x}{y}}{z} - \frac{1}{2} \cdot \frac{1}{y}\right)} \]
                                          3. Step-by-step derivation
                                            1. Applied rewrites38.0%

                                              \[\leadsto \left(\left(0 - \frac{0.5}{y}\right) \cdot z\right) \cdot \color{blue}{z} \]
                                            2. Taylor expanded in y around 0

                                              \[\leadsto \left(\frac{-1}{2} \cdot \frac{z}{y}\right) \cdot z \]
                                            3. Step-by-step derivation
                                              1. Applied rewrites38.0%

                                                \[\leadsto \left(-0.5 \cdot \frac{z}{y}\right) \cdot z \]

                                              if 0.0 < (/.f64 (-.f64 (+.f64 (*.f64 x x) (*.f64 y y)) (*.f64 z z)) (*.f64 y #s(literal 2 binary64))) < +inf.0

                                              1. Initial program 76.9%

                                                \[\frac{\left(x \cdot x + y \cdot y\right) - z \cdot z}{y \cdot 2} \]
                                              2. Add Preprocessing
                                              3. Taylor expanded in z around 0

                                                \[\leadsto \color{blue}{\frac{1}{2} \cdot \frac{{x}^{2} + {y}^{2}}{y}} \]
                                              4. Step-by-step derivation
                                                1. div-addN/A

                                                  \[\leadsto \frac{1}{2} \cdot \color{blue}{\left(\frac{{x}^{2}}{y} + \frac{{y}^{2}}{y}\right)} \]
                                                2. distribute-lft-inN/A

                                                  \[\leadsto \color{blue}{\frac{1}{2} \cdot \frac{{x}^{2}}{y} + \frac{1}{2} \cdot \frac{{y}^{2}}{y}} \]
                                                3. unpow2N/A

                                                  \[\leadsto \frac{1}{2} \cdot \frac{{x}^{2}}{y} + \frac{1}{2} \cdot \frac{\color{blue}{y \cdot y}}{y} \]
                                                4. associate-/l*N/A

                                                  \[\leadsto \frac{1}{2} \cdot \frac{{x}^{2}}{y} + \frac{1}{2} \cdot \color{blue}{\left(y \cdot \frac{y}{y}\right)} \]
                                                5. *-inversesN/A

                                                  \[\leadsto \frac{1}{2} \cdot \frac{{x}^{2}}{y} + \frac{1}{2} \cdot \left(y \cdot \color{blue}{1}\right) \]
                                                6. *-rgt-identityN/A

                                                  \[\leadsto \frac{1}{2} \cdot \frac{{x}^{2}}{y} + \frac{1}{2} \cdot \color{blue}{y} \]
                                                7. distribute-lft-inN/A

                                                  \[\leadsto \color{blue}{\frac{1}{2} \cdot \left(\frac{{x}^{2}}{y} + y\right)} \]
                                                8. +-commutativeN/A

                                                  \[\leadsto \frac{1}{2} \cdot \color{blue}{\left(y + \frac{{x}^{2}}{y}\right)} \]
                                                9. *-commutativeN/A

                                                  \[\leadsto \color{blue}{\left(y + \frac{{x}^{2}}{y}\right) \cdot \frac{1}{2}} \]
                                                10. lower-*.f64N/A

                                                  \[\leadsto \color{blue}{\left(y + \frac{{x}^{2}}{y}\right) \cdot \frac{1}{2}} \]
                                                11. +-commutativeN/A

                                                  \[\leadsto \color{blue}{\left(\frac{{x}^{2}}{y} + y\right)} \cdot \frac{1}{2} \]
                                                12. unpow2N/A

                                                  \[\leadsto \left(\frac{\color{blue}{x \cdot x}}{y} + y\right) \cdot \frac{1}{2} \]
                                                13. associate-/l*N/A

                                                  \[\leadsto \left(\color{blue}{x \cdot \frac{x}{y}} + y\right) \cdot \frac{1}{2} \]
                                                14. *-commutativeN/A

                                                  \[\leadsto \left(\color{blue}{\frac{x}{y} \cdot x} + y\right) \cdot \frac{1}{2} \]
                                                15. lower-fma.f64N/A

                                                  \[\leadsto \color{blue}{\mathsf{fma}\left(\frac{x}{y}, x, y\right)} \cdot \frac{1}{2} \]
                                                16. lower-/.f6470.1

                                                  \[\leadsto \mathsf{fma}\left(\color{blue}{\frac{x}{y}}, x, y\right) \cdot 0.5 \]
                                              5. Applied rewrites70.1%

                                                \[\leadsto \color{blue}{\mathsf{fma}\left(\frac{x}{y}, x, y\right) \cdot 0.5} \]
                                            4. Recombined 2 regimes into one program.
                                            5. Final simplification50.8%

                                              \[\leadsto \begin{array}{l} \mathbf{if}\;\frac{\left(x \cdot x + y \cdot y\right) - z \cdot z}{y \cdot 2} \leq 0 \lor \neg \left(\frac{\left(x \cdot x + y \cdot y\right) - z \cdot z}{y \cdot 2} \leq \infty\right):\\ \;\;\;\;\left(-0.5 \cdot \frac{z}{y}\right) \cdot z\\ \mathbf{else}:\\ \;\;\;\;\mathsf{fma}\left(\frac{x}{y}, x, y\right) \cdot 0.5\\ \end{array} \]
                                            6. Add Preprocessing

                                            Alternative 6: 51.2% accurate, 1.2× speedup?

                                            \[\begin{array}{l} x_m = \left|x\right| \\ y\_m = \left|y\right| \\ y\_s = \mathsf{copysign}\left(1, y\right) \\ y\_s \cdot \begin{array}{l} \mathbf{if}\;x\_m \cdot x\_m \leq 8.6 \cdot 10^{+176}:\\ \;\;\;\;0.5 \cdot y\_m\\ \mathbf{else}:\\ \;\;\;\;\frac{x\_m \cdot x\_m}{y\_m + y\_m}\\ \end{array} \end{array} \]
                                            x_m = (fabs.f64 x)
                                            y\_m = (fabs.f64 y)
                                            y\_s = (copysign.f64 #s(literal 1 binary64) y)
                                            (FPCore (y_s x_m y_m z)
                                             :precision binary64
                                             (*
                                              y_s
                                              (if (<= (* x_m x_m) 8.6e+176) (* 0.5 y_m) (/ (* x_m x_m) (+ y_m y_m)))))
                                            x_m = fabs(x);
                                            y\_m = fabs(y);
                                            y\_s = copysign(1.0, y);
                                            double code(double y_s, double x_m, double y_m, double z) {
                                            	double tmp;
                                            	if ((x_m * x_m) <= 8.6e+176) {
                                            		tmp = 0.5 * y_m;
                                            	} else {
                                            		tmp = (x_m * x_m) / (y_m + y_m);
                                            	}
                                            	return y_s * tmp;
                                            }
                                            
                                            x_m = abs(x)
                                            y\_m = abs(y)
                                            y\_s = copysign(1.0d0, y)
                                            real(8) function code(y_s, x_m, y_m, z)
                                                real(8), intent (in) :: y_s
                                                real(8), intent (in) :: x_m
                                                real(8), intent (in) :: y_m
                                                real(8), intent (in) :: z
                                                real(8) :: tmp
                                                if ((x_m * x_m) <= 8.6d+176) then
                                                    tmp = 0.5d0 * y_m
                                                else
                                                    tmp = (x_m * x_m) / (y_m + y_m)
                                                end if
                                                code = y_s * tmp
                                            end function
                                            
                                            x_m = Math.abs(x);
                                            y\_m = Math.abs(y);
                                            y\_s = Math.copySign(1.0, y);
                                            public static double code(double y_s, double x_m, double y_m, double z) {
                                            	double tmp;
                                            	if ((x_m * x_m) <= 8.6e+176) {
                                            		tmp = 0.5 * y_m;
                                            	} else {
                                            		tmp = (x_m * x_m) / (y_m + y_m);
                                            	}
                                            	return y_s * tmp;
                                            }
                                            
                                            x_m = math.fabs(x)
                                            y\_m = math.fabs(y)
                                            y\_s = math.copysign(1.0, y)
                                            def code(y_s, x_m, y_m, z):
                                            	tmp = 0
                                            	if (x_m * x_m) <= 8.6e+176:
                                            		tmp = 0.5 * y_m
                                            	else:
                                            		tmp = (x_m * x_m) / (y_m + y_m)
                                            	return y_s * tmp
                                            
                                            x_m = abs(x)
                                            y\_m = abs(y)
                                            y\_s = copysign(1.0, y)
                                            function code(y_s, x_m, y_m, z)
                                            	tmp = 0.0
                                            	if (Float64(x_m * x_m) <= 8.6e+176)
                                            		tmp = Float64(0.5 * y_m);
                                            	else
                                            		tmp = Float64(Float64(x_m * x_m) / Float64(y_m + y_m));
                                            	end
                                            	return Float64(y_s * tmp)
                                            end
                                            
                                            x_m = abs(x);
                                            y\_m = abs(y);
                                            y\_s = sign(y) * abs(1.0);
                                            function tmp_2 = code(y_s, x_m, y_m, z)
                                            	tmp = 0.0;
                                            	if ((x_m * x_m) <= 8.6e+176)
                                            		tmp = 0.5 * y_m;
                                            	else
                                            		tmp = (x_m * x_m) / (y_m + y_m);
                                            	end
                                            	tmp_2 = y_s * tmp;
                                            end
                                            
                                            x_m = N[Abs[x], $MachinePrecision]
                                            y\_m = N[Abs[y], $MachinePrecision]
                                            y\_s = N[With[{TMP1 = Abs[1.0], TMP2 = Sign[y]}, TMP1 * If[TMP2 == 0, 1, TMP2]], $MachinePrecision]
                                            code[y$95$s_, x$95$m_, y$95$m_, z_] := N[(y$95$s * If[LessEqual[N[(x$95$m * x$95$m), $MachinePrecision], 8.6e+176], N[(0.5 * y$95$m), $MachinePrecision], N[(N[(x$95$m * x$95$m), $MachinePrecision] / N[(y$95$m + y$95$m), $MachinePrecision]), $MachinePrecision]]), $MachinePrecision]
                                            
                                            \begin{array}{l}
                                            x_m = \left|x\right|
                                            \\
                                            y\_m = \left|y\right|
                                            \\
                                            y\_s = \mathsf{copysign}\left(1, y\right)
                                            
                                            \\
                                            y\_s \cdot \begin{array}{l}
                                            \mathbf{if}\;x\_m \cdot x\_m \leq 8.6 \cdot 10^{+176}:\\
                                            \;\;\;\;0.5 \cdot y\_m\\
                                            
                                            \mathbf{else}:\\
                                            \;\;\;\;\frac{x\_m \cdot x\_m}{y\_m + y\_m}\\
                                            
                                            
                                            \end{array}
                                            \end{array}
                                            
                                            Derivation
                                            1. Split input into 2 regimes
                                            2. if (*.f64 x x) < 8.60000000000000051e176

                                              1. Initial program 71.1%

                                                \[\frac{\left(x \cdot x + y \cdot y\right) - z \cdot z}{y \cdot 2} \]
                                              2. Add Preprocessing
                                              3. Taylor expanded in y around inf

                                                \[\leadsto \color{blue}{\frac{1}{2} \cdot y} \]
                                              4. Step-by-step derivation
                                                1. lower-*.f6447.3

                                                  \[\leadsto \color{blue}{0.5 \cdot y} \]
                                              5. Applied rewrites47.3%

                                                \[\leadsto \color{blue}{0.5 \cdot y} \]

                                              if 8.60000000000000051e176 < (*.f64 x x)

                                              1. Initial program 73.0%

                                                \[\frac{\left(x \cdot x + y \cdot y\right) - z \cdot z}{y \cdot 2} \]
                                              2. Add Preprocessing
                                              3. Taylor expanded in z around inf

                                                \[\leadsto \frac{\color{blue}{-1 \cdot {z}^{2}}}{y \cdot 2} \]
                                              4. Step-by-step derivation
                                                1. mul-1-negN/A

                                                  \[\leadsto \frac{\color{blue}{\mathsf{neg}\left({z}^{2}\right)}}{y \cdot 2} \]
                                                2. unpow2N/A

                                                  \[\leadsto \frac{\mathsf{neg}\left(\color{blue}{z \cdot z}\right)}{y \cdot 2} \]
                                                3. distribute-lft-neg-inN/A

                                                  \[\leadsto \frac{\color{blue}{\left(\mathsf{neg}\left(z\right)\right) \cdot z}}{y \cdot 2} \]
                                                4. lower-*.f64N/A

                                                  \[\leadsto \frac{\color{blue}{\left(\mathsf{neg}\left(z\right)\right) \cdot z}}{y \cdot 2} \]
                                                5. lower-neg.f6414.9

                                                  \[\leadsto \frac{\color{blue}{\left(-z\right)} \cdot z}{y \cdot 2} \]
                                              5. Applied rewrites14.9%

                                                \[\leadsto \frac{\color{blue}{\left(-z\right) \cdot z}}{y \cdot 2} \]
                                              6. Step-by-step derivation
                                                1. lift-*.f64N/A

                                                  \[\leadsto \frac{\left(-z\right) \cdot z}{\color{blue}{y \cdot 2}} \]
                                                2. *-commutativeN/A

                                                  \[\leadsto \frac{\left(-z\right) \cdot z}{\color{blue}{2 \cdot y}} \]
                                                3. count-2-revN/A

                                                  \[\leadsto \frac{\left(-z\right) \cdot z}{\color{blue}{y + y}} \]
                                                4. lower-+.f6414.9

                                                  \[\leadsto \frac{\left(-z\right) \cdot z}{\color{blue}{y + y}} \]
                                              7. Applied rewrites14.9%

                                                \[\leadsto \frac{\left(-z\right) \cdot z}{\color{blue}{y + y}} \]
                                              8. Taylor expanded in x around inf

                                                \[\leadsto \frac{\color{blue}{{x}^{2}}}{y + y} \]
                                              9. Step-by-step derivation
                                                1. unpow2N/A

                                                  \[\leadsto \frac{\color{blue}{x \cdot x}}{y + y} \]
                                                2. lower-*.f6469.2

                                                  \[\leadsto \frac{\color{blue}{x \cdot x}}{y + y} \]
                                              10. Applied rewrites69.2%

                                                \[\leadsto \frac{\color{blue}{x \cdot x}}{y + y} \]
                                            3. Recombined 2 regimes into one program.
                                            4. Add Preprocessing

                                            Alternative 7: 99.9% accurate, 1.2× speedup?

                                            \[\begin{array}{l} x_m = \left|x\right| \\ y\_m = \left|y\right| \\ y\_s = \mathsf{copysign}\left(1, y\right) \\ y\_s \cdot \left(\left(\frac{x\_m - z}{y\_m} \cdot \left(z + x\_m\right) + y\_m\right) \cdot 0.5\right) \end{array} \]
                                            x_m = (fabs.f64 x)
                                            y\_m = (fabs.f64 y)
                                            y\_s = (copysign.f64 #s(literal 1 binary64) y)
                                            (FPCore (y_s x_m y_m z)
                                             :precision binary64
                                             (* y_s (* (+ (* (/ (- x_m z) y_m) (+ z x_m)) y_m) 0.5)))
                                            x_m = fabs(x);
                                            y\_m = fabs(y);
                                            y\_s = copysign(1.0, y);
                                            double code(double y_s, double x_m, double y_m, double z) {
                                            	return y_s * (((((x_m - z) / y_m) * (z + x_m)) + y_m) * 0.5);
                                            }
                                            
                                            x_m = abs(x)
                                            y\_m = abs(y)
                                            y\_s = copysign(1.0d0, y)
                                            real(8) function code(y_s, x_m, y_m, z)
                                                real(8), intent (in) :: y_s
                                                real(8), intent (in) :: x_m
                                                real(8), intent (in) :: y_m
                                                real(8), intent (in) :: z
                                                code = y_s * (((((x_m - z) / y_m) * (z + x_m)) + y_m) * 0.5d0)
                                            end function
                                            
                                            x_m = Math.abs(x);
                                            y\_m = Math.abs(y);
                                            y\_s = Math.copySign(1.0, y);
                                            public static double code(double y_s, double x_m, double y_m, double z) {
                                            	return y_s * (((((x_m - z) / y_m) * (z + x_m)) + y_m) * 0.5);
                                            }
                                            
                                            x_m = math.fabs(x)
                                            y\_m = math.fabs(y)
                                            y\_s = math.copysign(1.0, y)
                                            def code(y_s, x_m, y_m, z):
                                            	return y_s * (((((x_m - z) / y_m) * (z + x_m)) + y_m) * 0.5)
                                            
                                            x_m = abs(x)
                                            y\_m = abs(y)
                                            y\_s = copysign(1.0, y)
                                            function code(y_s, x_m, y_m, z)
                                            	return Float64(y_s * Float64(Float64(Float64(Float64(Float64(x_m - z) / y_m) * Float64(z + x_m)) + y_m) * 0.5))
                                            end
                                            
                                            x_m = abs(x);
                                            y\_m = abs(y);
                                            y\_s = sign(y) * abs(1.0);
                                            function tmp = code(y_s, x_m, y_m, z)
                                            	tmp = y_s * (((((x_m - z) / y_m) * (z + x_m)) + y_m) * 0.5);
                                            end
                                            
                                            x_m = N[Abs[x], $MachinePrecision]
                                            y\_m = N[Abs[y], $MachinePrecision]
                                            y\_s = N[With[{TMP1 = Abs[1.0], TMP2 = Sign[y]}, TMP1 * If[TMP2 == 0, 1, TMP2]], $MachinePrecision]
                                            code[y$95$s_, x$95$m_, y$95$m_, z_] := N[(y$95$s * N[(N[(N[(N[(N[(x$95$m - z), $MachinePrecision] / y$95$m), $MachinePrecision] * N[(z + x$95$m), $MachinePrecision]), $MachinePrecision] + y$95$m), $MachinePrecision] * 0.5), $MachinePrecision]), $MachinePrecision]
                                            
                                            \begin{array}{l}
                                            x_m = \left|x\right|
                                            \\
                                            y\_m = \left|y\right|
                                            \\
                                            y\_s = \mathsf{copysign}\left(1, y\right)
                                            
                                            \\
                                            y\_s \cdot \left(\left(\frac{x\_m - z}{y\_m} \cdot \left(z + x\_m\right) + y\_m\right) \cdot 0.5\right)
                                            \end{array}
                                            
                                            Derivation
                                            1. Initial program 71.8%

                                              \[\frac{\left(x \cdot x + y \cdot y\right) - z \cdot z}{y \cdot 2} \]
                                            2. Add Preprocessing
                                            3. Taylor expanded in x around 0

                                              \[\leadsto \color{blue}{\frac{1}{2} \cdot \frac{{y}^{2} - {z}^{2}}{y} + \frac{1}{2} \cdot \frac{{x}^{2}}{y}} \]
                                            4. Step-by-step derivation
                                              1. distribute-lft-outN/A

                                                \[\leadsto \color{blue}{\frac{1}{2} \cdot \left(\frac{{y}^{2} - {z}^{2}}{y} + \frac{{x}^{2}}{y}\right)} \]
                                              2. *-commutativeN/A

                                                \[\leadsto \color{blue}{\left(\frac{{y}^{2} - {z}^{2}}{y} + \frac{{x}^{2}}{y}\right) \cdot \frac{1}{2}} \]
                                              3. lower-*.f64N/A

                                                \[\leadsto \color{blue}{\left(\frac{{y}^{2} - {z}^{2}}{y} + \frac{{x}^{2}}{y}\right) \cdot \frac{1}{2}} \]
                                            5. Applied rewrites99.6%

                                              \[\leadsto \color{blue}{\mathsf{fma}\left(z + x, \frac{x - z}{y}, y\right) \cdot 0.5} \]
                                            6. Step-by-step derivation
                                              1. Applied rewrites99.6%

                                                \[\leadsto \left(\frac{x - z}{y} \cdot \left(z + x\right) + y\right) \cdot 0.5 \]
                                              2. Add Preprocessing

                                              Alternative 8: 99.9% accurate, 1.3× speedup?

                                              \[\begin{array}{l} x_m = \left|x\right| \\ y\_m = \left|y\right| \\ y\_s = \mathsf{copysign}\left(1, y\right) \\ y\_s \cdot \left(\mathsf{fma}\left(z + x\_m, \frac{x\_m - z}{y\_m}, y\_m\right) \cdot 0.5\right) \end{array} \]
                                              x_m = (fabs.f64 x)
                                              y\_m = (fabs.f64 y)
                                              y\_s = (copysign.f64 #s(literal 1 binary64) y)
                                              (FPCore (y_s x_m y_m z)
                                               :precision binary64
                                               (* y_s (* (fma (+ z x_m) (/ (- x_m z) y_m) y_m) 0.5)))
                                              x_m = fabs(x);
                                              y\_m = fabs(y);
                                              y\_s = copysign(1.0, y);
                                              double code(double y_s, double x_m, double y_m, double z) {
                                              	return y_s * (fma((z + x_m), ((x_m - z) / y_m), y_m) * 0.5);
                                              }
                                              
                                              x_m = abs(x)
                                              y\_m = abs(y)
                                              y\_s = copysign(1.0, y)
                                              function code(y_s, x_m, y_m, z)
                                              	return Float64(y_s * Float64(fma(Float64(z + x_m), Float64(Float64(x_m - z) / y_m), y_m) * 0.5))
                                              end
                                              
                                              x_m = N[Abs[x], $MachinePrecision]
                                              y\_m = N[Abs[y], $MachinePrecision]
                                              y\_s = N[With[{TMP1 = Abs[1.0], TMP2 = Sign[y]}, TMP1 * If[TMP2 == 0, 1, TMP2]], $MachinePrecision]
                                              code[y$95$s_, x$95$m_, y$95$m_, z_] := N[(y$95$s * N[(N[(N[(z + x$95$m), $MachinePrecision] * N[(N[(x$95$m - z), $MachinePrecision] / y$95$m), $MachinePrecision] + y$95$m), $MachinePrecision] * 0.5), $MachinePrecision]), $MachinePrecision]
                                              
                                              \begin{array}{l}
                                              x_m = \left|x\right|
                                              \\
                                              y\_m = \left|y\right|
                                              \\
                                              y\_s = \mathsf{copysign}\left(1, y\right)
                                              
                                              \\
                                              y\_s \cdot \left(\mathsf{fma}\left(z + x\_m, \frac{x\_m - z}{y\_m}, y\_m\right) \cdot 0.5\right)
                                              \end{array}
                                              
                                              Derivation
                                              1. Initial program 71.8%

                                                \[\frac{\left(x \cdot x + y \cdot y\right) - z \cdot z}{y \cdot 2} \]
                                              2. Add Preprocessing
                                              3. Taylor expanded in x around 0

                                                \[\leadsto \color{blue}{\frac{1}{2} \cdot \frac{{y}^{2} - {z}^{2}}{y} + \frac{1}{2} \cdot \frac{{x}^{2}}{y}} \]
                                              4. Step-by-step derivation
                                                1. distribute-lft-outN/A

                                                  \[\leadsto \color{blue}{\frac{1}{2} \cdot \left(\frac{{y}^{2} - {z}^{2}}{y} + \frac{{x}^{2}}{y}\right)} \]
                                                2. *-commutativeN/A

                                                  \[\leadsto \color{blue}{\left(\frac{{y}^{2} - {z}^{2}}{y} + \frac{{x}^{2}}{y}\right) \cdot \frac{1}{2}} \]
                                                3. lower-*.f64N/A

                                                  \[\leadsto \color{blue}{\left(\frac{{y}^{2} - {z}^{2}}{y} + \frac{{x}^{2}}{y}\right) \cdot \frac{1}{2}} \]
                                              5. Applied rewrites99.6%

                                                \[\leadsto \color{blue}{\mathsf{fma}\left(z + x, \frac{x - z}{y}, y\right) \cdot 0.5} \]
                                              6. Add Preprocessing

                                              Alternative 9: 34.6% accurate, 6.3× speedup?

                                              \[\begin{array}{l} x_m = \left|x\right| \\ y\_m = \left|y\right| \\ y\_s = \mathsf{copysign}\left(1, y\right) \\ y\_s \cdot \left(0.5 \cdot y\_m\right) \end{array} \]
                                              x_m = (fabs.f64 x)
                                              y\_m = (fabs.f64 y)
                                              y\_s = (copysign.f64 #s(literal 1 binary64) y)
                                              (FPCore (y_s x_m y_m z) :precision binary64 (* y_s (* 0.5 y_m)))
                                              x_m = fabs(x);
                                              y\_m = fabs(y);
                                              y\_s = copysign(1.0, y);
                                              double code(double y_s, double x_m, double y_m, double z) {
                                              	return y_s * (0.5 * y_m);
                                              }
                                              
                                              x_m = abs(x)
                                              y\_m = abs(y)
                                              y\_s = copysign(1.0d0, y)
                                              real(8) function code(y_s, x_m, y_m, z)
                                                  real(8), intent (in) :: y_s
                                                  real(8), intent (in) :: x_m
                                                  real(8), intent (in) :: y_m
                                                  real(8), intent (in) :: z
                                                  code = y_s * (0.5d0 * y_m)
                                              end function
                                              
                                              x_m = Math.abs(x);
                                              y\_m = Math.abs(y);
                                              y\_s = Math.copySign(1.0, y);
                                              public static double code(double y_s, double x_m, double y_m, double z) {
                                              	return y_s * (0.5 * y_m);
                                              }
                                              
                                              x_m = math.fabs(x)
                                              y\_m = math.fabs(y)
                                              y\_s = math.copysign(1.0, y)
                                              def code(y_s, x_m, y_m, z):
                                              	return y_s * (0.5 * y_m)
                                              
                                              x_m = abs(x)
                                              y\_m = abs(y)
                                              y\_s = copysign(1.0, y)
                                              function code(y_s, x_m, y_m, z)
                                              	return Float64(y_s * Float64(0.5 * y_m))
                                              end
                                              
                                              x_m = abs(x);
                                              y\_m = abs(y);
                                              y\_s = sign(y) * abs(1.0);
                                              function tmp = code(y_s, x_m, y_m, z)
                                              	tmp = y_s * (0.5 * y_m);
                                              end
                                              
                                              x_m = N[Abs[x], $MachinePrecision]
                                              y\_m = N[Abs[y], $MachinePrecision]
                                              y\_s = N[With[{TMP1 = Abs[1.0], TMP2 = Sign[y]}, TMP1 * If[TMP2 == 0, 1, TMP2]], $MachinePrecision]
                                              code[y$95$s_, x$95$m_, y$95$m_, z_] := N[(y$95$s * N[(0.5 * y$95$m), $MachinePrecision]), $MachinePrecision]
                                              
                                              \begin{array}{l}
                                              x_m = \left|x\right|
                                              \\
                                              y\_m = \left|y\right|
                                              \\
                                              y\_s = \mathsf{copysign}\left(1, y\right)
                                              
                                              \\
                                              y\_s \cdot \left(0.5 \cdot y\_m\right)
                                              \end{array}
                                              
                                              Derivation
                                              1. Initial program 71.8%

                                                \[\frac{\left(x \cdot x + y \cdot y\right) - z \cdot z}{y \cdot 2} \]
                                              2. Add Preprocessing
                                              3. Taylor expanded in y around inf

                                                \[\leadsto \color{blue}{\frac{1}{2} \cdot y} \]
                                              4. Step-by-step derivation
                                                1. lower-*.f6433.7

                                                  \[\leadsto \color{blue}{0.5 \cdot y} \]
                                              5. Applied rewrites33.7%

                                                \[\leadsto \color{blue}{0.5 \cdot y} \]
                                              6. Add Preprocessing

                                              Developer Target 1: 99.9% accurate, 1.1× speedup?

                                              \[\begin{array}{l} \\ y \cdot 0.5 - \left(\frac{0.5}{y} \cdot \left(z + x\right)\right) \cdot \left(z - x\right) \end{array} \]
                                              (FPCore (x y z)
                                               :precision binary64
                                               (- (* y 0.5) (* (* (/ 0.5 y) (+ z x)) (- z x))))
                                              double code(double x, double y, double z) {
                                              	return (y * 0.5) - (((0.5 / y) * (z + x)) * (z - x));
                                              }
                                              
                                              real(8) function code(x, y, z)
                                                  real(8), intent (in) :: x
                                                  real(8), intent (in) :: y
                                                  real(8), intent (in) :: z
                                                  code = (y * 0.5d0) - (((0.5d0 / y) * (z + x)) * (z - x))
                                              end function
                                              
                                              public static double code(double x, double y, double z) {
                                              	return (y * 0.5) - (((0.5 / y) * (z + x)) * (z - x));
                                              }
                                              
                                              def code(x, y, z):
                                              	return (y * 0.5) - (((0.5 / y) * (z + x)) * (z - x))
                                              
                                              function code(x, y, z)
                                              	return Float64(Float64(y * 0.5) - Float64(Float64(Float64(0.5 / y) * Float64(z + x)) * Float64(z - x)))
                                              end
                                              
                                              function tmp = code(x, y, z)
                                              	tmp = (y * 0.5) - (((0.5 / y) * (z + x)) * (z - x));
                                              end
                                              
                                              code[x_, y_, z_] := N[(N[(y * 0.5), $MachinePrecision] - N[(N[(N[(0.5 / y), $MachinePrecision] * N[(z + x), $MachinePrecision]), $MachinePrecision] * N[(z - x), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]
                                              
                                              \begin{array}{l}
                                              
                                              \\
                                              y \cdot 0.5 - \left(\frac{0.5}{y} \cdot \left(z + x\right)\right) \cdot \left(z - x\right)
                                              \end{array}
                                              

                                              Reproduce

                                              ?
                                              herbie shell --seed 2024320 
                                              (FPCore (x y z)
                                                :name "Diagrams.TwoD.Apollonian:initialConfig from diagrams-contrib-1.3.0.5, A"
                                                :precision binary64
                                              
                                                :alt
                                                (! :herbie-platform default (- (* y 1/2) (* (* (/ 1/2 y) (+ z x)) (- z x))))
                                              
                                                (/ (- (+ (* x x) (* y y)) (* z z)) (* y 2.0)))