Migdal et al, Equation (64)

Percentage Accurate: 99.6% → 99.6%
Time: 8.7s
Alternatives: 10
Speedup: 1.9×

Specification

?
\[\begin{array}{l} \\ \begin{array}{l} t_1 := \frac{\cos th}{\sqrt{2}}\\ t\_1 \cdot \left(a1 \cdot a1\right) + t\_1 \cdot \left(a2 \cdot a2\right) \end{array} \end{array} \]
(FPCore (a1 a2 th)
 :precision binary64
 (let* ((t_1 (/ (cos th) (sqrt 2.0))))
   (+ (* t_1 (* a1 a1)) (* t_1 (* a2 a2)))))
double code(double a1, double a2, double th) {
	double t_1 = cos(th) / sqrt(2.0);
	return (t_1 * (a1 * a1)) + (t_1 * (a2 * a2));
}
module fmin_fmax_functions
    implicit none
    private
    public fmax
    public fmin

    interface fmax
        module procedure fmax88
        module procedure fmax44
        module procedure fmax84
        module procedure fmax48
    end interface
    interface fmin
        module procedure fmin88
        module procedure fmin44
        module procedure fmin84
        module procedure fmin48
    end interface
contains
    real(8) function fmax88(x, y) result (res)
        real(8), intent (in) :: x
        real(8), intent (in) :: y
        res = merge(y, merge(x, max(x, y), y /= y), x /= x)
    end function
    real(4) function fmax44(x, y) result (res)
        real(4), intent (in) :: x
        real(4), intent (in) :: y
        res = merge(y, merge(x, max(x, y), y /= y), x /= x)
    end function
    real(8) function fmax84(x, y) result(res)
        real(8), intent (in) :: x
        real(4), intent (in) :: y
        res = merge(dble(y), merge(x, max(x, dble(y)), y /= y), x /= x)
    end function
    real(8) function fmax48(x, y) result(res)
        real(4), intent (in) :: x
        real(8), intent (in) :: y
        res = merge(y, merge(dble(x), max(dble(x), y), y /= y), x /= x)
    end function
    real(8) function fmin88(x, y) result (res)
        real(8), intent (in) :: x
        real(8), intent (in) :: y
        res = merge(y, merge(x, min(x, y), y /= y), x /= x)
    end function
    real(4) function fmin44(x, y) result (res)
        real(4), intent (in) :: x
        real(4), intent (in) :: y
        res = merge(y, merge(x, min(x, y), y /= y), x /= x)
    end function
    real(8) function fmin84(x, y) result(res)
        real(8), intent (in) :: x
        real(4), intent (in) :: y
        res = merge(dble(y), merge(x, min(x, dble(y)), y /= y), x /= x)
    end function
    real(8) function fmin48(x, y) result(res)
        real(4), intent (in) :: x
        real(8), intent (in) :: y
        res = merge(y, merge(dble(x), min(dble(x), y), y /= y), x /= x)
    end function
end module

real(8) function code(a1, a2, th)
use fmin_fmax_functions
    real(8), intent (in) :: a1
    real(8), intent (in) :: a2
    real(8), intent (in) :: th
    real(8) :: t_1
    t_1 = cos(th) / sqrt(2.0d0)
    code = (t_1 * (a1 * a1)) + (t_1 * (a2 * a2))
end function
public static double code(double a1, double a2, double th) {
	double t_1 = Math.cos(th) / Math.sqrt(2.0);
	return (t_1 * (a1 * a1)) + (t_1 * (a2 * a2));
}
def code(a1, a2, th):
	t_1 = math.cos(th) / math.sqrt(2.0)
	return (t_1 * (a1 * a1)) + (t_1 * (a2 * a2))
function code(a1, a2, th)
	t_1 = Float64(cos(th) / sqrt(2.0))
	return Float64(Float64(t_1 * Float64(a1 * a1)) + Float64(t_1 * Float64(a2 * a2)))
end
function tmp = code(a1, a2, th)
	t_1 = cos(th) / sqrt(2.0);
	tmp = (t_1 * (a1 * a1)) + (t_1 * (a2 * a2));
end
code[a1_, a2_, th_] := Block[{t$95$1 = N[(N[Cos[th], $MachinePrecision] / N[Sqrt[2.0], $MachinePrecision]), $MachinePrecision]}, N[(N[(t$95$1 * N[(a1 * a1), $MachinePrecision]), $MachinePrecision] + N[(t$95$1 * N[(a2 * a2), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]]
\begin{array}{l}

\\
\begin{array}{l}
t_1 := \frac{\cos th}{\sqrt{2}}\\
t\_1 \cdot \left(a1 \cdot a1\right) + t\_1 \cdot \left(a2 \cdot a2\right)
\end{array}
\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 10 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: 99.6% accurate, 1.0× speedup?

\[\begin{array}{l} \\ \begin{array}{l} t_1 := \frac{\cos th}{\sqrt{2}}\\ t\_1 \cdot \left(a1 \cdot a1\right) + t\_1 \cdot \left(a2 \cdot a2\right) \end{array} \end{array} \]
(FPCore (a1 a2 th)
 :precision binary64
 (let* ((t_1 (/ (cos th) (sqrt 2.0))))
   (+ (* t_1 (* a1 a1)) (* t_1 (* a2 a2)))))
double code(double a1, double a2, double th) {
	double t_1 = cos(th) / sqrt(2.0);
	return (t_1 * (a1 * a1)) + (t_1 * (a2 * a2));
}
module fmin_fmax_functions
    implicit none
    private
    public fmax
    public fmin

    interface fmax
        module procedure fmax88
        module procedure fmax44
        module procedure fmax84
        module procedure fmax48
    end interface
    interface fmin
        module procedure fmin88
        module procedure fmin44
        module procedure fmin84
        module procedure fmin48
    end interface
contains
    real(8) function fmax88(x, y) result (res)
        real(8), intent (in) :: x
        real(8), intent (in) :: y
        res = merge(y, merge(x, max(x, y), y /= y), x /= x)
    end function
    real(4) function fmax44(x, y) result (res)
        real(4), intent (in) :: x
        real(4), intent (in) :: y
        res = merge(y, merge(x, max(x, y), y /= y), x /= x)
    end function
    real(8) function fmax84(x, y) result(res)
        real(8), intent (in) :: x
        real(4), intent (in) :: y
        res = merge(dble(y), merge(x, max(x, dble(y)), y /= y), x /= x)
    end function
    real(8) function fmax48(x, y) result(res)
        real(4), intent (in) :: x
        real(8), intent (in) :: y
        res = merge(y, merge(dble(x), max(dble(x), y), y /= y), x /= x)
    end function
    real(8) function fmin88(x, y) result (res)
        real(8), intent (in) :: x
        real(8), intent (in) :: y
        res = merge(y, merge(x, min(x, y), y /= y), x /= x)
    end function
    real(4) function fmin44(x, y) result (res)
        real(4), intent (in) :: x
        real(4), intent (in) :: y
        res = merge(y, merge(x, min(x, y), y /= y), x /= x)
    end function
    real(8) function fmin84(x, y) result(res)
        real(8), intent (in) :: x
        real(4), intent (in) :: y
        res = merge(dble(y), merge(x, min(x, dble(y)), y /= y), x /= x)
    end function
    real(8) function fmin48(x, y) result(res)
        real(4), intent (in) :: x
        real(8), intent (in) :: y
        res = merge(y, merge(dble(x), min(dble(x), y), y /= y), x /= x)
    end function
end module

real(8) function code(a1, a2, th)
use fmin_fmax_functions
    real(8), intent (in) :: a1
    real(8), intent (in) :: a2
    real(8), intent (in) :: th
    real(8) :: t_1
    t_1 = cos(th) / sqrt(2.0d0)
    code = (t_1 * (a1 * a1)) + (t_1 * (a2 * a2))
end function
public static double code(double a1, double a2, double th) {
	double t_1 = Math.cos(th) / Math.sqrt(2.0);
	return (t_1 * (a1 * a1)) + (t_1 * (a2 * a2));
}
def code(a1, a2, th):
	t_1 = math.cos(th) / math.sqrt(2.0)
	return (t_1 * (a1 * a1)) + (t_1 * (a2 * a2))
function code(a1, a2, th)
	t_1 = Float64(cos(th) / sqrt(2.0))
	return Float64(Float64(t_1 * Float64(a1 * a1)) + Float64(t_1 * Float64(a2 * a2)))
end
function tmp = code(a1, a2, th)
	t_1 = cos(th) / sqrt(2.0);
	tmp = (t_1 * (a1 * a1)) + (t_1 * (a2 * a2));
end
code[a1_, a2_, th_] := Block[{t$95$1 = N[(N[Cos[th], $MachinePrecision] / N[Sqrt[2.0], $MachinePrecision]), $MachinePrecision]}, N[(N[(t$95$1 * N[(a1 * a1), $MachinePrecision]), $MachinePrecision] + N[(t$95$1 * N[(a2 * a2), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]]
\begin{array}{l}

\\
\begin{array}{l}
t_1 := \frac{\cos th}{\sqrt{2}}\\
t\_1 \cdot \left(a1 \cdot a1\right) + t\_1 \cdot \left(a2 \cdot a2\right)
\end{array}
\end{array}

Alternative 1: 99.6% accurate, 1.0× speedup?

\[\begin{array}{l} a2_m = \left|a2\right| \\ a1_m = \left|a1\right| \\ [a1_m, a2_m, th] = \mathsf{sort}([a1_m, a2_m, th])\\ \\ \mathsf{fma}\left(\frac{\cos th}{\sqrt{2}}, a1\_m \cdot a1\_m, \left(a2\_m \cdot \cos th\right) \cdot \frac{a2\_m}{\sqrt{2}}\right) \end{array} \]
a2_m = (fabs.f64 a2)
a1_m = (fabs.f64 a1)
NOTE: a1_m, a2_m, and th should be sorted in increasing order before calling this function.
(FPCore (a1_m a2_m th)
 :precision binary64
 (fma
  (/ (cos th) (sqrt 2.0))
  (* a1_m a1_m)
  (* (* a2_m (cos th)) (/ a2_m (sqrt 2.0)))))
a2_m = fabs(a2);
a1_m = fabs(a1);
assert(a1_m < a2_m && a2_m < th);
double code(double a1_m, double a2_m, double th) {
	return fma((cos(th) / sqrt(2.0)), (a1_m * a1_m), ((a2_m * cos(th)) * (a2_m / sqrt(2.0))));
}
a2_m = abs(a2)
a1_m = abs(a1)
a1_m, a2_m, th = sort([a1_m, a2_m, th])
function code(a1_m, a2_m, th)
	return fma(Float64(cos(th) / sqrt(2.0)), Float64(a1_m * a1_m), Float64(Float64(a2_m * cos(th)) * Float64(a2_m / sqrt(2.0))))
end
a2_m = N[Abs[a2], $MachinePrecision]
a1_m = N[Abs[a1], $MachinePrecision]
NOTE: a1_m, a2_m, and th should be sorted in increasing order before calling this function.
code[a1$95$m_, a2$95$m_, th_] := N[(N[(N[Cos[th], $MachinePrecision] / N[Sqrt[2.0], $MachinePrecision]), $MachinePrecision] * N[(a1$95$m * a1$95$m), $MachinePrecision] + N[(N[(a2$95$m * N[Cos[th], $MachinePrecision]), $MachinePrecision] * N[(a2$95$m / N[Sqrt[2.0], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
a2_m = \left|a2\right|
\\
a1_m = \left|a1\right|
\\
[a1_m, a2_m, th] = \mathsf{sort}([a1_m, a2_m, th])\\
\\
\mathsf{fma}\left(\frac{\cos th}{\sqrt{2}}, a1\_m \cdot a1\_m, \left(a2\_m \cdot \cos th\right) \cdot \frac{a2\_m}{\sqrt{2}}\right)
\end{array}
Derivation
  1. Initial program 99.5%

    \[\frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right) + \frac{\cos th}{\sqrt{2}} \cdot \left(a2 \cdot a2\right) \]
  2. Add Preprocessing
  3. Applied rewrites99.6%

    \[\leadsto \color{blue}{\mathsf{fma}\left(\frac{\cos th}{\sqrt{2}}, a1 \cdot a1, \left(a2 \cdot \cos th\right) \cdot \frac{a2}{\sqrt{2}}\right)} \]
  4. Add Preprocessing

Alternative 2: 80.1% accurate, 1.6× speedup?

\[\begin{array}{l} a2_m = \left|a2\right| \\ a1_m = \left|a1\right| \\ [a1_m, a2_m, th] = \mathsf{sort}([a1_m, a2_m, th])\\ \\ \begin{array}{l} t_1 := \mathsf{fma}\left(a1\_m, a1\_m, a2\_m \cdot a2\_m\right)\\ \mathbf{if}\;\frac{\cos th}{\sqrt{2}} \leq -0.005:\\ \;\;\;\;\frac{\sqrt{2} \cdot t\_1}{-2}\\ \mathbf{else}:\\ \;\;\;\;\left(0.5 \cdot \sqrt{2}\right) \cdot t\_1\\ \end{array} \end{array} \]
a2_m = (fabs.f64 a2)
a1_m = (fabs.f64 a1)
NOTE: a1_m, a2_m, and th should be sorted in increasing order before calling this function.
(FPCore (a1_m a2_m th)
 :precision binary64
 (let* ((t_1 (fma a1_m a1_m (* a2_m a2_m))))
   (if (<= (/ (cos th) (sqrt 2.0)) -0.005)
     (/ (* (sqrt 2.0) t_1) (- 2.0))
     (* (* 0.5 (sqrt 2.0)) t_1))))
a2_m = fabs(a2);
a1_m = fabs(a1);
assert(a1_m < a2_m && a2_m < th);
double code(double a1_m, double a2_m, double th) {
	double t_1 = fma(a1_m, a1_m, (a2_m * a2_m));
	double tmp;
	if ((cos(th) / sqrt(2.0)) <= -0.005) {
		tmp = (sqrt(2.0) * t_1) / -2.0;
	} else {
		tmp = (0.5 * sqrt(2.0)) * t_1;
	}
	return tmp;
}
a2_m = abs(a2)
a1_m = abs(a1)
a1_m, a2_m, th = sort([a1_m, a2_m, th])
function code(a1_m, a2_m, th)
	t_1 = fma(a1_m, a1_m, Float64(a2_m * a2_m))
	tmp = 0.0
	if (Float64(cos(th) / sqrt(2.0)) <= -0.005)
		tmp = Float64(Float64(sqrt(2.0) * t_1) / Float64(-2.0));
	else
		tmp = Float64(Float64(0.5 * sqrt(2.0)) * t_1);
	end
	return tmp
end
a2_m = N[Abs[a2], $MachinePrecision]
a1_m = N[Abs[a1], $MachinePrecision]
NOTE: a1_m, a2_m, and th should be sorted in increasing order before calling this function.
code[a1$95$m_, a2$95$m_, th_] := Block[{t$95$1 = N[(a1$95$m * a1$95$m + N[(a2$95$m * a2$95$m), $MachinePrecision]), $MachinePrecision]}, If[LessEqual[N[(N[Cos[th], $MachinePrecision] / N[Sqrt[2.0], $MachinePrecision]), $MachinePrecision], -0.005], N[(N[(N[Sqrt[2.0], $MachinePrecision] * t$95$1), $MachinePrecision] / (-2.0)), $MachinePrecision], N[(N[(0.5 * N[Sqrt[2.0], $MachinePrecision]), $MachinePrecision] * t$95$1), $MachinePrecision]]]
\begin{array}{l}
a2_m = \left|a2\right|
\\
a1_m = \left|a1\right|
\\
[a1_m, a2_m, th] = \mathsf{sort}([a1_m, a2_m, th])\\
\\
\begin{array}{l}
t_1 := \mathsf{fma}\left(a1\_m, a1\_m, a2\_m \cdot a2\_m\right)\\
\mathbf{if}\;\frac{\cos th}{\sqrt{2}} \leq -0.005:\\
\;\;\;\;\frac{\sqrt{2} \cdot t\_1}{-2}\\

\mathbf{else}:\\
\;\;\;\;\left(0.5 \cdot \sqrt{2}\right) \cdot t\_1\\


\end{array}
\end{array}
Derivation
  1. Split input into 2 regimes
  2. if (/.f64 (cos.f64 th) (sqrt.f64 #s(literal 2 binary64))) < -0.0050000000000000001

    1. Initial program 99.6%

      \[\frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right) + \frac{\cos th}{\sqrt{2}} \cdot \left(a2 \cdot a2\right) \]
    2. Add Preprocessing
    3. Taylor expanded in th around 0

      \[\leadsto \color{blue}{\frac{{a1}^{2}}{\sqrt{2}} + \frac{{a2}^{2}}{\sqrt{2}}} \]
    4. Step-by-step derivation
      1. div-add-revN/A

        \[\leadsto \color{blue}{\frac{{a1}^{2} + {a2}^{2}}{\sqrt{2}}} \]
      2. lower-/.f64N/A

        \[\leadsto \color{blue}{\frac{{a1}^{2} + {a2}^{2}}{\sqrt{2}}} \]
      3. +-commutativeN/A

        \[\leadsto \frac{\color{blue}{{a2}^{2} + {a1}^{2}}}{\sqrt{2}} \]
      4. unpow2N/A

        \[\leadsto \frac{\color{blue}{a2 \cdot a2} + {a1}^{2}}{\sqrt{2}} \]
      5. lower-fma.f64N/A

        \[\leadsto \frac{\color{blue}{\mathsf{fma}\left(a2, a2, {a1}^{2}\right)}}{\sqrt{2}} \]
      6. unpow2N/A

        \[\leadsto \frac{\mathsf{fma}\left(a2, a2, \color{blue}{a1 \cdot a1}\right)}{\sqrt{2}} \]
      7. lower-*.f64N/A

        \[\leadsto \frac{\mathsf{fma}\left(a2, a2, \color{blue}{a1 \cdot a1}\right)}{\sqrt{2}} \]
      8. lower-sqrt.f645.6

        \[\leadsto \frac{\mathsf{fma}\left(a2, a2, a1 \cdot a1\right)}{\color{blue}{\sqrt{2}}} \]
    5. Applied rewrites5.6%

      \[\leadsto \color{blue}{\frac{\mathsf{fma}\left(a2, a2, a1 \cdot a1\right)}{\sqrt{2}}} \]
    6. Step-by-step derivation
      1. Applied rewrites5.6%

        \[\leadsto \frac{\mathsf{fma}\left(a2 \cdot a2, \sqrt{2}, \sqrt{2} \cdot \left(a1 \cdot a1\right)\right)}{\color{blue}{2}} \]
      2. Applied rewrites64.7%

        \[\leadsto \frac{-\sqrt{2} \cdot \mathsf{fma}\left(a1, a1, a2 \cdot a2\right)}{2} \]

      if -0.0050000000000000001 < (/.f64 (cos.f64 th) (sqrt.f64 #s(literal 2 binary64)))

      1. Initial program 99.5%

        \[\frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right) + \frac{\cos th}{\sqrt{2}} \cdot \left(a2 \cdot a2\right) \]
      2. Add Preprocessing
      3. Taylor expanded in th around 0

        \[\leadsto \color{blue}{\frac{{a1}^{2}}{\sqrt{2}} + \frac{{a2}^{2}}{\sqrt{2}}} \]
      4. Step-by-step derivation
        1. div-add-revN/A

          \[\leadsto \color{blue}{\frac{{a1}^{2} + {a2}^{2}}{\sqrt{2}}} \]
        2. lower-/.f64N/A

          \[\leadsto \color{blue}{\frac{{a1}^{2} + {a2}^{2}}{\sqrt{2}}} \]
        3. +-commutativeN/A

          \[\leadsto \frac{\color{blue}{{a2}^{2} + {a1}^{2}}}{\sqrt{2}} \]
        4. unpow2N/A

          \[\leadsto \frac{\color{blue}{a2 \cdot a2} + {a1}^{2}}{\sqrt{2}} \]
        5. lower-fma.f64N/A

          \[\leadsto \frac{\color{blue}{\mathsf{fma}\left(a2, a2, {a1}^{2}\right)}}{\sqrt{2}} \]
        6. unpow2N/A

          \[\leadsto \frac{\mathsf{fma}\left(a2, a2, \color{blue}{a1 \cdot a1}\right)}{\sqrt{2}} \]
        7. lower-*.f64N/A

          \[\leadsto \frac{\mathsf{fma}\left(a2, a2, \color{blue}{a1 \cdot a1}\right)}{\sqrt{2}} \]
        8. lower-sqrt.f6485.8

          \[\leadsto \frac{\mathsf{fma}\left(a2, a2, a1 \cdot a1\right)}{\color{blue}{\sqrt{2}}} \]
      5. Applied rewrites85.8%

        \[\leadsto \color{blue}{\frac{\mathsf{fma}\left(a2, a2, a1 \cdot a1\right)}{\sqrt{2}}} \]
      6. Step-by-step derivation
        1. Applied rewrites85.8%

          \[\leadsto \frac{\mathsf{fma}\left(a2 \cdot a2, \sqrt{2}, \sqrt{2} \cdot \left(a1 \cdot a1\right)\right)}{\color{blue}{2}} \]
        2. Taylor expanded in a1 around 0

          \[\leadsto \frac{1}{2} \cdot \left({a1}^{2} \cdot \sqrt{2}\right) + \color{blue}{\frac{1}{2} \cdot \left({a2}^{2} \cdot \sqrt{2}\right)} \]
        3. Step-by-step derivation
          1. Applied rewrites85.8%

            \[\leadsto \left(0.5 \cdot \sqrt{2}\right) \cdot \color{blue}{\mathsf{fma}\left(a1, a1, a2 \cdot a2\right)} \]
        4. Recombined 2 regimes into one program.
        5. Final simplification80.5%

          \[\leadsto \begin{array}{l} \mathbf{if}\;\frac{\cos th}{\sqrt{2}} \leq -0.005:\\ \;\;\;\;\frac{\sqrt{2} \cdot \mathsf{fma}\left(a1, a1, a2 \cdot a2\right)}{-2}\\ \mathbf{else}:\\ \;\;\;\;\left(0.5 \cdot \sqrt{2}\right) \cdot \mathsf{fma}\left(a1, a1, a2 \cdot a2\right)\\ \end{array} \]
        6. Add Preprocessing

        Alternative 3: 99.6% accurate, 1.7× speedup?

        \[\begin{array}{l} a2_m = \left|a2\right| \\ a1_m = \left|a1\right| \\ [a1_m, a2_m, th] = \mathsf{sort}([a1_m, a2_m, th])\\ \\ \cos th \cdot \mathsf{fma}\left(a2\_m, \frac{a2\_m}{\sqrt{2}}, a1\_m \cdot \frac{a1\_m}{\sqrt{2}}\right) \end{array} \]
        a2_m = (fabs.f64 a2)
        a1_m = (fabs.f64 a1)
        NOTE: a1_m, a2_m, and th should be sorted in increasing order before calling this function.
        (FPCore (a1_m a2_m th)
         :precision binary64
         (* (cos th) (fma a2_m (/ a2_m (sqrt 2.0)) (* a1_m (/ a1_m (sqrt 2.0))))))
        a2_m = fabs(a2);
        a1_m = fabs(a1);
        assert(a1_m < a2_m && a2_m < th);
        double code(double a1_m, double a2_m, double th) {
        	return cos(th) * fma(a2_m, (a2_m / sqrt(2.0)), (a1_m * (a1_m / sqrt(2.0))));
        }
        
        a2_m = abs(a2)
        a1_m = abs(a1)
        a1_m, a2_m, th = sort([a1_m, a2_m, th])
        function code(a1_m, a2_m, th)
        	return Float64(cos(th) * fma(a2_m, Float64(a2_m / sqrt(2.0)), Float64(a1_m * Float64(a1_m / sqrt(2.0)))))
        end
        
        a2_m = N[Abs[a2], $MachinePrecision]
        a1_m = N[Abs[a1], $MachinePrecision]
        NOTE: a1_m, a2_m, and th should be sorted in increasing order before calling this function.
        code[a1$95$m_, a2$95$m_, th_] := N[(N[Cos[th], $MachinePrecision] * N[(a2$95$m * N[(a2$95$m / N[Sqrt[2.0], $MachinePrecision]), $MachinePrecision] + N[(a1$95$m * N[(a1$95$m / N[Sqrt[2.0], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]
        
        \begin{array}{l}
        a2_m = \left|a2\right|
        \\
        a1_m = \left|a1\right|
        \\
        [a1_m, a2_m, th] = \mathsf{sort}([a1_m, a2_m, th])\\
        \\
        \cos th \cdot \mathsf{fma}\left(a2\_m, \frac{a2\_m}{\sqrt{2}}, a1\_m \cdot \frac{a1\_m}{\sqrt{2}}\right)
        \end{array}
        
        Derivation
        1. Initial program 99.5%

          \[\frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right) + \frac{\cos th}{\sqrt{2}} \cdot \left(a2 \cdot a2\right) \]
        2. Add Preprocessing
        3. Step-by-step derivation
          1. lift-+.f64N/A

            \[\leadsto \color{blue}{\frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right) + \frac{\cos th}{\sqrt{2}} \cdot \left(a2 \cdot a2\right)} \]
          2. +-commutativeN/A

            \[\leadsto \color{blue}{\frac{\cos th}{\sqrt{2}} \cdot \left(a2 \cdot a2\right) + \frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right)} \]
          3. lift-*.f64N/A

            \[\leadsto \color{blue}{\frac{\cos th}{\sqrt{2}} \cdot \left(a2 \cdot a2\right)} + \frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right) \]
          4. lift-/.f64N/A

            \[\leadsto \color{blue}{\frac{\cos th}{\sqrt{2}}} \cdot \left(a2 \cdot a2\right) + \frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right) \]
          5. associate-*l/N/A

            \[\leadsto \color{blue}{\frac{\cos th \cdot \left(a2 \cdot a2\right)}{\sqrt{2}}} + \frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right) \]
          6. associate-/l*N/A

            \[\leadsto \color{blue}{\cos th \cdot \frac{a2 \cdot a2}{\sqrt{2}}} + \frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right) \]
          7. lift-*.f64N/A

            \[\leadsto \cos th \cdot \frac{a2 \cdot a2}{\sqrt{2}} + \color{blue}{\frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right)} \]
          8. lift-/.f64N/A

            \[\leadsto \cos th \cdot \frac{a2 \cdot a2}{\sqrt{2}} + \color{blue}{\frac{\cos th}{\sqrt{2}}} \cdot \left(a1 \cdot a1\right) \]
          9. associate-*l/N/A

            \[\leadsto \cos th \cdot \frac{a2 \cdot a2}{\sqrt{2}} + \color{blue}{\frac{\cos th \cdot \left(a1 \cdot a1\right)}{\sqrt{2}}} \]
          10. associate-/l*N/A

            \[\leadsto \cos th \cdot \frac{a2 \cdot a2}{\sqrt{2}} + \color{blue}{\cos th \cdot \frac{a1 \cdot a1}{\sqrt{2}}} \]
          11. distribute-lft-outN/A

            \[\leadsto \color{blue}{\cos th \cdot \left(\frac{a2 \cdot a2}{\sqrt{2}} + \frac{a1 \cdot a1}{\sqrt{2}}\right)} \]
          12. lower-*.f64N/A

            \[\leadsto \color{blue}{\cos th \cdot \left(\frac{a2 \cdot a2}{\sqrt{2}} + \frac{a1 \cdot a1}{\sqrt{2}}\right)} \]
          13. lift-*.f64N/A

            \[\leadsto \cos th \cdot \left(\frac{\color{blue}{a2 \cdot a2}}{\sqrt{2}} + \frac{a1 \cdot a1}{\sqrt{2}}\right) \]
          14. associate-/l*N/A

            \[\leadsto \cos th \cdot \left(\color{blue}{a2 \cdot \frac{a2}{\sqrt{2}}} + \frac{a1 \cdot a1}{\sqrt{2}}\right) \]
          15. lower-fma.f64N/A

            \[\leadsto \cos th \cdot \color{blue}{\mathsf{fma}\left(a2, \frac{a2}{\sqrt{2}}, \frac{a1 \cdot a1}{\sqrt{2}}\right)} \]
          16. lower-/.f64N/A

            \[\leadsto \cos th \cdot \mathsf{fma}\left(a2, \color{blue}{\frac{a2}{\sqrt{2}}}, \frac{a1 \cdot a1}{\sqrt{2}}\right) \]
          17. lift-*.f64N/A

            \[\leadsto \cos th \cdot \mathsf{fma}\left(a2, \frac{a2}{\sqrt{2}}, \frac{\color{blue}{a1 \cdot a1}}{\sqrt{2}}\right) \]
        4. Applied rewrites99.6%

          \[\leadsto \color{blue}{\cos th \cdot \mathsf{fma}\left(a2, \frac{a2}{\sqrt{2}}, a1 \cdot \frac{a1}{\sqrt{2}}\right)} \]
        5. Add Preprocessing

        Alternative 4: 99.6% accurate, 1.9× speedup?

        \[\begin{array}{l} a2_m = \left|a2\right| \\ a1_m = \left|a1\right| \\ [a1_m, a2_m, th] = \mathsf{sort}([a1_m, a2_m, th])\\ \\ \left(\left(0.5 \cdot \cos th\right) \cdot \mathsf{fma}\left(a2\_m, a2\_m, a1\_m \cdot a1\_m\right)\right) \cdot \sqrt{2} \end{array} \]
        a2_m = (fabs.f64 a2)
        a1_m = (fabs.f64 a1)
        NOTE: a1_m, a2_m, and th should be sorted in increasing order before calling this function.
        (FPCore (a1_m a2_m th)
         :precision binary64
         (* (* (* 0.5 (cos th)) (fma a2_m a2_m (* a1_m a1_m))) (sqrt 2.0)))
        a2_m = fabs(a2);
        a1_m = fabs(a1);
        assert(a1_m < a2_m && a2_m < th);
        double code(double a1_m, double a2_m, double th) {
        	return ((0.5 * cos(th)) * fma(a2_m, a2_m, (a1_m * a1_m))) * sqrt(2.0);
        }
        
        a2_m = abs(a2)
        a1_m = abs(a1)
        a1_m, a2_m, th = sort([a1_m, a2_m, th])
        function code(a1_m, a2_m, th)
        	return Float64(Float64(Float64(0.5 * cos(th)) * fma(a2_m, a2_m, Float64(a1_m * a1_m))) * sqrt(2.0))
        end
        
        a2_m = N[Abs[a2], $MachinePrecision]
        a1_m = N[Abs[a1], $MachinePrecision]
        NOTE: a1_m, a2_m, and th should be sorted in increasing order before calling this function.
        code[a1$95$m_, a2$95$m_, th_] := N[(N[(N[(0.5 * N[Cos[th], $MachinePrecision]), $MachinePrecision] * N[(a2$95$m * a2$95$m + N[(a1$95$m * a1$95$m), $MachinePrecision]), $MachinePrecision]), $MachinePrecision] * N[Sqrt[2.0], $MachinePrecision]), $MachinePrecision]
        
        \begin{array}{l}
        a2_m = \left|a2\right|
        \\
        a1_m = \left|a1\right|
        \\
        [a1_m, a2_m, th] = \mathsf{sort}([a1_m, a2_m, th])\\
        \\
        \left(\left(0.5 \cdot \cos th\right) \cdot \mathsf{fma}\left(a2\_m, a2\_m, a1\_m \cdot a1\_m\right)\right) \cdot \sqrt{2}
        \end{array}
        
        Derivation
        1. Initial program 99.5%

          \[\frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right) + \frac{\cos th}{\sqrt{2}} \cdot \left(a2 \cdot a2\right) \]
        2. Add Preprocessing
        3. Step-by-step derivation
          1. lift-+.f64N/A

            \[\leadsto \color{blue}{\frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right) + \frac{\cos th}{\sqrt{2}} \cdot \left(a2 \cdot a2\right)} \]
          2. +-commutativeN/A

            \[\leadsto \color{blue}{\frac{\cos th}{\sqrt{2}} \cdot \left(a2 \cdot a2\right) + \frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right)} \]
          3. lift-*.f64N/A

            \[\leadsto \color{blue}{\frac{\cos th}{\sqrt{2}} \cdot \left(a2 \cdot a2\right)} + \frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right) \]
          4. lift-/.f64N/A

            \[\leadsto \color{blue}{\frac{\cos th}{\sqrt{2}}} \cdot \left(a2 \cdot a2\right) + \frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right) \]
          5. associate-*l/N/A

            \[\leadsto \color{blue}{\frac{\cos th \cdot \left(a2 \cdot a2\right)}{\sqrt{2}}} + \frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right) \]
          6. associate-/l*N/A

            \[\leadsto \color{blue}{\cos th \cdot \frac{a2 \cdot a2}{\sqrt{2}}} + \frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right) \]
          7. lift-*.f64N/A

            \[\leadsto \cos th \cdot \frac{a2 \cdot a2}{\sqrt{2}} + \color{blue}{\frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right)} \]
          8. lift-/.f64N/A

            \[\leadsto \cos th \cdot \frac{a2 \cdot a2}{\sqrt{2}} + \color{blue}{\frac{\cos th}{\sqrt{2}}} \cdot \left(a1 \cdot a1\right) \]
          9. associate-*l/N/A

            \[\leadsto \cos th \cdot \frac{a2 \cdot a2}{\sqrt{2}} + \color{blue}{\frac{\cos th \cdot \left(a1 \cdot a1\right)}{\sqrt{2}}} \]
          10. associate-/l*N/A

            \[\leadsto \cos th \cdot \frac{a2 \cdot a2}{\sqrt{2}} + \color{blue}{\cos th \cdot \frac{a1 \cdot a1}{\sqrt{2}}} \]
          11. distribute-lft-outN/A

            \[\leadsto \color{blue}{\cos th \cdot \left(\frac{a2 \cdot a2}{\sqrt{2}} + \frac{a1 \cdot a1}{\sqrt{2}}\right)} \]
          12. lower-*.f64N/A

            \[\leadsto \color{blue}{\cos th \cdot \left(\frac{a2 \cdot a2}{\sqrt{2}} + \frac{a1 \cdot a1}{\sqrt{2}}\right)} \]
          13. lift-*.f64N/A

            \[\leadsto \cos th \cdot \left(\frac{\color{blue}{a2 \cdot a2}}{\sqrt{2}} + \frac{a1 \cdot a1}{\sqrt{2}}\right) \]
          14. associate-/l*N/A

            \[\leadsto \cos th \cdot \left(\color{blue}{a2 \cdot \frac{a2}{\sqrt{2}}} + \frac{a1 \cdot a1}{\sqrt{2}}\right) \]
          15. lower-fma.f64N/A

            \[\leadsto \cos th \cdot \color{blue}{\mathsf{fma}\left(a2, \frac{a2}{\sqrt{2}}, \frac{a1 \cdot a1}{\sqrt{2}}\right)} \]
          16. lower-/.f64N/A

            \[\leadsto \cos th \cdot \mathsf{fma}\left(a2, \color{blue}{\frac{a2}{\sqrt{2}}}, \frac{a1 \cdot a1}{\sqrt{2}}\right) \]
          17. lift-*.f64N/A

            \[\leadsto \cos th \cdot \mathsf{fma}\left(a2, \frac{a2}{\sqrt{2}}, \frac{\color{blue}{a1 \cdot a1}}{\sqrt{2}}\right) \]
        4. Applied rewrites99.6%

          \[\leadsto \color{blue}{\cos th \cdot \mathsf{fma}\left(a2, \frac{a2}{\sqrt{2}}, a1 \cdot \frac{a1}{\sqrt{2}}\right)} \]
        5. Step-by-step derivation
          1. lift-fma.f64N/A

            \[\leadsto \cos th \cdot \color{blue}{\left(a2 \cdot \frac{a2}{\sqrt{2}} + a1 \cdot \frac{a1}{\sqrt{2}}\right)} \]
          2. lift-/.f64N/A

            \[\leadsto \cos th \cdot \left(a2 \cdot \color{blue}{\frac{a2}{\sqrt{2}}} + a1 \cdot \frac{a1}{\sqrt{2}}\right) \]
          3. associate-*r/N/A

            \[\leadsto \cos th \cdot \left(\color{blue}{\frac{a2 \cdot a2}{\sqrt{2}}} + a1 \cdot \frac{a1}{\sqrt{2}}\right) \]
          4. lift-*.f64N/A

            \[\leadsto \cos th \cdot \left(\frac{\color{blue}{a2 \cdot a2}}{\sqrt{2}} + a1 \cdot \frac{a1}{\sqrt{2}}\right) \]
          5. lift-*.f64N/A

            \[\leadsto \cos th \cdot \left(\frac{a2 \cdot a2}{\sqrt{2}} + \color{blue}{a1 \cdot \frac{a1}{\sqrt{2}}}\right) \]
          6. lift-/.f64N/A

            \[\leadsto \cos th \cdot \left(\frac{a2 \cdot a2}{\sqrt{2}} + a1 \cdot \color{blue}{\frac{a1}{\sqrt{2}}}\right) \]
          7. associate-*r/N/A

            \[\leadsto \cos th \cdot \left(\frac{a2 \cdot a2}{\sqrt{2}} + \color{blue}{\frac{a1 \cdot a1}{\sqrt{2}}}\right) \]
          8. lift-*.f64N/A

            \[\leadsto \cos th \cdot \left(\frac{a2 \cdot a2}{\sqrt{2}} + \frac{\color{blue}{a1 \cdot a1}}{\sqrt{2}}\right) \]
          9. frac-addN/A

            \[\leadsto \cos th \cdot \color{blue}{\frac{\left(a2 \cdot a2\right) \cdot \sqrt{2} + \sqrt{2} \cdot \left(a1 \cdot a1\right)}{\sqrt{2} \cdot \sqrt{2}}} \]
          10. lift-*.f64N/A

            \[\leadsto \cos th \cdot \frac{\left(a2 \cdot a2\right) \cdot \sqrt{2} + \color{blue}{\sqrt{2} \cdot \left(a1 \cdot a1\right)}}{\sqrt{2} \cdot \sqrt{2}} \]
          11. lift-fma.f64N/A

            \[\leadsto \cos th \cdot \frac{\color{blue}{\mathsf{fma}\left(a2 \cdot a2, \sqrt{2}, \sqrt{2} \cdot \left(a1 \cdot a1\right)\right)}}{\sqrt{2} \cdot \sqrt{2}} \]
          12. lift-sqrt.f64N/A

            \[\leadsto \cos th \cdot \frac{\mathsf{fma}\left(a2 \cdot a2, \sqrt{2}, \sqrt{2} \cdot \left(a1 \cdot a1\right)\right)}{\color{blue}{\sqrt{2}} \cdot \sqrt{2}} \]
          13. lift-sqrt.f64N/A

            \[\leadsto \cos th \cdot \frac{\mathsf{fma}\left(a2 \cdot a2, \sqrt{2}, \sqrt{2} \cdot \left(a1 \cdot a1\right)\right)}{\sqrt{2} \cdot \color{blue}{\sqrt{2}}} \]
          14. rem-square-sqrtN/A

            \[\leadsto \cos th \cdot \frac{\mathsf{fma}\left(a2 \cdot a2, \sqrt{2}, \sqrt{2} \cdot \left(a1 \cdot a1\right)\right)}{\color{blue}{2}} \]
          15. lift-/.f6499.6

            \[\leadsto \cos th \cdot \color{blue}{\frac{\mathsf{fma}\left(a2 \cdot a2, \sqrt{2}, \sqrt{2} \cdot \left(a1 \cdot a1\right)\right)}{2}} \]
        6. Applied rewrites99.6%

          \[\leadsto \cos th \cdot \color{blue}{\frac{\sqrt{2} \cdot \mathsf{fma}\left(a1, a1, a2 \cdot a2\right)}{2}} \]
        7. Taylor expanded in a1 around 0

          \[\leadsto \color{blue}{\frac{1}{2} \cdot \left({a1}^{2} \cdot \left(\cos th \cdot \sqrt{2}\right)\right) + \frac{1}{2} \cdot \left({a2}^{2} \cdot \left(\cos th \cdot \sqrt{2}\right)\right)} \]
        8. Step-by-step derivation
          1. +-commutativeN/A

            \[\leadsto \color{blue}{\frac{1}{2} \cdot \left({a2}^{2} \cdot \left(\cos th \cdot \sqrt{2}\right)\right) + \frac{1}{2} \cdot \left({a1}^{2} \cdot \left(\cos th \cdot \sqrt{2}\right)\right)} \]
          2. distribute-lft-outN/A

            \[\leadsto \color{blue}{\frac{1}{2} \cdot \left({a2}^{2} \cdot \left(\cos th \cdot \sqrt{2}\right) + {a1}^{2} \cdot \left(\cos th \cdot \sqrt{2}\right)\right)} \]
          3. distribute-rgt-inN/A

            \[\leadsto \frac{1}{2} \cdot \color{blue}{\left(\left(\cos th \cdot \sqrt{2}\right) \cdot \left({a2}^{2} + {a1}^{2}\right)\right)} \]
          4. +-commutativeN/A

            \[\leadsto \frac{1}{2} \cdot \left(\left(\cos th \cdot \sqrt{2}\right) \cdot \color{blue}{\left({a1}^{2} + {a2}^{2}\right)}\right) \]
          5. associate-*r*N/A

            \[\leadsto \frac{1}{2} \cdot \color{blue}{\left(\cos th \cdot \left(\sqrt{2} \cdot \left({a1}^{2} + {a2}^{2}\right)\right)\right)} \]
          6. associate-*r*N/A

            \[\leadsto \color{blue}{\left(\frac{1}{2} \cdot \cos th\right) \cdot \left(\sqrt{2} \cdot \left({a1}^{2} + {a2}^{2}\right)\right)} \]
          7. *-commutativeN/A

            \[\leadsto \left(\frac{1}{2} \cdot \cos th\right) \cdot \color{blue}{\left(\left({a1}^{2} + {a2}^{2}\right) \cdot \sqrt{2}\right)} \]
          8. associate-*r*N/A

            \[\leadsto \color{blue}{\left(\left(\frac{1}{2} \cdot \cos th\right) \cdot \left({a1}^{2} + {a2}^{2}\right)\right) \cdot \sqrt{2}} \]
          9. lower-*.f64N/A

            \[\leadsto \color{blue}{\left(\left(\frac{1}{2} \cdot \cos th\right) \cdot \left({a1}^{2} + {a2}^{2}\right)\right) \cdot \sqrt{2}} \]
        9. Applied rewrites99.6%

          \[\leadsto \color{blue}{\left(\left(0.5 \cdot \cos th\right) \cdot \mathsf{fma}\left(a2, a2, a1 \cdot a1\right)\right) \cdot \sqrt{2}} \]
        10. Add Preprocessing

        Alternative 5: 99.1% accurate, 2.0× speedup?

        \[\begin{array}{l} a2_m = \left|a2\right| \\ a1_m = \left|a1\right| \\ [a1_m, a2_m, th] = \mathsf{sort}([a1_m, a2_m, th])\\ \\ \left(\cos th \cdot a2\_m\right) \cdot \frac{a2\_m}{\sqrt{2}} \end{array} \]
        a2_m = (fabs.f64 a2)
        a1_m = (fabs.f64 a1)
        NOTE: a1_m, a2_m, and th should be sorted in increasing order before calling this function.
        (FPCore (a1_m a2_m th)
         :precision binary64
         (* (* (cos th) a2_m) (/ a2_m (sqrt 2.0))))
        a2_m = fabs(a2);
        a1_m = fabs(a1);
        assert(a1_m < a2_m && a2_m < th);
        double code(double a1_m, double a2_m, double th) {
        	return (cos(th) * a2_m) * (a2_m / sqrt(2.0));
        }
        
        a2_m =     private
        a1_m =     private
        NOTE: a1_m, a2_m, and th should be sorted in increasing order before calling this function.
        module fmin_fmax_functions
            implicit none
            private
            public fmax
            public fmin
        
            interface fmax
                module procedure fmax88
                module procedure fmax44
                module procedure fmax84
                module procedure fmax48
            end interface
            interface fmin
                module procedure fmin88
                module procedure fmin44
                module procedure fmin84
                module procedure fmin48
            end interface
        contains
            real(8) function fmax88(x, y) result (res)
                real(8), intent (in) :: x
                real(8), intent (in) :: y
                res = merge(y, merge(x, max(x, y), y /= y), x /= x)
            end function
            real(4) function fmax44(x, y) result (res)
                real(4), intent (in) :: x
                real(4), intent (in) :: y
                res = merge(y, merge(x, max(x, y), y /= y), x /= x)
            end function
            real(8) function fmax84(x, y) result(res)
                real(8), intent (in) :: x
                real(4), intent (in) :: y
                res = merge(dble(y), merge(x, max(x, dble(y)), y /= y), x /= x)
            end function
            real(8) function fmax48(x, y) result(res)
                real(4), intent (in) :: x
                real(8), intent (in) :: y
                res = merge(y, merge(dble(x), max(dble(x), y), y /= y), x /= x)
            end function
            real(8) function fmin88(x, y) result (res)
                real(8), intent (in) :: x
                real(8), intent (in) :: y
                res = merge(y, merge(x, min(x, y), y /= y), x /= x)
            end function
            real(4) function fmin44(x, y) result (res)
                real(4), intent (in) :: x
                real(4), intent (in) :: y
                res = merge(y, merge(x, min(x, y), y /= y), x /= x)
            end function
            real(8) function fmin84(x, y) result(res)
                real(8), intent (in) :: x
                real(4), intent (in) :: y
                res = merge(dble(y), merge(x, min(x, dble(y)), y /= y), x /= x)
            end function
            real(8) function fmin48(x, y) result(res)
                real(4), intent (in) :: x
                real(8), intent (in) :: y
                res = merge(y, merge(dble(x), min(dble(x), y), y /= y), x /= x)
            end function
        end module
        
        real(8) function code(a1_m, a2_m, th)
        use fmin_fmax_functions
            real(8), intent (in) :: a1_m
            real(8), intent (in) :: a2_m
            real(8), intent (in) :: th
            code = (cos(th) * a2_m) * (a2_m / sqrt(2.0d0))
        end function
        
        a2_m = Math.abs(a2);
        a1_m = Math.abs(a1);
        assert a1_m < a2_m && a2_m < th;
        public static double code(double a1_m, double a2_m, double th) {
        	return (Math.cos(th) * a2_m) * (a2_m / Math.sqrt(2.0));
        }
        
        a2_m = math.fabs(a2)
        a1_m = math.fabs(a1)
        [a1_m, a2_m, th] = sort([a1_m, a2_m, th])
        def code(a1_m, a2_m, th):
        	return (math.cos(th) * a2_m) * (a2_m / math.sqrt(2.0))
        
        a2_m = abs(a2)
        a1_m = abs(a1)
        a1_m, a2_m, th = sort([a1_m, a2_m, th])
        function code(a1_m, a2_m, th)
        	return Float64(Float64(cos(th) * a2_m) * Float64(a2_m / sqrt(2.0)))
        end
        
        a2_m = abs(a2);
        a1_m = abs(a1);
        a1_m, a2_m, th = num2cell(sort([a1_m, a2_m, th])){:}
        function tmp = code(a1_m, a2_m, th)
        	tmp = (cos(th) * a2_m) * (a2_m / sqrt(2.0));
        end
        
        a2_m = N[Abs[a2], $MachinePrecision]
        a1_m = N[Abs[a1], $MachinePrecision]
        NOTE: a1_m, a2_m, and th should be sorted in increasing order before calling this function.
        code[a1$95$m_, a2$95$m_, th_] := N[(N[(N[Cos[th], $MachinePrecision] * a2$95$m), $MachinePrecision] * N[(a2$95$m / N[Sqrt[2.0], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]
        
        \begin{array}{l}
        a2_m = \left|a2\right|
        \\
        a1_m = \left|a1\right|
        \\
        [a1_m, a2_m, th] = \mathsf{sort}([a1_m, a2_m, th])\\
        \\
        \left(\cos th \cdot a2\_m\right) \cdot \frac{a2\_m}{\sqrt{2}}
        \end{array}
        
        Derivation
        1. Initial program 99.5%

          \[\frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right) + \frac{\cos th}{\sqrt{2}} \cdot \left(a2 \cdot a2\right) \]
        2. Add Preprocessing
        3. Taylor expanded in a1 around 0

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

            \[\leadsto \frac{\color{blue}{\cos th \cdot {a2}^{2}}}{\sqrt{2}} \]
          2. unpow2N/A

            \[\leadsto \frac{\cos th \cdot \color{blue}{\left(a2 \cdot a2\right)}}{\sqrt{2}} \]
          3. associate-*r*N/A

            \[\leadsto \frac{\color{blue}{\left(\cos th \cdot a2\right) \cdot a2}}{\sqrt{2}} \]
          4. associate-/l*N/A

            \[\leadsto \color{blue}{\left(\cos th \cdot a2\right) \cdot \frac{a2}{\sqrt{2}}} \]
          5. lower-*.f64N/A

            \[\leadsto \color{blue}{\left(\cos th \cdot a2\right) \cdot \frac{a2}{\sqrt{2}}} \]
          6. lower-*.f64N/A

            \[\leadsto \color{blue}{\left(\cos th \cdot a2\right)} \cdot \frac{a2}{\sqrt{2}} \]
          7. lower-cos.f64N/A

            \[\leadsto \left(\color{blue}{\cos th} \cdot a2\right) \cdot \frac{a2}{\sqrt{2}} \]
          8. lower-/.f64N/A

            \[\leadsto \left(\cos th \cdot a2\right) \cdot \color{blue}{\frac{a2}{\sqrt{2}}} \]
          9. lower-sqrt.f6458.9

            \[\leadsto \left(\cos th \cdot a2\right) \cdot \frac{a2}{\color{blue}{\sqrt{2}}} \]
        5. Applied rewrites58.9%

          \[\leadsto \color{blue}{\left(\cos th \cdot a2\right) \cdot \frac{a2}{\sqrt{2}}} \]
        6. Add Preprocessing

        Alternative 6: 99.1% accurate, 2.0× speedup?

        \[\begin{array}{l} a2_m = \left|a2\right| \\ a1_m = \left|a1\right| \\ [a1_m, a2_m, th] = \mathsf{sort}([a1_m, a2_m, th])\\ \\ \left(\left(a2\_m \cdot a2\_m\right) \cdot 0.5\right) \cdot \left(\sqrt{2} \cdot \cos th\right) \end{array} \]
        a2_m = (fabs.f64 a2)
        a1_m = (fabs.f64 a1)
        NOTE: a1_m, a2_m, and th should be sorted in increasing order before calling this function.
        (FPCore (a1_m a2_m th)
         :precision binary64
         (* (* (* a2_m a2_m) 0.5) (* (sqrt 2.0) (cos th))))
        a2_m = fabs(a2);
        a1_m = fabs(a1);
        assert(a1_m < a2_m && a2_m < th);
        double code(double a1_m, double a2_m, double th) {
        	return ((a2_m * a2_m) * 0.5) * (sqrt(2.0) * cos(th));
        }
        
        a2_m =     private
        a1_m =     private
        NOTE: a1_m, a2_m, and th should be sorted in increasing order before calling this function.
        module fmin_fmax_functions
            implicit none
            private
            public fmax
            public fmin
        
            interface fmax
                module procedure fmax88
                module procedure fmax44
                module procedure fmax84
                module procedure fmax48
            end interface
            interface fmin
                module procedure fmin88
                module procedure fmin44
                module procedure fmin84
                module procedure fmin48
            end interface
        contains
            real(8) function fmax88(x, y) result (res)
                real(8), intent (in) :: x
                real(8), intent (in) :: y
                res = merge(y, merge(x, max(x, y), y /= y), x /= x)
            end function
            real(4) function fmax44(x, y) result (res)
                real(4), intent (in) :: x
                real(4), intent (in) :: y
                res = merge(y, merge(x, max(x, y), y /= y), x /= x)
            end function
            real(8) function fmax84(x, y) result(res)
                real(8), intent (in) :: x
                real(4), intent (in) :: y
                res = merge(dble(y), merge(x, max(x, dble(y)), y /= y), x /= x)
            end function
            real(8) function fmax48(x, y) result(res)
                real(4), intent (in) :: x
                real(8), intent (in) :: y
                res = merge(y, merge(dble(x), max(dble(x), y), y /= y), x /= x)
            end function
            real(8) function fmin88(x, y) result (res)
                real(8), intent (in) :: x
                real(8), intent (in) :: y
                res = merge(y, merge(x, min(x, y), y /= y), x /= x)
            end function
            real(4) function fmin44(x, y) result (res)
                real(4), intent (in) :: x
                real(4), intent (in) :: y
                res = merge(y, merge(x, min(x, y), y /= y), x /= x)
            end function
            real(8) function fmin84(x, y) result(res)
                real(8), intent (in) :: x
                real(4), intent (in) :: y
                res = merge(dble(y), merge(x, min(x, dble(y)), y /= y), x /= x)
            end function
            real(8) function fmin48(x, y) result(res)
                real(4), intent (in) :: x
                real(8), intent (in) :: y
                res = merge(y, merge(dble(x), min(dble(x), y), y /= y), x /= x)
            end function
        end module
        
        real(8) function code(a1_m, a2_m, th)
        use fmin_fmax_functions
            real(8), intent (in) :: a1_m
            real(8), intent (in) :: a2_m
            real(8), intent (in) :: th
            code = ((a2_m * a2_m) * 0.5d0) * (sqrt(2.0d0) * cos(th))
        end function
        
        a2_m = Math.abs(a2);
        a1_m = Math.abs(a1);
        assert a1_m < a2_m && a2_m < th;
        public static double code(double a1_m, double a2_m, double th) {
        	return ((a2_m * a2_m) * 0.5) * (Math.sqrt(2.0) * Math.cos(th));
        }
        
        a2_m = math.fabs(a2)
        a1_m = math.fabs(a1)
        [a1_m, a2_m, th] = sort([a1_m, a2_m, th])
        def code(a1_m, a2_m, th):
        	return ((a2_m * a2_m) * 0.5) * (math.sqrt(2.0) * math.cos(th))
        
        a2_m = abs(a2)
        a1_m = abs(a1)
        a1_m, a2_m, th = sort([a1_m, a2_m, th])
        function code(a1_m, a2_m, th)
        	return Float64(Float64(Float64(a2_m * a2_m) * 0.5) * Float64(sqrt(2.0) * cos(th)))
        end
        
        a2_m = abs(a2);
        a1_m = abs(a1);
        a1_m, a2_m, th = num2cell(sort([a1_m, a2_m, th])){:}
        function tmp = code(a1_m, a2_m, th)
        	tmp = ((a2_m * a2_m) * 0.5) * (sqrt(2.0) * cos(th));
        end
        
        a2_m = N[Abs[a2], $MachinePrecision]
        a1_m = N[Abs[a1], $MachinePrecision]
        NOTE: a1_m, a2_m, and th should be sorted in increasing order before calling this function.
        code[a1$95$m_, a2$95$m_, th_] := N[(N[(N[(a2$95$m * a2$95$m), $MachinePrecision] * 0.5), $MachinePrecision] * N[(N[Sqrt[2.0], $MachinePrecision] * N[Cos[th], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]
        
        \begin{array}{l}
        a2_m = \left|a2\right|
        \\
        a1_m = \left|a1\right|
        \\
        [a1_m, a2_m, th] = \mathsf{sort}([a1_m, a2_m, th])\\
        \\
        \left(\left(a2\_m \cdot a2\_m\right) \cdot 0.5\right) \cdot \left(\sqrt{2} \cdot \cos th\right)
        \end{array}
        
        Derivation
        1. Initial program 99.5%

          \[\frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right) + \frac{\cos th}{\sqrt{2}} \cdot \left(a2 \cdot a2\right) \]
        2. Add Preprocessing
        3. Step-by-step derivation
          1. lift-+.f64N/A

            \[\leadsto \color{blue}{\frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right) + \frac{\cos th}{\sqrt{2}} \cdot \left(a2 \cdot a2\right)} \]
          2. +-commutativeN/A

            \[\leadsto \color{blue}{\frac{\cos th}{\sqrt{2}} \cdot \left(a2 \cdot a2\right) + \frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right)} \]
          3. lift-*.f64N/A

            \[\leadsto \color{blue}{\frac{\cos th}{\sqrt{2}} \cdot \left(a2 \cdot a2\right)} + \frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right) \]
          4. lift-/.f64N/A

            \[\leadsto \color{blue}{\frac{\cos th}{\sqrt{2}}} \cdot \left(a2 \cdot a2\right) + \frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right) \]
          5. associate-*l/N/A

            \[\leadsto \color{blue}{\frac{\cos th \cdot \left(a2 \cdot a2\right)}{\sqrt{2}}} + \frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right) \]
          6. associate-/l*N/A

            \[\leadsto \color{blue}{\cos th \cdot \frac{a2 \cdot a2}{\sqrt{2}}} + \frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right) \]
          7. lift-*.f64N/A

            \[\leadsto \cos th \cdot \frac{a2 \cdot a2}{\sqrt{2}} + \color{blue}{\frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right)} \]
          8. lift-/.f64N/A

            \[\leadsto \cos th \cdot \frac{a2 \cdot a2}{\sqrt{2}} + \color{blue}{\frac{\cos th}{\sqrt{2}}} \cdot \left(a1 \cdot a1\right) \]
          9. associate-*l/N/A

            \[\leadsto \cos th \cdot \frac{a2 \cdot a2}{\sqrt{2}} + \color{blue}{\frac{\cos th \cdot \left(a1 \cdot a1\right)}{\sqrt{2}}} \]
          10. associate-/l*N/A

            \[\leadsto \cos th \cdot \frac{a2 \cdot a2}{\sqrt{2}} + \color{blue}{\cos th \cdot \frac{a1 \cdot a1}{\sqrt{2}}} \]
          11. distribute-lft-outN/A

            \[\leadsto \color{blue}{\cos th \cdot \left(\frac{a2 \cdot a2}{\sqrt{2}} + \frac{a1 \cdot a1}{\sqrt{2}}\right)} \]
          12. lower-*.f64N/A

            \[\leadsto \color{blue}{\cos th \cdot \left(\frac{a2 \cdot a2}{\sqrt{2}} + \frac{a1 \cdot a1}{\sqrt{2}}\right)} \]
          13. lift-*.f64N/A

            \[\leadsto \cos th \cdot \left(\frac{\color{blue}{a2 \cdot a2}}{\sqrt{2}} + \frac{a1 \cdot a1}{\sqrt{2}}\right) \]
          14. associate-/l*N/A

            \[\leadsto \cos th \cdot \left(\color{blue}{a2 \cdot \frac{a2}{\sqrt{2}}} + \frac{a1 \cdot a1}{\sqrt{2}}\right) \]
          15. lower-fma.f64N/A

            \[\leadsto \cos th \cdot \color{blue}{\mathsf{fma}\left(a2, \frac{a2}{\sqrt{2}}, \frac{a1 \cdot a1}{\sqrt{2}}\right)} \]
          16. lower-/.f64N/A

            \[\leadsto \cos th \cdot \mathsf{fma}\left(a2, \color{blue}{\frac{a2}{\sqrt{2}}}, \frac{a1 \cdot a1}{\sqrt{2}}\right) \]
          17. lift-*.f64N/A

            \[\leadsto \cos th \cdot \mathsf{fma}\left(a2, \frac{a2}{\sqrt{2}}, \frac{\color{blue}{a1 \cdot a1}}{\sqrt{2}}\right) \]
        4. Applied rewrites99.6%

          \[\leadsto \color{blue}{\cos th \cdot \mathsf{fma}\left(a2, \frac{a2}{\sqrt{2}}, a1 \cdot \frac{a1}{\sqrt{2}}\right)} \]
        5. Step-by-step derivation
          1. lift-fma.f64N/A

            \[\leadsto \cos th \cdot \color{blue}{\left(a2 \cdot \frac{a2}{\sqrt{2}} + a1 \cdot \frac{a1}{\sqrt{2}}\right)} \]
          2. lift-/.f64N/A

            \[\leadsto \cos th \cdot \left(a2 \cdot \color{blue}{\frac{a2}{\sqrt{2}}} + a1 \cdot \frac{a1}{\sqrt{2}}\right) \]
          3. associate-*r/N/A

            \[\leadsto \cos th \cdot \left(\color{blue}{\frac{a2 \cdot a2}{\sqrt{2}}} + a1 \cdot \frac{a1}{\sqrt{2}}\right) \]
          4. lift-*.f64N/A

            \[\leadsto \cos th \cdot \left(\frac{\color{blue}{a2 \cdot a2}}{\sqrt{2}} + a1 \cdot \frac{a1}{\sqrt{2}}\right) \]
          5. lift-*.f64N/A

            \[\leadsto \cos th \cdot \left(\frac{a2 \cdot a2}{\sqrt{2}} + \color{blue}{a1 \cdot \frac{a1}{\sqrt{2}}}\right) \]
          6. lift-/.f64N/A

            \[\leadsto \cos th \cdot \left(\frac{a2 \cdot a2}{\sqrt{2}} + a1 \cdot \color{blue}{\frac{a1}{\sqrt{2}}}\right) \]
          7. associate-*r/N/A

            \[\leadsto \cos th \cdot \left(\frac{a2 \cdot a2}{\sqrt{2}} + \color{blue}{\frac{a1 \cdot a1}{\sqrt{2}}}\right) \]
          8. lift-*.f64N/A

            \[\leadsto \cos th \cdot \left(\frac{a2 \cdot a2}{\sqrt{2}} + \frac{\color{blue}{a1 \cdot a1}}{\sqrt{2}}\right) \]
          9. frac-addN/A

            \[\leadsto \cos th \cdot \color{blue}{\frac{\left(a2 \cdot a2\right) \cdot \sqrt{2} + \sqrt{2} \cdot \left(a1 \cdot a1\right)}{\sqrt{2} \cdot \sqrt{2}}} \]
          10. lift-*.f64N/A

            \[\leadsto \cos th \cdot \frac{\left(a2 \cdot a2\right) \cdot \sqrt{2} + \color{blue}{\sqrt{2} \cdot \left(a1 \cdot a1\right)}}{\sqrt{2} \cdot \sqrt{2}} \]
          11. lift-fma.f64N/A

            \[\leadsto \cos th \cdot \frac{\color{blue}{\mathsf{fma}\left(a2 \cdot a2, \sqrt{2}, \sqrt{2} \cdot \left(a1 \cdot a1\right)\right)}}{\sqrt{2} \cdot \sqrt{2}} \]
          12. lift-sqrt.f64N/A

            \[\leadsto \cos th \cdot \frac{\mathsf{fma}\left(a2 \cdot a2, \sqrt{2}, \sqrt{2} \cdot \left(a1 \cdot a1\right)\right)}{\color{blue}{\sqrt{2}} \cdot \sqrt{2}} \]
          13. lift-sqrt.f64N/A

            \[\leadsto \cos th \cdot \frac{\mathsf{fma}\left(a2 \cdot a2, \sqrt{2}, \sqrt{2} \cdot \left(a1 \cdot a1\right)\right)}{\sqrt{2} \cdot \color{blue}{\sqrt{2}}} \]
          14. rem-square-sqrtN/A

            \[\leadsto \cos th \cdot \frac{\mathsf{fma}\left(a2 \cdot a2, \sqrt{2}, \sqrt{2} \cdot \left(a1 \cdot a1\right)\right)}{\color{blue}{2}} \]
          15. lift-/.f6499.6

            \[\leadsto \cos th \cdot \color{blue}{\frac{\mathsf{fma}\left(a2 \cdot a2, \sqrt{2}, \sqrt{2} \cdot \left(a1 \cdot a1\right)\right)}{2}} \]
        6. Applied rewrites99.6%

          \[\leadsto \cos th \cdot \color{blue}{\frac{\sqrt{2} \cdot \mathsf{fma}\left(a1, a1, a2 \cdot a2\right)}{2}} \]
        7. Taylor expanded in a1 around 0

          \[\leadsto \color{blue}{\frac{1}{2} \cdot \left({a2}^{2} \cdot \left(\cos th \cdot \sqrt{2}\right)\right)} \]
        8. Step-by-step derivation
          1. associate-*r*N/A

            \[\leadsto \color{blue}{\left(\frac{1}{2} \cdot {a2}^{2}\right) \cdot \left(\cos th \cdot \sqrt{2}\right)} \]
          2. lower-*.f64N/A

            \[\leadsto \color{blue}{\left(\frac{1}{2} \cdot {a2}^{2}\right) \cdot \left(\cos th \cdot \sqrt{2}\right)} \]
          3. *-commutativeN/A

            \[\leadsto \color{blue}{\left({a2}^{2} \cdot \frac{1}{2}\right)} \cdot \left(\cos th \cdot \sqrt{2}\right) \]
          4. lower-*.f64N/A

            \[\leadsto \color{blue}{\left({a2}^{2} \cdot \frac{1}{2}\right)} \cdot \left(\cos th \cdot \sqrt{2}\right) \]
          5. unpow2N/A

            \[\leadsto \left(\color{blue}{\left(a2 \cdot a2\right)} \cdot \frac{1}{2}\right) \cdot \left(\cos th \cdot \sqrt{2}\right) \]
          6. lower-*.f64N/A

            \[\leadsto \left(\color{blue}{\left(a2 \cdot a2\right)} \cdot \frac{1}{2}\right) \cdot \left(\cos th \cdot \sqrt{2}\right) \]
          7. *-commutativeN/A

            \[\leadsto \left(\left(a2 \cdot a2\right) \cdot \frac{1}{2}\right) \cdot \color{blue}{\left(\sqrt{2} \cdot \cos th\right)} \]
          8. lower-*.f64N/A

            \[\leadsto \left(\left(a2 \cdot a2\right) \cdot \frac{1}{2}\right) \cdot \color{blue}{\left(\sqrt{2} \cdot \cos th\right)} \]
          9. lower-sqrt.f64N/A

            \[\leadsto \left(\left(a2 \cdot a2\right) \cdot \frac{1}{2}\right) \cdot \left(\color{blue}{\sqrt{2}} \cdot \cos th\right) \]
          10. lower-cos.f6458.9

            \[\leadsto \left(\left(a2 \cdot a2\right) \cdot 0.5\right) \cdot \left(\sqrt{2} \cdot \color{blue}{\cos th}\right) \]
        9. Applied rewrites58.9%

          \[\leadsto \color{blue}{\left(\left(a2 \cdot a2\right) \cdot 0.5\right) \cdot \left(\sqrt{2} \cdot \cos th\right)} \]
        10. Add Preprocessing

        Alternative 7: 66.6% accurate, 8.3× speedup?

        \[\begin{array}{l} a2_m = \left|a2\right| \\ a1_m = \left|a1\right| \\ [a1_m, a2_m, th] = \mathsf{sort}([a1_m, a2_m, th])\\ \\ \left(0.5 \cdot \sqrt{2}\right) \cdot \mathsf{fma}\left(a1\_m, a1\_m, a2\_m \cdot a2\_m\right) \end{array} \]
        a2_m = (fabs.f64 a2)
        a1_m = (fabs.f64 a1)
        NOTE: a1_m, a2_m, and th should be sorted in increasing order before calling this function.
        (FPCore (a1_m a2_m th)
         :precision binary64
         (* (* 0.5 (sqrt 2.0)) (fma a1_m a1_m (* a2_m a2_m))))
        a2_m = fabs(a2);
        a1_m = fabs(a1);
        assert(a1_m < a2_m && a2_m < th);
        double code(double a1_m, double a2_m, double th) {
        	return (0.5 * sqrt(2.0)) * fma(a1_m, a1_m, (a2_m * a2_m));
        }
        
        a2_m = abs(a2)
        a1_m = abs(a1)
        a1_m, a2_m, th = sort([a1_m, a2_m, th])
        function code(a1_m, a2_m, th)
        	return Float64(Float64(0.5 * sqrt(2.0)) * fma(a1_m, a1_m, Float64(a2_m * a2_m)))
        end
        
        a2_m = N[Abs[a2], $MachinePrecision]
        a1_m = N[Abs[a1], $MachinePrecision]
        NOTE: a1_m, a2_m, and th should be sorted in increasing order before calling this function.
        code[a1$95$m_, a2$95$m_, th_] := N[(N[(0.5 * N[Sqrt[2.0], $MachinePrecision]), $MachinePrecision] * N[(a1$95$m * a1$95$m + N[(a2$95$m * a2$95$m), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]
        
        \begin{array}{l}
        a2_m = \left|a2\right|
        \\
        a1_m = \left|a1\right|
        \\
        [a1_m, a2_m, th] = \mathsf{sort}([a1_m, a2_m, th])\\
        \\
        \left(0.5 \cdot \sqrt{2}\right) \cdot \mathsf{fma}\left(a1\_m, a1\_m, a2\_m \cdot a2\_m\right)
        \end{array}
        
        Derivation
        1. Initial program 99.5%

          \[\frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right) + \frac{\cos th}{\sqrt{2}} \cdot \left(a2 \cdot a2\right) \]
        2. Add Preprocessing
        3. Taylor expanded in th around 0

          \[\leadsto \color{blue}{\frac{{a1}^{2}}{\sqrt{2}} + \frac{{a2}^{2}}{\sqrt{2}}} \]
        4. Step-by-step derivation
          1. div-add-revN/A

            \[\leadsto \color{blue}{\frac{{a1}^{2} + {a2}^{2}}{\sqrt{2}}} \]
          2. lower-/.f64N/A

            \[\leadsto \color{blue}{\frac{{a1}^{2} + {a2}^{2}}{\sqrt{2}}} \]
          3. +-commutativeN/A

            \[\leadsto \frac{\color{blue}{{a2}^{2} + {a1}^{2}}}{\sqrt{2}} \]
          4. unpow2N/A

            \[\leadsto \frac{\color{blue}{a2 \cdot a2} + {a1}^{2}}{\sqrt{2}} \]
          5. lower-fma.f64N/A

            \[\leadsto \frac{\color{blue}{\mathsf{fma}\left(a2, a2, {a1}^{2}\right)}}{\sqrt{2}} \]
          6. unpow2N/A

            \[\leadsto \frac{\mathsf{fma}\left(a2, a2, \color{blue}{a1 \cdot a1}\right)}{\sqrt{2}} \]
          7. lower-*.f64N/A

            \[\leadsto \frac{\mathsf{fma}\left(a2, a2, \color{blue}{a1 \cdot a1}\right)}{\sqrt{2}} \]
          8. lower-sqrt.f6465.4

            \[\leadsto \frac{\mathsf{fma}\left(a2, a2, a1 \cdot a1\right)}{\color{blue}{\sqrt{2}}} \]
        5. Applied rewrites65.4%

          \[\leadsto \color{blue}{\frac{\mathsf{fma}\left(a2, a2, a1 \cdot a1\right)}{\sqrt{2}}} \]
        6. Step-by-step derivation
          1. Applied rewrites65.5%

            \[\leadsto \frac{\mathsf{fma}\left(a2 \cdot a2, \sqrt{2}, \sqrt{2} \cdot \left(a1 \cdot a1\right)\right)}{\color{blue}{2}} \]
          2. Taylor expanded in a1 around 0

            \[\leadsto \frac{1}{2} \cdot \left({a1}^{2} \cdot \sqrt{2}\right) + \color{blue}{\frac{1}{2} \cdot \left({a2}^{2} \cdot \sqrt{2}\right)} \]
          3. Step-by-step derivation
            1. Applied rewrites65.5%

              \[\leadsto \left(0.5 \cdot \sqrt{2}\right) \cdot \color{blue}{\mathsf{fma}\left(a1, a1, a2 \cdot a2\right)} \]
            2. Add Preprocessing

            Alternative 8: 66.4% accurate, 9.9× speedup?

            \[\begin{array}{l} a2_m = \left|a2\right| \\ a1_m = \left|a1\right| \\ [a1_m, a2_m, th] = \mathsf{sort}([a1_m, a2_m, th])\\ \\ a2\_m \cdot \frac{a2\_m}{\sqrt{2}} \end{array} \]
            a2_m = (fabs.f64 a2)
            a1_m = (fabs.f64 a1)
            NOTE: a1_m, a2_m, and th should be sorted in increasing order before calling this function.
            (FPCore (a1_m a2_m th) :precision binary64 (* a2_m (/ a2_m (sqrt 2.0))))
            a2_m = fabs(a2);
            a1_m = fabs(a1);
            assert(a1_m < a2_m && a2_m < th);
            double code(double a1_m, double a2_m, double th) {
            	return a2_m * (a2_m / sqrt(2.0));
            }
            
            a2_m =     private
            a1_m =     private
            NOTE: a1_m, a2_m, and th should be sorted in increasing order before calling this function.
            module fmin_fmax_functions
                implicit none
                private
                public fmax
                public fmin
            
                interface fmax
                    module procedure fmax88
                    module procedure fmax44
                    module procedure fmax84
                    module procedure fmax48
                end interface
                interface fmin
                    module procedure fmin88
                    module procedure fmin44
                    module procedure fmin84
                    module procedure fmin48
                end interface
            contains
                real(8) function fmax88(x, y) result (res)
                    real(8), intent (in) :: x
                    real(8), intent (in) :: y
                    res = merge(y, merge(x, max(x, y), y /= y), x /= x)
                end function
                real(4) function fmax44(x, y) result (res)
                    real(4), intent (in) :: x
                    real(4), intent (in) :: y
                    res = merge(y, merge(x, max(x, y), y /= y), x /= x)
                end function
                real(8) function fmax84(x, y) result(res)
                    real(8), intent (in) :: x
                    real(4), intent (in) :: y
                    res = merge(dble(y), merge(x, max(x, dble(y)), y /= y), x /= x)
                end function
                real(8) function fmax48(x, y) result(res)
                    real(4), intent (in) :: x
                    real(8), intent (in) :: y
                    res = merge(y, merge(dble(x), max(dble(x), y), y /= y), x /= x)
                end function
                real(8) function fmin88(x, y) result (res)
                    real(8), intent (in) :: x
                    real(8), intent (in) :: y
                    res = merge(y, merge(x, min(x, y), y /= y), x /= x)
                end function
                real(4) function fmin44(x, y) result (res)
                    real(4), intent (in) :: x
                    real(4), intent (in) :: y
                    res = merge(y, merge(x, min(x, y), y /= y), x /= x)
                end function
                real(8) function fmin84(x, y) result(res)
                    real(8), intent (in) :: x
                    real(4), intent (in) :: y
                    res = merge(dble(y), merge(x, min(x, dble(y)), y /= y), x /= x)
                end function
                real(8) function fmin48(x, y) result(res)
                    real(4), intent (in) :: x
                    real(8), intent (in) :: y
                    res = merge(y, merge(dble(x), min(dble(x), y), y /= y), x /= x)
                end function
            end module
            
            real(8) function code(a1_m, a2_m, th)
            use fmin_fmax_functions
                real(8), intent (in) :: a1_m
                real(8), intent (in) :: a2_m
                real(8), intent (in) :: th
                code = a2_m * (a2_m / sqrt(2.0d0))
            end function
            
            a2_m = Math.abs(a2);
            a1_m = Math.abs(a1);
            assert a1_m < a2_m && a2_m < th;
            public static double code(double a1_m, double a2_m, double th) {
            	return a2_m * (a2_m / Math.sqrt(2.0));
            }
            
            a2_m = math.fabs(a2)
            a1_m = math.fabs(a1)
            [a1_m, a2_m, th] = sort([a1_m, a2_m, th])
            def code(a1_m, a2_m, th):
            	return a2_m * (a2_m / math.sqrt(2.0))
            
            a2_m = abs(a2)
            a1_m = abs(a1)
            a1_m, a2_m, th = sort([a1_m, a2_m, th])
            function code(a1_m, a2_m, th)
            	return Float64(a2_m * Float64(a2_m / sqrt(2.0)))
            end
            
            a2_m = abs(a2);
            a1_m = abs(a1);
            a1_m, a2_m, th = num2cell(sort([a1_m, a2_m, th])){:}
            function tmp = code(a1_m, a2_m, th)
            	tmp = a2_m * (a2_m / sqrt(2.0));
            end
            
            a2_m = N[Abs[a2], $MachinePrecision]
            a1_m = N[Abs[a1], $MachinePrecision]
            NOTE: a1_m, a2_m, and th should be sorted in increasing order before calling this function.
            code[a1$95$m_, a2$95$m_, th_] := N[(a2$95$m * N[(a2$95$m / N[Sqrt[2.0], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]
            
            \begin{array}{l}
            a2_m = \left|a2\right|
            \\
            a1_m = \left|a1\right|
            \\
            [a1_m, a2_m, th] = \mathsf{sort}([a1_m, a2_m, th])\\
            \\
            a2\_m \cdot \frac{a2\_m}{\sqrt{2}}
            \end{array}
            
            Derivation
            1. Initial program 99.5%

              \[\frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right) + \frac{\cos th}{\sqrt{2}} \cdot \left(a2 \cdot a2\right) \]
            2. Add Preprocessing
            3. Taylor expanded in th around 0

              \[\leadsto \color{blue}{\frac{{a1}^{2}}{\sqrt{2}} + \frac{{a2}^{2}}{\sqrt{2}}} \]
            4. Step-by-step derivation
              1. div-add-revN/A

                \[\leadsto \color{blue}{\frac{{a1}^{2} + {a2}^{2}}{\sqrt{2}}} \]
              2. lower-/.f64N/A

                \[\leadsto \color{blue}{\frac{{a1}^{2} + {a2}^{2}}{\sqrt{2}}} \]
              3. +-commutativeN/A

                \[\leadsto \frac{\color{blue}{{a2}^{2} + {a1}^{2}}}{\sqrt{2}} \]
              4. unpow2N/A

                \[\leadsto \frac{\color{blue}{a2 \cdot a2} + {a1}^{2}}{\sqrt{2}} \]
              5. lower-fma.f64N/A

                \[\leadsto \frac{\color{blue}{\mathsf{fma}\left(a2, a2, {a1}^{2}\right)}}{\sqrt{2}} \]
              6. unpow2N/A

                \[\leadsto \frac{\mathsf{fma}\left(a2, a2, \color{blue}{a1 \cdot a1}\right)}{\sqrt{2}} \]
              7. lower-*.f64N/A

                \[\leadsto \frac{\mathsf{fma}\left(a2, a2, \color{blue}{a1 \cdot a1}\right)}{\sqrt{2}} \]
              8. lower-sqrt.f6465.4

                \[\leadsto \frac{\mathsf{fma}\left(a2, a2, a1 \cdot a1\right)}{\color{blue}{\sqrt{2}}} \]
            5. Applied rewrites65.4%

              \[\leadsto \color{blue}{\frac{\mathsf{fma}\left(a2, a2, a1 \cdot a1\right)}{\sqrt{2}}} \]
            6. Taylor expanded in a1 around 0

              \[\leadsto \frac{{a2}^{2}}{\color{blue}{\sqrt{2}}} \]
            7. Step-by-step derivation
              1. Applied rewrites40.4%

                \[\leadsto a2 \cdot \color{blue}{\frac{a2}{\sqrt{2}}} \]
              2. Add Preprocessing

              Alternative 9: 66.4% accurate, 10.2× speedup?

              \[\begin{array}{l} a2_m = \left|a2\right| \\ a1_m = \left|a1\right| \\ [a1_m, a2_m, th] = \mathsf{sort}([a1_m, a2_m, th])\\ \\ \left(\left(0.5 \cdot \sqrt{2}\right) \cdot a2\_m\right) \cdot a2\_m \end{array} \]
              a2_m = (fabs.f64 a2)
              a1_m = (fabs.f64 a1)
              NOTE: a1_m, a2_m, and th should be sorted in increasing order before calling this function.
              (FPCore (a1_m a2_m th)
               :precision binary64
               (* (* (* 0.5 (sqrt 2.0)) a2_m) a2_m))
              a2_m = fabs(a2);
              a1_m = fabs(a1);
              assert(a1_m < a2_m && a2_m < th);
              double code(double a1_m, double a2_m, double th) {
              	return ((0.5 * sqrt(2.0)) * a2_m) * a2_m;
              }
              
              a2_m =     private
              a1_m =     private
              NOTE: a1_m, a2_m, and th should be sorted in increasing order before calling this function.
              module fmin_fmax_functions
                  implicit none
                  private
                  public fmax
                  public fmin
              
                  interface fmax
                      module procedure fmax88
                      module procedure fmax44
                      module procedure fmax84
                      module procedure fmax48
                  end interface
                  interface fmin
                      module procedure fmin88
                      module procedure fmin44
                      module procedure fmin84
                      module procedure fmin48
                  end interface
              contains
                  real(8) function fmax88(x, y) result (res)
                      real(8), intent (in) :: x
                      real(8), intent (in) :: y
                      res = merge(y, merge(x, max(x, y), y /= y), x /= x)
                  end function
                  real(4) function fmax44(x, y) result (res)
                      real(4), intent (in) :: x
                      real(4), intent (in) :: y
                      res = merge(y, merge(x, max(x, y), y /= y), x /= x)
                  end function
                  real(8) function fmax84(x, y) result(res)
                      real(8), intent (in) :: x
                      real(4), intent (in) :: y
                      res = merge(dble(y), merge(x, max(x, dble(y)), y /= y), x /= x)
                  end function
                  real(8) function fmax48(x, y) result(res)
                      real(4), intent (in) :: x
                      real(8), intent (in) :: y
                      res = merge(y, merge(dble(x), max(dble(x), y), y /= y), x /= x)
                  end function
                  real(8) function fmin88(x, y) result (res)
                      real(8), intent (in) :: x
                      real(8), intent (in) :: y
                      res = merge(y, merge(x, min(x, y), y /= y), x /= x)
                  end function
                  real(4) function fmin44(x, y) result (res)
                      real(4), intent (in) :: x
                      real(4), intent (in) :: y
                      res = merge(y, merge(x, min(x, y), y /= y), x /= x)
                  end function
                  real(8) function fmin84(x, y) result(res)
                      real(8), intent (in) :: x
                      real(4), intent (in) :: y
                      res = merge(dble(y), merge(x, min(x, dble(y)), y /= y), x /= x)
                  end function
                  real(8) function fmin48(x, y) result(res)
                      real(4), intent (in) :: x
                      real(8), intent (in) :: y
                      res = merge(y, merge(dble(x), min(dble(x), y), y /= y), x /= x)
                  end function
              end module
              
              real(8) function code(a1_m, a2_m, th)
              use fmin_fmax_functions
                  real(8), intent (in) :: a1_m
                  real(8), intent (in) :: a2_m
                  real(8), intent (in) :: th
                  code = ((0.5d0 * sqrt(2.0d0)) * a2_m) * a2_m
              end function
              
              a2_m = Math.abs(a2);
              a1_m = Math.abs(a1);
              assert a1_m < a2_m && a2_m < th;
              public static double code(double a1_m, double a2_m, double th) {
              	return ((0.5 * Math.sqrt(2.0)) * a2_m) * a2_m;
              }
              
              a2_m = math.fabs(a2)
              a1_m = math.fabs(a1)
              [a1_m, a2_m, th] = sort([a1_m, a2_m, th])
              def code(a1_m, a2_m, th):
              	return ((0.5 * math.sqrt(2.0)) * a2_m) * a2_m
              
              a2_m = abs(a2)
              a1_m = abs(a1)
              a1_m, a2_m, th = sort([a1_m, a2_m, th])
              function code(a1_m, a2_m, th)
              	return Float64(Float64(Float64(0.5 * sqrt(2.0)) * a2_m) * a2_m)
              end
              
              a2_m = abs(a2);
              a1_m = abs(a1);
              a1_m, a2_m, th = num2cell(sort([a1_m, a2_m, th])){:}
              function tmp = code(a1_m, a2_m, th)
              	tmp = ((0.5 * sqrt(2.0)) * a2_m) * a2_m;
              end
              
              a2_m = N[Abs[a2], $MachinePrecision]
              a1_m = N[Abs[a1], $MachinePrecision]
              NOTE: a1_m, a2_m, and th should be sorted in increasing order before calling this function.
              code[a1$95$m_, a2$95$m_, th_] := N[(N[(N[(0.5 * N[Sqrt[2.0], $MachinePrecision]), $MachinePrecision] * a2$95$m), $MachinePrecision] * a2$95$m), $MachinePrecision]
              
              \begin{array}{l}
              a2_m = \left|a2\right|
              \\
              a1_m = \left|a1\right|
              \\
              [a1_m, a2_m, th] = \mathsf{sort}([a1_m, a2_m, th])\\
              \\
              \left(\left(0.5 \cdot \sqrt{2}\right) \cdot a2\_m\right) \cdot a2\_m
              \end{array}
              
              Derivation
              1. Initial program 99.5%

                \[\frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right) + \frac{\cos th}{\sqrt{2}} \cdot \left(a2 \cdot a2\right) \]
              2. Add Preprocessing
              3. Taylor expanded in th around 0

                \[\leadsto \color{blue}{\frac{{a1}^{2}}{\sqrt{2}} + \frac{{a2}^{2}}{\sqrt{2}}} \]
              4. Step-by-step derivation
                1. div-add-revN/A

                  \[\leadsto \color{blue}{\frac{{a1}^{2} + {a2}^{2}}{\sqrt{2}}} \]
                2. lower-/.f64N/A

                  \[\leadsto \color{blue}{\frac{{a1}^{2} + {a2}^{2}}{\sqrt{2}}} \]
                3. +-commutativeN/A

                  \[\leadsto \frac{\color{blue}{{a2}^{2} + {a1}^{2}}}{\sqrt{2}} \]
                4. unpow2N/A

                  \[\leadsto \frac{\color{blue}{a2 \cdot a2} + {a1}^{2}}{\sqrt{2}} \]
                5. lower-fma.f64N/A

                  \[\leadsto \frac{\color{blue}{\mathsf{fma}\left(a2, a2, {a1}^{2}\right)}}{\sqrt{2}} \]
                6. unpow2N/A

                  \[\leadsto \frac{\mathsf{fma}\left(a2, a2, \color{blue}{a1 \cdot a1}\right)}{\sqrt{2}} \]
                7. lower-*.f64N/A

                  \[\leadsto \frac{\mathsf{fma}\left(a2, a2, \color{blue}{a1 \cdot a1}\right)}{\sqrt{2}} \]
                8. lower-sqrt.f6465.4

                  \[\leadsto \frac{\mathsf{fma}\left(a2, a2, a1 \cdot a1\right)}{\color{blue}{\sqrt{2}}} \]
              5. Applied rewrites65.4%

                \[\leadsto \color{blue}{\frac{\mathsf{fma}\left(a2, a2, a1 \cdot a1\right)}{\sqrt{2}}} \]
              6. Step-by-step derivation
                1. Applied rewrites65.5%

                  \[\leadsto \frac{\mathsf{fma}\left(a2 \cdot a2, \sqrt{2}, \sqrt{2} \cdot \left(a1 \cdot a1\right)\right)}{\color{blue}{2}} \]
                2. Taylor expanded in a1 around 0

                  \[\leadsto \frac{1}{2} \cdot \color{blue}{\left({a2}^{2} \cdot \sqrt{2}\right)} \]
                3. Step-by-step derivation
                  1. Applied rewrites40.4%

                    \[\leadsto \left(\left(0.5 \cdot \sqrt{2}\right) \cdot a2\right) \cdot \color{blue}{a2} \]
                  2. Add Preprocessing

                  Alternative 10: 14.0% accurate, 10.2× speedup?

                  \[\begin{array}{l} a2_m = \left|a2\right| \\ a1_m = \left|a1\right| \\ [a1_m, a2_m, th] = \mathsf{sort}([a1_m, a2_m, th])\\ \\ \left(a1\_m \cdot a1\_m\right) \cdot \left(0.5 \cdot \sqrt{2}\right) \end{array} \]
                  a2_m = (fabs.f64 a2)
                  a1_m = (fabs.f64 a1)
                  NOTE: a1_m, a2_m, and th should be sorted in increasing order before calling this function.
                  (FPCore (a1_m a2_m th)
                   :precision binary64
                   (* (* a1_m a1_m) (* 0.5 (sqrt 2.0))))
                  a2_m = fabs(a2);
                  a1_m = fabs(a1);
                  assert(a1_m < a2_m && a2_m < th);
                  double code(double a1_m, double a2_m, double th) {
                  	return (a1_m * a1_m) * (0.5 * sqrt(2.0));
                  }
                  
                  a2_m =     private
                  a1_m =     private
                  NOTE: a1_m, a2_m, and th should be sorted in increasing order before calling this function.
                  module fmin_fmax_functions
                      implicit none
                      private
                      public fmax
                      public fmin
                  
                      interface fmax
                          module procedure fmax88
                          module procedure fmax44
                          module procedure fmax84
                          module procedure fmax48
                      end interface
                      interface fmin
                          module procedure fmin88
                          module procedure fmin44
                          module procedure fmin84
                          module procedure fmin48
                      end interface
                  contains
                      real(8) function fmax88(x, y) result (res)
                          real(8), intent (in) :: x
                          real(8), intent (in) :: y
                          res = merge(y, merge(x, max(x, y), y /= y), x /= x)
                      end function
                      real(4) function fmax44(x, y) result (res)
                          real(4), intent (in) :: x
                          real(4), intent (in) :: y
                          res = merge(y, merge(x, max(x, y), y /= y), x /= x)
                      end function
                      real(8) function fmax84(x, y) result(res)
                          real(8), intent (in) :: x
                          real(4), intent (in) :: y
                          res = merge(dble(y), merge(x, max(x, dble(y)), y /= y), x /= x)
                      end function
                      real(8) function fmax48(x, y) result(res)
                          real(4), intent (in) :: x
                          real(8), intent (in) :: y
                          res = merge(y, merge(dble(x), max(dble(x), y), y /= y), x /= x)
                      end function
                      real(8) function fmin88(x, y) result (res)
                          real(8), intent (in) :: x
                          real(8), intent (in) :: y
                          res = merge(y, merge(x, min(x, y), y /= y), x /= x)
                      end function
                      real(4) function fmin44(x, y) result (res)
                          real(4), intent (in) :: x
                          real(4), intent (in) :: y
                          res = merge(y, merge(x, min(x, y), y /= y), x /= x)
                      end function
                      real(8) function fmin84(x, y) result(res)
                          real(8), intent (in) :: x
                          real(4), intent (in) :: y
                          res = merge(dble(y), merge(x, min(x, dble(y)), y /= y), x /= x)
                      end function
                      real(8) function fmin48(x, y) result(res)
                          real(4), intent (in) :: x
                          real(8), intent (in) :: y
                          res = merge(y, merge(dble(x), min(dble(x), y), y /= y), x /= x)
                      end function
                  end module
                  
                  real(8) function code(a1_m, a2_m, th)
                  use fmin_fmax_functions
                      real(8), intent (in) :: a1_m
                      real(8), intent (in) :: a2_m
                      real(8), intent (in) :: th
                      code = (a1_m * a1_m) * (0.5d0 * sqrt(2.0d0))
                  end function
                  
                  a2_m = Math.abs(a2);
                  a1_m = Math.abs(a1);
                  assert a1_m < a2_m && a2_m < th;
                  public static double code(double a1_m, double a2_m, double th) {
                  	return (a1_m * a1_m) * (0.5 * Math.sqrt(2.0));
                  }
                  
                  a2_m = math.fabs(a2)
                  a1_m = math.fabs(a1)
                  [a1_m, a2_m, th] = sort([a1_m, a2_m, th])
                  def code(a1_m, a2_m, th):
                  	return (a1_m * a1_m) * (0.5 * math.sqrt(2.0))
                  
                  a2_m = abs(a2)
                  a1_m = abs(a1)
                  a1_m, a2_m, th = sort([a1_m, a2_m, th])
                  function code(a1_m, a2_m, th)
                  	return Float64(Float64(a1_m * a1_m) * Float64(0.5 * sqrt(2.0)))
                  end
                  
                  a2_m = abs(a2);
                  a1_m = abs(a1);
                  a1_m, a2_m, th = num2cell(sort([a1_m, a2_m, th])){:}
                  function tmp = code(a1_m, a2_m, th)
                  	tmp = (a1_m * a1_m) * (0.5 * sqrt(2.0));
                  end
                  
                  a2_m = N[Abs[a2], $MachinePrecision]
                  a1_m = N[Abs[a1], $MachinePrecision]
                  NOTE: a1_m, a2_m, and th should be sorted in increasing order before calling this function.
                  code[a1$95$m_, a2$95$m_, th_] := N[(N[(a1$95$m * a1$95$m), $MachinePrecision] * N[(0.5 * N[Sqrt[2.0], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]
                  
                  \begin{array}{l}
                  a2_m = \left|a2\right|
                  \\
                  a1_m = \left|a1\right|
                  \\
                  [a1_m, a2_m, th] = \mathsf{sort}([a1_m, a2_m, th])\\
                  \\
                  \left(a1\_m \cdot a1\_m\right) \cdot \left(0.5 \cdot \sqrt{2}\right)
                  \end{array}
                  
                  Derivation
                  1. Initial program 99.5%

                    \[\frac{\cos th}{\sqrt{2}} \cdot \left(a1 \cdot a1\right) + \frac{\cos th}{\sqrt{2}} \cdot \left(a2 \cdot a2\right) \]
                  2. Add Preprocessing
                  3. Taylor expanded in th around 0

                    \[\leadsto \color{blue}{\frac{{a1}^{2}}{\sqrt{2}} + \frac{{a2}^{2}}{\sqrt{2}}} \]
                  4. Step-by-step derivation
                    1. div-add-revN/A

                      \[\leadsto \color{blue}{\frac{{a1}^{2} + {a2}^{2}}{\sqrt{2}}} \]
                    2. lower-/.f64N/A

                      \[\leadsto \color{blue}{\frac{{a1}^{2} + {a2}^{2}}{\sqrt{2}}} \]
                    3. +-commutativeN/A

                      \[\leadsto \frac{\color{blue}{{a2}^{2} + {a1}^{2}}}{\sqrt{2}} \]
                    4. unpow2N/A

                      \[\leadsto \frac{\color{blue}{a2 \cdot a2} + {a1}^{2}}{\sqrt{2}} \]
                    5. lower-fma.f64N/A

                      \[\leadsto \frac{\color{blue}{\mathsf{fma}\left(a2, a2, {a1}^{2}\right)}}{\sqrt{2}} \]
                    6. unpow2N/A

                      \[\leadsto \frac{\mathsf{fma}\left(a2, a2, \color{blue}{a1 \cdot a1}\right)}{\sqrt{2}} \]
                    7. lower-*.f64N/A

                      \[\leadsto \frac{\mathsf{fma}\left(a2, a2, \color{blue}{a1 \cdot a1}\right)}{\sqrt{2}} \]
                    8. lower-sqrt.f6465.4

                      \[\leadsto \frac{\mathsf{fma}\left(a2, a2, a1 \cdot a1\right)}{\color{blue}{\sqrt{2}}} \]
                  5. Applied rewrites65.4%

                    \[\leadsto \color{blue}{\frac{\mathsf{fma}\left(a2, a2, a1 \cdot a1\right)}{\sqrt{2}}} \]
                  6. Step-by-step derivation
                    1. Applied rewrites65.5%

                      \[\leadsto \frac{\mathsf{fma}\left(a2 \cdot a2, \sqrt{2}, \sqrt{2} \cdot \left(a1 \cdot a1\right)\right)}{\color{blue}{2}} \]
                    2. Taylor expanded in a1 around inf

                      \[\leadsto \frac{1}{2} \cdot \color{blue}{\left({a1}^{2} \cdot \sqrt{2}\right)} \]
                    3. Step-by-step derivation
                      1. Applied rewrites36.5%

                        \[\leadsto \left(\left(0.5 \cdot \sqrt{2}\right) \cdot a1\right) \cdot \color{blue}{a1} \]
                      2. Step-by-step derivation
                        1. Applied rewrites36.5%

                          \[\leadsto \left(a1 \cdot a1\right) \cdot \left(0.5 \cdot \color{blue}{\sqrt{2}}\right) \]
                        2. Add Preprocessing

                        Reproduce

                        ?
                        herbie shell --seed 2024350 
                        (FPCore (a1 a2 th)
                          :name "Migdal et al, Equation (64)"
                          :precision binary64
                          (+ (* (/ (cos th) (sqrt 2.0)) (* a1 a1)) (* (/ (cos th) (sqrt 2.0)) (* a2 a2))))