Toniolo and Linder, Equation (3a)

Percentage Accurate: 98.4% → 99.4%
Time: 11.5s
Alternatives: 7
Speedup: 1.0×

Specification

?
\[\begin{array}{l} \\ \sqrt{\frac{1}{2} \cdot \left(1 + \frac{1}{\sqrt{1 + {\left(\frac{2 \cdot \ell}{Om}\right)}^{2} \cdot \left({\sin kx}^{2} + {\sin ky}^{2}\right)}}\right)} \end{array} \]
(FPCore (l Om kx ky)
 :precision binary64
 (sqrt
  (*
   (/ 1.0 2.0)
   (+
    1.0
    (/
     1.0
     (sqrt
      (+
       1.0
       (*
        (pow (/ (* 2.0 l) Om) 2.0)
        (+ (pow (sin kx) 2.0) (pow (sin ky) 2.0))))))))))
double code(double l, double Om, double kx, double ky) {
	return sqrt(((1.0 / 2.0) * (1.0 + (1.0 / sqrt((1.0 + (pow(((2.0 * l) / Om), 2.0) * (pow(sin(kx), 2.0) + pow(sin(ky), 2.0)))))))));
}
real(8) function code(l, om, kx, ky)
    real(8), intent (in) :: l
    real(8), intent (in) :: om
    real(8), intent (in) :: kx
    real(8), intent (in) :: ky
    code = sqrt(((1.0d0 / 2.0d0) * (1.0d0 + (1.0d0 / sqrt((1.0d0 + ((((2.0d0 * l) / om) ** 2.0d0) * ((sin(kx) ** 2.0d0) + (sin(ky) ** 2.0d0)))))))))
end function
public static double code(double l, double Om, double kx, double ky) {
	return Math.sqrt(((1.0 / 2.0) * (1.0 + (1.0 / Math.sqrt((1.0 + (Math.pow(((2.0 * l) / Om), 2.0) * (Math.pow(Math.sin(kx), 2.0) + Math.pow(Math.sin(ky), 2.0)))))))));
}
def code(l, Om, kx, ky):
	return math.sqrt(((1.0 / 2.0) * (1.0 + (1.0 / math.sqrt((1.0 + (math.pow(((2.0 * l) / Om), 2.0) * (math.pow(math.sin(kx), 2.0) + math.pow(math.sin(ky), 2.0)))))))))
function code(l, Om, kx, ky)
	return sqrt(Float64(Float64(1.0 / 2.0) * Float64(1.0 + Float64(1.0 / sqrt(Float64(1.0 + Float64((Float64(Float64(2.0 * l) / Om) ^ 2.0) * Float64((sin(kx) ^ 2.0) + (sin(ky) ^ 2.0)))))))))
end
function tmp = code(l, Om, kx, ky)
	tmp = sqrt(((1.0 / 2.0) * (1.0 + (1.0 / sqrt((1.0 + ((((2.0 * l) / Om) ^ 2.0) * ((sin(kx) ^ 2.0) + (sin(ky) ^ 2.0)))))))));
end
code[l_, Om_, kx_, ky_] := N[Sqrt[N[(N[(1.0 / 2.0), $MachinePrecision] * N[(1.0 + N[(1.0 / N[Sqrt[N[(1.0 + N[(N[Power[N[(N[(2.0 * l), $MachinePrecision] / Om), $MachinePrecision], 2.0], $MachinePrecision] * N[(N[Power[N[Sin[kx], $MachinePrecision], 2.0], $MachinePrecision] + N[Power[N[Sin[ky], $MachinePrecision], 2.0], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]
\begin{array}{l}

\\
\sqrt{\frac{1}{2} \cdot \left(1 + \frac{1}{\sqrt{1 + {\left(\frac{2 \cdot \ell}{Om}\right)}^{2} \cdot \left({\sin kx}^{2} + {\sin ky}^{2}\right)}}\right)}
\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 7 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: 98.4% accurate, 1.0× speedup?

\[\begin{array}{l} \\ \sqrt{\frac{1}{2} \cdot \left(1 + \frac{1}{\sqrt{1 + {\left(\frac{2 \cdot \ell}{Om}\right)}^{2} \cdot \left({\sin kx}^{2} + {\sin ky}^{2}\right)}}\right)} \end{array} \]
(FPCore (l Om kx ky)
 :precision binary64
 (sqrt
  (*
   (/ 1.0 2.0)
   (+
    1.0
    (/
     1.0
     (sqrt
      (+
       1.0
       (*
        (pow (/ (* 2.0 l) Om) 2.0)
        (+ (pow (sin kx) 2.0) (pow (sin ky) 2.0))))))))))
double code(double l, double Om, double kx, double ky) {
	return sqrt(((1.0 / 2.0) * (1.0 + (1.0 / sqrt((1.0 + (pow(((2.0 * l) / Om), 2.0) * (pow(sin(kx), 2.0) + pow(sin(ky), 2.0)))))))));
}
real(8) function code(l, om, kx, ky)
    real(8), intent (in) :: l
    real(8), intent (in) :: om
    real(8), intent (in) :: kx
    real(8), intent (in) :: ky
    code = sqrt(((1.0d0 / 2.0d0) * (1.0d0 + (1.0d0 / sqrt((1.0d0 + ((((2.0d0 * l) / om) ** 2.0d0) * ((sin(kx) ** 2.0d0) + (sin(ky) ** 2.0d0)))))))))
end function
public static double code(double l, double Om, double kx, double ky) {
	return Math.sqrt(((1.0 / 2.0) * (1.0 + (1.0 / Math.sqrt((1.0 + (Math.pow(((2.0 * l) / Om), 2.0) * (Math.pow(Math.sin(kx), 2.0) + Math.pow(Math.sin(ky), 2.0)))))))));
}
def code(l, Om, kx, ky):
	return math.sqrt(((1.0 / 2.0) * (1.0 + (1.0 / math.sqrt((1.0 + (math.pow(((2.0 * l) / Om), 2.0) * (math.pow(math.sin(kx), 2.0) + math.pow(math.sin(ky), 2.0)))))))))
function code(l, Om, kx, ky)
	return sqrt(Float64(Float64(1.0 / 2.0) * Float64(1.0 + Float64(1.0 / sqrt(Float64(1.0 + Float64((Float64(Float64(2.0 * l) / Om) ^ 2.0) * Float64((sin(kx) ^ 2.0) + (sin(ky) ^ 2.0)))))))))
end
function tmp = code(l, Om, kx, ky)
	tmp = sqrt(((1.0 / 2.0) * (1.0 + (1.0 / sqrt((1.0 + ((((2.0 * l) / Om) ^ 2.0) * ((sin(kx) ^ 2.0) + (sin(ky) ^ 2.0)))))))));
end
code[l_, Om_, kx_, ky_] := N[Sqrt[N[(N[(1.0 / 2.0), $MachinePrecision] * N[(1.0 + N[(1.0 / N[Sqrt[N[(1.0 + N[(N[Power[N[(N[(2.0 * l), $MachinePrecision] / Om), $MachinePrecision], 2.0], $MachinePrecision] * N[(N[Power[N[Sin[kx], $MachinePrecision], 2.0], $MachinePrecision] + N[Power[N[Sin[ky], $MachinePrecision], 2.0], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]
\begin{array}{l}

\\
\sqrt{\frac{1}{2} \cdot \left(1 + \frac{1}{\sqrt{1 + {\left(\frac{2 \cdot \ell}{Om}\right)}^{2} \cdot \left({\sin kx}^{2} + {\sin ky}^{2}\right)}}\right)}
\end{array}

Alternative 1: 99.4% accurate, 1.0× speedup?

\[\begin{array}{l} Om_m = \left|Om\right| \\ l_m = \left|\ell\right| \\ \begin{array}{l} t_0 := \frac{l\_m \cdot 2}{Om\_m}\\ \mathbf{if}\;t\_0 \leq 10^{+133}:\\ \;\;\;\;\sqrt{\left(\frac{1}{\sqrt{\left({\sin ky}^{2} + {\sin kx}^{2}\right) \cdot {t\_0}^{2} + 1}} + 1\right) \cdot \frac{1}{2}}\\ \mathbf{else}:\\ \;\;\;\;\sqrt{0.5}\\ \end{array} \end{array} \]
Om_m = (fabs.f64 Om)
l_m = (fabs.f64 l)
(FPCore (l_m Om_m kx ky)
 :precision binary64
 (let* ((t_0 (/ (* l_m 2.0) Om_m)))
   (if (<= t_0 1e+133)
     (sqrt
      (*
       (+
        (/
         1.0
         (sqrt
          (+ (* (+ (pow (sin ky) 2.0) (pow (sin kx) 2.0)) (pow t_0 2.0)) 1.0)))
        1.0)
       (/ 1.0 2.0)))
     (sqrt 0.5))))
Om_m = fabs(Om);
l_m = fabs(l);
double code(double l_m, double Om_m, double kx, double ky) {
	double t_0 = (l_m * 2.0) / Om_m;
	double tmp;
	if (t_0 <= 1e+133) {
		tmp = sqrt((((1.0 / sqrt((((pow(sin(ky), 2.0) + pow(sin(kx), 2.0)) * pow(t_0, 2.0)) + 1.0))) + 1.0) * (1.0 / 2.0)));
	} else {
		tmp = sqrt(0.5);
	}
	return tmp;
}
Om_m = abs(om)
l_m = abs(l)
real(8) function code(l_m, om_m, kx, ky)
    real(8), intent (in) :: l_m
    real(8), intent (in) :: om_m
    real(8), intent (in) :: kx
    real(8), intent (in) :: ky
    real(8) :: t_0
    real(8) :: tmp
    t_0 = (l_m * 2.0d0) / om_m
    if (t_0 <= 1d+133) then
        tmp = sqrt((((1.0d0 / sqrt(((((sin(ky) ** 2.0d0) + (sin(kx) ** 2.0d0)) * (t_0 ** 2.0d0)) + 1.0d0))) + 1.0d0) * (1.0d0 / 2.0d0)))
    else
        tmp = sqrt(0.5d0)
    end if
    code = tmp
end function
Om_m = Math.abs(Om);
l_m = Math.abs(l);
public static double code(double l_m, double Om_m, double kx, double ky) {
	double t_0 = (l_m * 2.0) / Om_m;
	double tmp;
	if (t_0 <= 1e+133) {
		tmp = Math.sqrt((((1.0 / Math.sqrt((((Math.pow(Math.sin(ky), 2.0) + Math.pow(Math.sin(kx), 2.0)) * Math.pow(t_0, 2.0)) + 1.0))) + 1.0) * (1.0 / 2.0)));
	} else {
		tmp = Math.sqrt(0.5);
	}
	return tmp;
}
Om_m = math.fabs(Om)
l_m = math.fabs(l)
def code(l_m, Om_m, kx, ky):
	t_0 = (l_m * 2.0) / Om_m
	tmp = 0
	if t_0 <= 1e+133:
		tmp = math.sqrt((((1.0 / math.sqrt((((math.pow(math.sin(ky), 2.0) + math.pow(math.sin(kx), 2.0)) * math.pow(t_0, 2.0)) + 1.0))) + 1.0) * (1.0 / 2.0)))
	else:
		tmp = math.sqrt(0.5)
	return tmp
Om_m = abs(Om)
l_m = abs(l)
function code(l_m, Om_m, kx, ky)
	t_0 = Float64(Float64(l_m * 2.0) / Om_m)
	tmp = 0.0
	if (t_0 <= 1e+133)
		tmp = sqrt(Float64(Float64(Float64(1.0 / sqrt(Float64(Float64(Float64((sin(ky) ^ 2.0) + (sin(kx) ^ 2.0)) * (t_0 ^ 2.0)) + 1.0))) + 1.0) * Float64(1.0 / 2.0)));
	else
		tmp = sqrt(0.5);
	end
	return tmp
end
Om_m = abs(Om);
l_m = abs(l);
function tmp_2 = code(l_m, Om_m, kx, ky)
	t_0 = (l_m * 2.0) / Om_m;
	tmp = 0.0;
	if (t_0 <= 1e+133)
		tmp = sqrt((((1.0 / sqrt(((((sin(ky) ^ 2.0) + (sin(kx) ^ 2.0)) * (t_0 ^ 2.0)) + 1.0))) + 1.0) * (1.0 / 2.0)));
	else
		tmp = sqrt(0.5);
	end
	tmp_2 = tmp;
end
Om_m = N[Abs[Om], $MachinePrecision]
l_m = N[Abs[l], $MachinePrecision]
code[l$95$m_, Om$95$m_, kx_, ky_] := Block[{t$95$0 = N[(N[(l$95$m * 2.0), $MachinePrecision] / Om$95$m), $MachinePrecision]}, If[LessEqual[t$95$0, 1e+133], N[Sqrt[N[(N[(N[(1.0 / N[Sqrt[N[(N[(N[(N[Power[N[Sin[ky], $MachinePrecision], 2.0], $MachinePrecision] + N[Power[N[Sin[kx], $MachinePrecision], 2.0], $MachinePrecision]), $MachinePrecision] * N[Power[t$95$0, 2.0], $MachinePrecision]), $MachinePrecision] + 1.0), $MachinePrecision]], $MachinePrecision]), $MachinePrecision] + 1.0), $MachinePrecision] * N[(1.0 / 2.0), $MachinePrecision]), $MachinePrecision]], $MachinePrecision], N[Sqrt[0.5], $MachinePrecision]]]
\begin{array}{l}
Om_m = \left|Om\right|
\\
l_m = \left|\ell\right|

\\
\begin{array}{l}
t_0 := \frac{l\_m \cdot 2}{Om\_m}\\
\mathbf{if}\;t\_0 \leq 10^{+133}:\\
\;\;\;\;\sqrt{\left(\frac{1}{\sqrt{\left({\sin ky}^{2} + {\sin kx}^{2}\right) \cdot {t\_0}^{2} + 1}} + 1\right) \cdot \frac{1}{2}}\\

\mathbf{else}:\\
\;\;\;\;\sqrt{0.5}\\


\end{array}
\end{array}
Derivation
  1. Split input into 2 regimes
  2. if (/.f64 (*.f64 #s(literal 2 binary64) l) Om) < 1e133

    1. Initial program 99.5%

      \[\sqrt{\frac{1}{2} \cdot \left(1 + \frac{1}{\sqrt{1 + {\left(\frac{2 \cdot \ell}{Om}\right)}^{2} \cdot \left({\sin kx}^{2} + {\sin ky}^{2}\right)}}\right)} \]
    2. Add Preprocessing

    if 1e133 < (/.f64 (*.f64 #s(literal 2 binary64) l) Om)

    1. Initial program 91.3%

      \[\sqrt{\frac{1}{2} \cdot \left(1 + \frac{1}{\sqrt{1 + {\left(\frac{2 \cdot \ell}{Om}\right)}^{2} \cdot \left({\sin kx}^{2} + {\sin ky}^{2}\right)}}\right)} \]
    2. Add Preprocessing
    3. Taylor expanded in Om around 0

      \[\leadsto \sqrt{\color{blue}{\frac{1}{2}}} \]
    4. Step-by-step derivation
      1. Applied rewrites98.3%

        \[\leadsto \sqrt{\color{blue}{0.5}} \]
    5. Recombined 2 regimes into one program.
    6. Final simplification99.3%

      \[\leadsto \begin{array}{l} \mathbf{if}\;\frac{\ell \cdot 2}{Om} \leq 10^{+133}:\\ \;\;\;\;\sqrt{\left(\frac{1}{\sqrt{\left({\sin ky}^{2} + {\sin kx}^{2}\right) \cdot {\left(\frac{\ell \cdot 2}{Om}\right)}^{2} + 1}} + 1\right) \cdot \frac{1}{2}}\\ \mathbf{else}:\\ \;\;\;\;\sqrt{0.5}\\ \end{array} \]
    7. Add Preprocessing

    Alternative 2: 95.7% accurate, 0.7× speedup?

    \[\begin{array}{l} Om_m = \left|Om\right| \\ l_m = \left|\ell\right| \\ \begin{array}{l} \mathbf{if}\;\left({\sin ky}^{2} + {\sin kx}^{2}\right) \cdot {\left(\frac{l\_m \cdot 2}{Om\_m}\right)}^{2} \leq 0.0005:\\ \;\;\;\;\sqrt{1}\\ \mathbf{else}:\\ \;\;\;\;\sqrt{\frac{0.5}{\mathsf{fma}\left(\frac{0.25}{\sin ky}, \frac{Om\_m}{l\_m \cdot l\_m}, \frac{\sin ky}{Om\_m} \cdot 2\right) \cdot l\_m} + 0.5}\\ \end{array} \end{array} \]
    Om_m = (fabs.f64 Om)
    l_m = (fabs.f64 l)
    (FPCore (l_m Om_m kx ky)
     :precision binary64
     (if (<=
          (*
           (+ (pow (sin ky) 2.0) (pow (sin kx) 2.0))
           (pow (/ (* l_m 2.0) Om_m) 2.0))
          0.0005)
       (sqrt 1.0)
       (sqrt
        (+
         (/
          0.5
          (*
           (fma (/ 0.25 (sin ky)) (/ Om_m (* l_m l_m)) (* (/ (sin ky) Om_m) 2.0))
           l_m))
         0.5))))
    Om_m = fabs(Om);
    l_m = fabs(l);
    double code(double l_m, double Om_m, double kx, double ky) {
    	double tmp;
    	if (((pow(sin(ky), 2.0) + pow(sin(kx), 2.0)) * pow(((l_m * 2.0) / Om_m), 2.0)) <= 0.0005) {
    		tmp = sqrt(1.0);
    	} else {
    		tmp = sqrt(((0.5 / (fma((0.25 / sin(ky)), (Om_m / (l_m * l_m)), ((sin(ky) / Om_m) * 2.0)) * l_m)) + 0.5));
    	}
    	return tmp;
    }
    
    Om_m = abs(Om)
    l_m = abs(l)
    function code(l_m, Om_m, kx, ky)
    	tmp = 0.0
    	if (Float64(Float64((sin(ky) ^ 2.0) + (sin(kx) ^ 2.0)) * (Float64(Float64(l_m * 2.0) / Om_m) ^ 2.0)) <= 0.0005)
    		tmp = sqrt(1.0);
    	else
    		tmp = sqrt(Float64(Float64(0.5 / Float64(fma(Float64(0.25 / sin(ky)), Float64(Om_m / Float64(l_m * l_m)), Float64(Float64(sin(ky) / Om_m) * 2.0)) * l_m)) + 0.5));
    	end
    	return tmp
    end
    
    Om_m = N[Abs[Om], $MachinePrecision]
    l_m = N[Abs[l], $MachinePrecision]
    code[l$95$m_, Om$95$m_, kx_, ky_] := If[LessEqual[N[(N[(N[Power[N[Sin[ky], $MachinePrecision], 2.0], $MachinePrecision] + N[Power[N[Sin[kx], $MachinePrecision], 2.0], $MachinePrecision]), $MachinePrecision] * N[Power[N[(N[(l$95$m * 2.0), $MachinePrecision] / Om$95$m), $MachinePrecision], 2.0], $MachinePrecision]), $MachinePrecision], 0.0005], N[Sqrt[1.0], $MachinePrecision], N[Sqrt[N[(N[(0.5 / N[(N[(N[(0.25 / N[Sin[ky], $MachinePrecision]), $MachinePrecision] * N[(Om$95$m / N[(l$95$m * l$95$m), $MachinePrecision]), $MachinePrecision] + N[(N[(N[Sin[ky], $MachinePrecision] / Om$95$m), $MachinePrecision] * 2.0), $MachinePrecision]), $MachinePrecision] * l$95$m), $MachinePrecision]), $MachinePrecision] + 0.5), $MachinePrecision]], $MachinePrecision]]
    
    \begin{array}{l}
    Om_m = \left|Om\right|
    \\
    l_m = \left|\ell\right|
    
    \\
    \begin{array}{l}
    \mathbf{if}\;\left({\sin ky}^{2} + {\sin kx}^{2}\right) \cdot {\left(\frac{l\_m \cdot 2}{Om\_m}\right)}^{2} \leq 0.0005:\\
    \;\;\;\;\sqrt{1}\\
    
    \mathbf{else}:\\
    \;\;\;\;\sqrt{\frac{0.5}{\mathsf{fma}\left(\frac{0.25}{\sin ky}, \frac{Om\_m}{l\_m \cdot l\_m}, \frac{\sin ky}{Om\_m} \cdot 2\right) \cdot l\_m} + 0.5}\\
    
    
    \end{array}
    \end{array}
    
    Derivation
    1. Split input into 2 regimes
    2. if (*.f64 (pow.f64 (/.f64 (*.f64 #s(literal 2 binary64) l) Om) #s(literal 2 binary64)) (+.f64 (pow.f64 (sin.f64 kx) #s(literal 2 binary64)) (pow.f64 (sin.f64 ky) #s(literal 2 binary64)))) < 5.0000000000000001e-4

      1. Initial program 100.0%

        \[\sqrt{\frac{1}{2} \cdot \left(1 + \frac{1}{\sqrt{1 + {\left(\frac{2 \cdot \ell}{Om}\right)}^{2} \cdot \left({\sin kx}^{2} + {\sin ky}^{2}\right)}}\right)} \]
      2. Add Preprocessing
      3. Taylor expanded in Om around inf

        \[\leadsto \sqrt{\color{blue}{1}} \]
      4. Step-by-step derivation
        1. Applied rewrites99.5%

          \[\leadsto \sqrt{\color{blue}{1}} \]

        if 5.0000000000000001e-4 < (*.f64 (pow.f64 (/.f64 (*.f64 #s(literal 2 binary64) l) Om) #s(literal 2 binary64)) (+.f64 (pow.f64 (sin.f64 kx) #s(literal 2 binary64)) (pow.f64 (sin.f64 ky) #s(literal 2 binary64))))

        1. Initial program 95.9%

          \[\sqrt{\frac{1}{2} \cdot \left(1 + \frac{1}{\sqrt{1 + {\left(\frac{2 \cdot \ell}{Om}\right)}^{2} \cdot \left({\sin kx}^{2} + {\sin ky}^{2}\right)}}\right)} \]
        2. Add Preprocessing
        3. Taylor expanded in kx around 0

          \[\leadsto \sqrt{\color{blue}{\frac{1}{2} \cdot \left(1 + \sqrt{\frac{1}{1 + 4 \cdot \frac{{\ell}^{2} \cdot {\sin ky}^{2}}{{Om}^{2}}}}\right)}} \]
        4. Step-by-step derivation
          1. +-commutativeN/A

            \[\leadsto \sqrt{\frac{1}{2} \cdot \color{blue}{\left(\sqrt{\frac{1}{1 + 4 \cdot \frac{{\ell}^{2} \cdot {\sin ky}^{2}}{{Om}^{2}}}} + 1\right)}} \]
          2. distribute-rgt-inN/A

            \[\leadsto \sqrt{\color{blue}{\sqrt{\frac{1}{1 + 4 \cdot \frac{{\ell}^{2} \cdot {\sin ky}^{2}}{{Om}^{2}}}} \cdot \frac{1}{2} + 1 \cdot \frac{1}{2}}} \]
          3. metadata-evalN/A

            \[\leadsto \sqrt{\sqrt{\frac{1}{1 + 4 \cdot \frac{{\ell}^{2} \cdot {\sin ky}^{2}}{{Om}^{2}}}} \cdot \frac{1}{2} + \color{blue}{\frac{1}{2}}} \]
          4. lower-fma.f64N/A

            \[\leadsto \sqrt{\color{blue}{\mathsf{fma}\left(\sqrt{\frac{1}{1 + 4 \cdot \frac{{\ell}^{2} \cdot {\sin ky}^{2}}{{Om}^{2}}}}, \frac{1}{2}, \frac{1}{2}\right)}} \]
        5. Applied rewrites74.6%

          \[\leadsto \sqrt{\color{blue}{\mathsf{fma}\left(\sqrt{\frac{1}{\mathsf{fma}\left(\frac{{\sin ky}^{2}}{Om} \cdot \frac{\ell \cdot \ell}{Om}, 4, 1\right)}}, 0.5, 0.5\right)}} \]
        6. Step-by-step derivation
          1. Applied rewrites76.4%

            \[\leadsto \sqrt{\frac{0.5}{\sqrt{\mathsf{fma}\left(4, {\sin ky}^{2} \cdot {\left(\frac{Om}{\ell}\right)}^{-2}, 1\right)}} + \color{blue}{0.5}} \]
          2. Taylor expanded in l around inf

            \[\leadsto \sqrt{\frac{\frac{1}{2}}{\ell \cdot \left(\frac{1}{4} \cdot \frac{Om}{{\ell}^{2} \cdot \sin ky} + 2 \cdot \frac{\sin ky}{Om}\right)} + \frac{1}{2}} \]
          3. Step-by-step derivation
            1. Applied rewrites92.4%

              \[\leadsto \sqrt{\frac{0.5}{\mathsf{fma}\left(\frac{0.25}{\sin ky}, \frac{Om}{\ell \cdot \ell}, \frac{\sin ky}{Om} \cdot 2\right) \cdot \ell} + 0.5} \]
          4. Recombined 2 regimes into one program.
          5. Final simplification96.1%

            \[\leadsto \begin{array}{l} \mathbf{if}\;\left({\sin ky}^{2} + {\sin kx}^{2}\right) \cdot {\left(\frac{\ell \cdot 2}{Om}\right)}^{2} \leq 0.0005:\\ \;\;\;\;\sqrt{1}\\ \mathbf{else}:\\ \;\;\;\;\sqrt{\frac{0.5}{\mathsf{fma}\left(\frac{0.25}{\sin ky}, \frac{Om}{\ell \cdot \ell}, \frac{\sin ky}{Om} \cdot 2\right) \cdot \ell} + 0.5}\\ \end{array} \]
          6. Add Preprocessing

          Alternative 3: 91.8% accurate, 0.9× speedup?

          \[\begin{array}{l} Om_m = \left|Om\right| \\ l_m = \left|\ell\right| \\ \begin{array}{l} \mathbf{if}\;\left({\sin ky}^{2} + {\sin kx}^{2}\right) \cdot {\left(\frac{l\_m \cdot 2}{Om\_m}\right)}^{2} \leq 2:\\ \;\;\;\;\sqrt{1}\\ \mathbf{else}:\\ \;\;\;\;\sqrt{\mathsf{fma}\left(\frac{Om\_m}{\sin ky \cdot l\_m}, 0.25, 0.5\right)}\\ \end{array} \end{array} \]
          Om_m = (fabs.f64 Om)
          l_m = (fabs.f64 l)
          (FPCore (l_m Om_m kx ky)
           :precision binary64
           (if (<=
                (*
                 (+ (pow (sin ky) 2.0) (pow (sin kx) 2.0))
                 (pow (/ (* l_m 2.0) Om_m) 2.0))
                2.0)
             (sqrt 1.0)
             (sqrt (fma (/ Om_m (* (sin ky) l_m)) 0.25 0.5))))
          Om_m = fabs(Om);
          l_m = fabs(l);
          double code(double l_m, double Om_m, double kx, double ky) {
          	double tmp;
          	if (((pow(sin(ky), 2.0) + pow(sin(kx), 2.0)) * pow(((l_m * 2.0) / Om_m), 2.0)) <= 2.0) {
          		tmp = sqrt(1.0);
          	} else {
          		tmp = sqrt(fma((Om_m / (sin(ky) * l_m)), 0.25, 0.5));
          	}
          	return tmp;
          }
          
          Om_m = abs(Om)
          l_m = abs(l)
          function code(l_m, Om_m, kx, ky)
          	tmp = 0.0
          	if (Float64(Float64((sin(ky) ^ 2.0) + (sin(kx) ^ 2.0)) * (Float64(Float64(l_m * 2.0) / Om_m) ^ 2.0)) <= 2.0)
          		tmp = sqrt(1.0);
          	else
          		tmp = sqrt(fma(Float64(Om_m / Float64(sin(ky) * l_m)), 0.25, 0.5));
          	end
          	return tmp
          end
          
          Om_m = N[Abs[Om], $MachinePrecision]
          l_m = N[Abs[l], $MachinePrecision]
          code[l$95$m_, Om$95$m_, kx_, ky_] := If[LessEqual[N[(N[(N[Power[N[Sin[ky], $MachinePrecision], 2.0], $MachinePrecision] + N[Power[N[Sin[kx], $MachinePrecision], 2.0], $MachinePrecision]), $MachinePrecision] * N[Power[N[(N[(l$95$m * 2.0), $MachinePrecision] / Om$95$m), $MachinePrecision], 2.0], $MachinePrecision]), $MachinePrecision], 2.0], N[Sqrt[1.0], $MachinePrecision], N[Sqrt[N[(N[(Om$95$m / N[(N[Sin[ky], $MachinePrecision] * l$95$m), $MachinePrecision]), $MachinePrecision] * 0.25 + 0.5), $MachinePrecision]], $MachinePrecision]]
          
          \begin{array}{l}
          Om_m = \left|Om\right|
          \\
          l_m = \left|\ell\right|
          
          \\
          \begin{array}{l}
          \mathbf{if}\;\left({\sin ky}^{2} + {\sin kx}^{2}\right) \cdot {\left(\frac{l\_m \cdot 2}{Om\_m}\right)}^{2} \leq 2:\\
          \;\;\;\;\sqrt{1}\\
          
          \mathbf{else}:\\
          \;\;\;\;\sqrt{\mathsf{fma}\left(\frac{Om\_m}{\sin ky \cdot l\_m}, 0.25, 0.5\right)}\\
          
          
          \end{array}
          \end{array}
          
          Derivation
          1. Split input into 2 regimes
          2. if (*.f64 (pow.f64 (/.f64 (*.f64 #s(literal 2 binary64) l) Om) #s(literal 2 binary64)) (+.f64 (pow.f64 (sin.f64 kx) #s(literal 2 binary64)) (pow.f64 (sin.f64 ky) #s(literal 2 binary64)))) < 2

            1. Initial program 100.0%

              \[\sqrt{\frac{1}{2} \cdot \left(1 + \frac{1}{\sqrt{1 + {\left(\frac{2 \cdot \ell}{Om}\right)}^{2} \cdot \left({\sin kx}^{2} + {\sin ky}^{2}\right)}}\right)} \]
            2. Add Preprocessing
            3. Taylor expanded in Om around inf

              \[\leadsto \sqrt{\color{blue}{1}} \]
            4. Step-by-step derivation
              1. Applied rewrites99.0%

                \[\leadsto \sqrt{\color{blue}{1}} \]

              if 2 < (*.f64 (pow.f64 (/.f64 (*.f64 #s(literal 2 binary64) l) Om) #s(literal 2 binary64)) (+.f64 (pow.f64 (sin.f64 kx) #s(literal 2 binary64)) (pow.f64 (sin.f64 ky) #s(literal 2 binary64))))

              1. Initial program 95.8%

                \[\sqrt{\frac{1}{2} \cdot \left(1 + \frac{1}{\sqrt{1 + {\left(\frac{2 \cdot \ell}{Om}\right)}^{2} \cdot \left({\sin kx}^{2} + {\sin ky}^{2}\right)}}\right)} \]
              2. Add Preprocessing
              3. Taylor expanded in kx around 0

                \[\leadsto \sqrt{\color{blue}{\frac{1}{2} \cdot \left(1 + \sqrt{\frac{1}{1 + 4 \cdot \frac{{\ell}^{2} \cdot {\sin ky}^{2}}{{Om}^{2}}}}\right)}} \]
              4. Step-by-step derivation
                1. +-commutativeN/A

                  \[\leadsto \sqrt{\frac{1}{2} \cdot \color{blue}{\left(\sqrt{\frac{1}{1 + 4 \cdot \frac{{\ell}^{2} \cdot {\sin ky}^{2}}{{Om}^{2}}}} + 1\right)}} \]
                2. distribute-rgt-inN/A

                  \[\leadsto \sqrt{\color{blue}{\sqrt{\frac{1}{1 + 4 \cdot \frac{{\ell}^{2} \cdot {\sin ky}^{2}}{{Om}^{2}}}} \cdot \frac{1}{2} + 1 \cdot \frac{1}{2}}} \]
                3. metadata-evalN/A

                  \[\leadsto \sqrt{\sqrt{\frac{1}{1 + 4 \cdot \frac{{\ell}^{2} \cdot {\sin ky}^{2}}{{Om}^{2}}}} \cdot \frac{1}{2} + \color{blue}{\frac{1}{2}}} \]
                4. lower-fma.f64N/A

                  \[\leadsto \sqrt{\color{blue}{\mathsf{fma}\left(\sqrt{\frac{1}{1 + 4 \cdot \frac{{\ell}^{2} \cdot {\sin ky}^{2}}{{Om}^{2}}}}, \frac{1}{2}, \frac{1}{2}\right)}} \]
              5. Applied rewrites75.0%

                \[\leadsto \sqrt{\color{blue}{\mathsf{fma}\left(\sqrt{\frac{1}{\mathsf{fma}\left(\frac{{\sin ky}^{2}}{Om} \cdot \frac{\ell \cdot \ell}{Om}, 4, 1\right)}}, 0.5, 0.5\right)}} \]
              6. Taylor expanded in Om around 0

                \[\leadsto \sqrt{\frac{1}{2} + \color{blue}{\frac{1}{4} \cdot \frac{Om}{\ell \cdot \sin ky}}} \]
              7. Step-by-step derivation
                1. Applied rewrites85.4%

                  \[\leadsto \sqrt{\mathsf{fma}\left(\frac{Om}{\sin ky \cdot \ell}, \color{blue}{0.25}, 0.5\right)} \]
              8. Recombined 2 regimes into one program.
              9. Final simplification92.6%

                \[\leadsto \begin{array}{l} \mathbf{if}\;\left({\sin ky}^{2} + {\sin kx}^{2}\right) \cdot {\left(\frac{\ell \cdot 2}{Om}\right)}^{2} \leq 2:\\ \;\;\;\;\sqrt{1}\\ \mathbf{else}:\\ \;\;\;\;\sqrt{\mathsf{fma}\left(\frac{Om}{\sin ky \cdot \ell}, 0.25, 0.5\right)}\\ \end{array} \]
              10. Add Preprocessing

              Alternative 4: 91.7% accurate, 0.9× speedup?

              \[\begin{array}{l} Om_m = \left|Om\right| \\ l_m = \left|\ell\right| \\ \begin{array}{l} \mathbf{if}\;\left({\sin ky}^{2} + {\sin kx}^{2}\right) \cdot {\left(\frac{l\_m \cdot 2}{Om\_m}\right)}^{2} \leq 2:\\ \;\;\;\;\sqrt{1}\\ \mathbf{else}:\\ \;\;\;\;\sqrt{\mathsf{fma}\left(\frac{Om\_m}{\sin ky \cdot l\_m}, -0.25, 0.5\right)}\\ \end{array} \end{array} \]
              Om_m = (fabs.f64 Om)
              l_m = (fabs.f64 l)
              (FPCore (l_m Om_m kx ky)
               :precision binary64
               (if (<=
                    (*
                     (+ (pow (sin ky) 2.0) (pow (sin kx) 2.0))
                     (pow (/ (* l_m 2.0) Om_m) 2.0))
                    2.0)
                 (sqrt 1.0)
                 (sqrt (fma (/ Om_m (* (sin ky) l_m)) -0.25 0.5))))
              Om_m = fabs(Om);
              l_m = fabs(l);
              double code(double l_m, double Om_m, double kx, double ky) {
              	double tmp;
              	if (((pow(sin(ky), 2.0) + pow(sin(kx), 2.0)) * pow(((l_m * 2.0) / Om_m), 2.0)) <= 2.0) {
              		tmp = sqrt(1.0);
              	} else {
              		tmp = sqrt(fma((Om_m / (sin(ky) * l_m)), -0.25, 0.5));
              	}
              	return tmp;
              }
              
              Om_m = abs(Om)
              l_m = abs(l)
              function code(l_m, Om_m, kx, ky)
              	tmp = 0.0
              	if (Float64(Float64((sin(ky) ^ 2.0) + (sin(kx) ^ 2.0)) * (Float64(Float64(l_m * 2.0) / Om_m) ^ 2.0)) <= 2.0)
              		tmp = sqrt(1.0);
              	else
              		tmp = sqrt(fma(Float64(Om_m / Float64(sin(ky) * l_m)), -0.25, 0.5));
              	end
              	return tmp
              end
              
              Om_m = N[Abs[Om], $MachinePrecision]
              l_m = N[Abs[l], $MachinePrecision]
              code[l$95$m_, Om$95$m_, kx_, ky_] := If[LessEqual[N[(N[(N[Power[N[Sin[ky], $MachinePrecision], 2.0], $MachinePrecision] + N[Power[N[Sin[kx], $MachinePrecision], 2.0], $MachinePrecision]), $MachinePrecision] * N[Power[N[(N[(l$95$m * 2.0), $MachinePrecision] / Om$95$m), $MachinePrecision], 2.0], $MachinePrecision]), $MachinePrecision], 2.0], N[Sqrt[1.0], $MachinePrecision], N[Sqrt[N[(N[(Om$95$m / N[(N[Sin[ky], $MachinePrecision] * l$95$m), $MachinePrecision]), $MachinePrecision] * -0.25 + 0.5), $MachinePrecision]], $MachinePrecision]]
              
              \begin{array}{l}
              Om_m = \left|Om\right|
              \\
              l_m = \left|\ell\right|
              
              \\
              \begin{array}{l}
              \mathbf{if}\;\left({\sin ky}^{2} + {\sin kx}^{2}\right) \cdot {\left(\frac{l\_m \cdot 2}{Om\_m}\right)}^{2} \leq 2:\\
              \;\;\;\;\sqrt{1}\\
              
              \mathbf{else}:\\
              \;\;\;\;\sqrt{\mathsf{fma}\left(\frac{Om\_m}{\sin ky \cdot l\_m}, -0.25, 0.5\right)}\\
              
              
              \end{array}
              \end{array}
              
              Derivation
              1. Split input into 2 regimes
              2. if (*.f64 (pow.f64 (/.f64 (*.f64 #s(literal 2 binary64) l) Om) #s(literal 2 binary64)) (+.f64 (pow.f64 (sin.f64 kx) #s(literal 2 binary64)) (pow.f64 (sin.f64 ky) #s(literal 2 binary64)))) < 2

                1. Initial program 100.0%

                  \[\sqrt{\frac{1}{2} \cdot \left(1 + \frac{1}{\sqrt{1 + {\left(\frac{2 \cdot \ell}{Om}\right)}^{2} \cdot \left({\sin kx}^{2} + {\sin ky}^{2}\right)}}\right)} \]
                2. Add Preprocessing
                3. Taylor expanded in Om around inf

                  \[\leadsto \sqrt{\color{blue}{1}} \]
                4. Step-by-step derivation
                  1. Applied rewrites99.0%

                    \[\leadsto \sqrt{\color{blue}{1}} \]

                  if 2 < (*.f64 (pow.f64 (/.f64 (*.f64 #s(literal 2 binary64) l) Om) #s(literal 2 binary64)) (+.f64 (pow.f64 (sin.f64 kx) #s(literal 2 binary64)) (pow.f64 (sin.f64 ky) #s(literal 2 binary64))))

                  1. Initial program 95.8%

                    \[\sqrt{\frac{1}{2} \cdot \left(1 + \frac{1}{\sqrt{1 + {\left(\frac{2 \cdot \ell}{Om}\right)}^{2} \cdot \left({\sin kx}^{2} + {\sin ky}^{2}\right)}}\right)} \]
                  2. Add Preprocessing
                  3. Taylor expanded in kx around 0

                    \[\leadsto \sqrt{\color{blue}{\frac{1}{2} \cdot \left(1 + \sqrt{\frac{1}{1 + 4 \cdot \frac{{\ell}^{2} \cdot {\sin ky}^{2}}{{Om}^{2}}}}\right)}} \]
                  4. Step-by-step derivation
                    1. +-commutativeN/A

                      \[\leadsto \sqrt{\frac{1}{2} \cdot \color{blue}{\left(\sqrt{\frac{1}{1 + 4 \cdot \frac{{\ell}^{2} \cdot {\sin ky}^{2}}{{Om}^{2}}}} + 1\right)}} \]
                    2. distribute-rgt-inN/A

                      \[\leadsto \sqrt{\color{blue}{\sqrt{\frac{1}{1 + 4 \cdot \frac{{\ell}^{2} \cdot {\sin ky}^{2}}{{Om}^{2}}}} \cdot \frac{1}{2} + 1 \cdot \frac{1}{2}}} \]
                    3. metadata-evalN/A

                      \[\leadsto \sqrt{\sqrt{\frac{1}{1 + 4 \cdot \frac{{\ell}^{2} \cdot {\sin ky}^{2}}{{Om}^{2}}}} \cdot \frac{1}{2} + \color{blue}{\frac{1}{2}}} \]
                    4. lower-fma.f64N/A

                      \[\leadsto \sqrt{\color{blue}{\mathsf{fma}\left(\sqrt{\frac{1}{1 + 4 \cdot \frac{{\ell}^{2} \cdot {\sin ky}^{2}}{{Om}^{2}}}}, \frac{1}{2}, \frac{1}{2}\right)}} \]
                  5. Applied rewrites75.0%

                    \[\leadsto \sqrt{\color{blue}{\mathsf{fma}\left(\sqrt{\frac{1}{\mathsf{fma}\left(\frac{{\sin ky}^{2}}{Om} \cdot \frac{\ell \cdot \ell}{Om}, 4, 1\right)}}, 0.5, 0.5\right)}} \]
                  6. Taylor expanded in l around -inf

                    \[\leadsto \sqrt{\frac{1}{2} + \color{blue}{\frac{-1}{4} \cdot \frac{Om}{\ell \cdot \sin ky}}} \]
                  7. Step-by-step derivation
                    1. Applied rewrites86.1%

                      \[\leadsto \sqrt{\mathsf{fma}\left(\frac{Om}{\sin ky \cdot \ell}, \color{blue}{-0.25}, 0.5\right)} \]
                  8. Recombined 2 regimes into one program.
                  9. Final simplification92.9%

                    \[\leadsto \begin{array}{l} \mathbf{if}\;\left({\sin ky}^{2} + {\sin kx}^{2}\right) \cdot {\left(\frac{\ell \cdot 2}{Om}\right)}^{2} \leq 2:\\ \;\;\;\;\sqrt{1}\\ \mathbf{else}:\\ \;\;\;\;\sqrt{\mathsf{fma}\left(\frac{Om}{\sin ky \cdot \ell}, -0.25, 0.5\right)}\\ \end{array} \]
                  10. Add Preprocessing

                  Alternative 5: 98.2% accurate, 1.1× speedup?

                  \[\begin{array}{l} Om_m = \left|Om\right| \\ l_m = \left|\ell\right| \\ \begin{array}{l} \mathbf{if}\;\left({\sin ky}^{2} + {\sin kx}^{2}\right) \cdot {\left(\frac{l\_m \cdot 2}{Om\_m}\right)}^{2} \leq 2:\\ \;\;\;\;\sqrt{1}\\ \mathbf{else}:\\ \;\;\;\;\sqrt{0.5}\\ \end{array} \end{array} \]
                  Om_m = (fabs.f64 Om)
                  l_m = (fabs.f64 l)
                  (FPCore (l_m Om_m kx ky)
                   :precision binary64
                   (if (<=
                        (*
                         (+ (pow (sin ky) 2.0) (pow (sin kx) 2.0))
                         (pow (/ (* l_m 2.0) Om_m) 2.0))
                        2.0)
                     (sqrt 1.0)
                     (sqrt 0.5)))
                  Om_m = fabs(Om);
                  l_m = fabs(l);
                  double code(double l_m, double Om_m, double kx, double ky) {
                  	double tmp;
                  	if (((pow(sin(ky), 2.0) + pow(sin(kx), 2.0)) * pow(((l_m * 2.0) / Om_m), 2.0)) <= 2.0) {
                  		tmp = sqrt(1.0);
                  	} else {
                  		tmp = sqrt(0.5);
                  	}
                  	return tmp;
                  }
                  
                  Om_m = abs(om)
                  l_m = abs(l)
                  real(8) function code(l_m, om_m, kx, ky)
                      real(8), intent (in) :: l_m
                      real(8), intent (in) :: om_m
                      real(8), intent (in) :: kx
                      real(8), intent (in) :: ky
                      real(8) :: tmp
                      if ((((sin(ky) ** 2.0d0) + (sin(kx) ** 2.0d0)) * (((l_m * 2.0d0) / om_m) ** 2.0d0)) <= 2.0d0) then
                          tmp = sqrt(1.0d0)
                      else
                          tmp = sqrt(0.5d0)
                      end if
                      code = tmp
                  end function
                  
                  Om_m = Math.abs(Om);
                  l_m = Math.abs(l);
                  public static double code(double l_m, double Om_m, double kx, double ky) {
                  	double tmp;
                  	if (((Math.pow(Math.sin(ky), 2.0) + Math.pow(Math.sin(kx), 2.0)) * Math.pow(((l_m * 2.0) / Om_m), 2.0)) <= 2.0) {
                  		tmp = Math.sqrt(1.0);
                  	} else {
                  		tmp = Math.sqrt(0.5);
                  	}
                  	return tmp;
                  }
                  
                  Om_m = math.fabs(Om)
                  l_m = math.fabs(l)
                  def code(l_m, Om_m, kx, ky):
                  	tmp = 0
                  	if ((math.pow(math.sin(ky), 2.0) + math.pow(math.sin(kx), 2.0)) * math.pow(((l_m * 2.0) / Om_m), 2.0)) <= 2.0:
                  		tmp = math.sqrt(1.0)
                  	else:
                  		tmp = math.sqrt(0.5)
                  	return tmp
                  
                  Om_m = abs(Om)
                  l_m = abs(l)
                  function code(l_m, Om_m, kx, ky)
                  	tmp = 0.0
                  	if (Float64(Float64((sin(ky) ^ 2.0) + (sin(kx) ^ 2.0)) * (Float64(Float64(l_m * 2.0) / Om_m) ^ 2.0)) <= 2.0)
                  		tmp = sqrt(1.0);
                  	else
                  		tmp = sqrt(0.5);
                  	end
                  	return tmp
                  end
                  
                  Om_m = abs(Om);
                  l_m = abs(l);
                  function tmp_2 = code(l_m, Om_m, kx, ky)
                  	tmp = 0.0;
                  	if ((((sin(ky) ^ 2.0) + (sin(kx) ^ 2.0)) * (((l_m * 2.0) / Om_m) ^ 2.0)) <= 2.0)
                  		tmp = sqrt(1.0);
                  	else
                  		tmp = sqrt(0.5);
                  	end
                  	tmp_2 = tmp;
                  end
                  
                  Om_m = N[Abs[Om], $MachinePrecision]
                  l_m = N[Abs[l], $MachinePrecision]
                  code[l$95$m_, Om$95$m_, kx_, ky_] := If[LessEqual[N[(N[(N[Power[N[Sin[ky], $MachinePrecision], 2.0], $MachinePrecision] + N[Power[N[Sin[kx], $MachinePrecision], 2.0], $MachinePrecision]), $MachinePrecision] * N[Power[N[(N[(l$95$m * 2.0), $MachinePrecision] / Om$95$m), $MachinePrecision], 2.0], $MachinePrecision]), $MachinePrecision], 2.0], N[Sqrt[1.0], $MachinePrecision], N[Sqrt[0.5], $MachinePrecision]]
                  
                  \begin{array}{l}
                  Om_m = \left|Om\right|
                  \\
                  l_m = \left|\ell\right|
                  
                  \\
                  \begin{array}{l}
                  \mathbf{if}\;\left({\sin ky}^{2} + {\sin kx}^{2}\right) \cdot {\left(\frac{l\_m \cdot 2}{Om\_m}\right)}^{2} \leq 2:\\
                  \;\;\;\;\sqrt{1}\\
                  
                  \mathbf{else}:\\
                  \;\;\;\;\sqrt{0.5}\\
                  
                  
                  \end{array}
                  \end{array}
                  
                  Derivation
                  1. Split input into 2 regimes
                  2. if (*.f64 (pow.f64 (/.f64 (*.f64 #s(literal 2 binary64) l) Om) #s(literal 2 binary64)) (+.f64 (pow.f64 (sin.f64 kx) #s(literal 2 binary64)) (pow.f64 (sin.f64 ky) #s(literal 2 binary64)))) < 2

                    1. Initial program 100.0%

                      \[\sqrt{\frac{1}{2} \cdot \left(1 + \frac{1}{\sqrt{1 + {\left(\frac{2 \cdot \ell}{Om}\right)}^{2} \cdot \left({\sin kx}^{2} + {\sin ky}^{2}\right)}}\right)} \]
                    2. Add Preprocessing
                    3. Taylor expanded in Om around inf

                      \[\leadsto \sqrt{\color{blue}{1}} \]
                    4. Step-by-step derivation
                      1. Applied rewrites99.0%

                        \[\leadsto \sqrt{\color{blue}{1}} \]

                      if 2 < (*.f64 (pow.f64 (/.f64 (*.f64 #s(literal 2 binary64) l) Om) #s(literal 2 binary64)) (+.f64 (pow.f64 (sin.f64 kx) #s(literal 2 binary64)) (pow.f64 (sin.f64 ky) #s(literal 2 binary64))))

                      1. Initial program 95.8%

                        \[\sqrt{\frac{1}{2} \cdot \left(1 + \frac{1}{\sqrt{1 + {\left(\frac{2 \cdot \ell}{Om}\right)}^{2} \cdot \left({\sin kx}^{2} + {\sin ky}^{2}\right)}}\right)} \]
                      2. Add Preprocessing
                      3. Taylor expanded in Om around 0

                        \[\leadsto \sqrt{\color{blue}{\frac{1}{2}}} \]
                      4. Step-by-step derivation
                        1. Applied rewrites98.6%

                          \[\leadsto \sqrt{\color{blue}{0.5}} \]
                      5. Recombined 2 regimes into one program.
                      6. Final simplification98.8%

                        \[\leadsto \begin{array}{l} \mathbf{if}\;\left({\sin ky}^{2} + {\sin kx}^{2}\right) \cdot {\left(\frac{\ell \cdot 2}{Om}\right)}^{2} \leq 2:\\ \;\;\;\;\sqrt{1}\\ \mathbf{else}:\\ \;\;\;\;\sqrt{0.5}\\ \end{array} \]
                      7. Add Preprocessing

                      Alternative 6: 95.3% accurate, 1.5× speedup?

                      \[\begin{array}{l} Om_m = \left|Om\right| \\ l_m = \left|\ell\right| \\ \begin{array}{l} \mathbf{if}\;\frac{l\_m \cdot 2}{Om\_m} \leq 10^{+133}:\\ \;\;\;\;\sqrt{\frac{0.5}{\sqrt{\mathsf{fma}\left(4, {\left(\frac{Om\_m}{l\_m}\right)}^{-2} \cdot {\sin ky}^{2}, 1\right)}} + 0.5}\\ \mathbf{else}:\\ \;\;\;\;\sqrt{0.5}\\ \end{array} \end{array} \]
                      Om_m = (fabs.f64 Om)
                      l_m = (fabs.f64 l)
                      (FPCore (l_m Om_m kx ky)
                       :precision binary64
                       (if (<= (/ (* l_m 2.0) Om_m) 1e+133)
                         (sqrt
                          (+
                           (/
                            0.5
                            (sqrt (fma 4.0 (* (pow (/ Om_m l_m) -2.0) (pow (sin ky) 2.0)) 1.0)))
                           0.5))
                         (sqrt 0.5)))
                      Om_m = fabs(Om);
                      l_m = fabs(l);
                      double code(double l_m, double Om_m, double kx, double ky) {
                      	double tmp;
                      	if (((l_m * 2.0) / Om_m) <= 1e+133) {
                      		tmp = sqrt(((0.5 / sqrt(fma(4.0, (pow((Om_m / l_m), -2.0) * pow(sin(ky), 2.0)), 1.0))) + 0.5));
                      	} else {
                      		tmp = sqrt(0.5);
                      	}
                      	return tmp;
                      }
                      
                      Om_m = abs(Om)
                      l_m = abs(l)
                      function code(l_m, Om_m, kx, ky)
                      	tmp = 0.0
                      	if (Float64(Float64(l_m * 2.0) / Om_m) <= 1e+133)
                      		tmp = sqrt(Float64(Float64(0.5 / sqrt(fma(4.0, Float64((Float64(Om_m / l_m) ^ -2.0) * (sin(ky) ^ 2.0)), 1.0))) + 0.5));
                      	else
                      		tmp = sqrt(0.5);
                      	end
                      	return tmp
                      end
                      
                      Om_m = N[Abs[Om], $MachinePrecision]
                      l_m = N[Abs[l], $MachinePrecision]
                      code[l$95$m_, Om$95$m_, kx_, ky_] := If[LessEqual[N[(N[(l$95$m * 2.0), $MachinePrecision] / Om$95$m), $MachinePrecision], 1e+133], N[Sqrt[N[(N[(0.5 / N[Sqrt[N[(4.0 * N[(N[Power[N[(Om$95$m / l$95$m), $MachinePrecision], -2.0], $MachinePrecision] * N[Power[N[Sin[ky], $MachinePrecision], 2.0], $MachinePrecision]), $MachinePrecision] + 1.0), $MachinePrecision]], $MachinePrecision]), $MachinePrecision] + 0.5), $MachinePrecision]], $MachinePrecision], N[Sqrt[0.5], $MachinePrecision]]
                      
                      \begin{array}{l}
                      Om_m = \left|Om\right|
                      \\
                      l_m = \left|\ell\right|
                      
                      \\
                      \begin{array}{l}
                      \mathbf{if}\;\frac{l\_m \cdot 2}{Om\_m} \leq 10^{+133}:\\
                      \;\;\;\;\sqrt{\frac{0.5}{\sqrt{\mathsf{fma}\left(4, {\left(\frac{Om\_m}{l\_m}\right)}^{-2} \cdot {\sin ky}^{2}, 1\right)}} + 0.5}\\
                      
                      \mathbf{else}:\\
                      \;\;\;\;\sqrt{0.5}\\
                      
                      
                      \end{array}
                      \end{array}
                      
                      Derivation
                      1. Split input into 2 regimes
                      2. if (/.f64 (*.f64 #s(literal 2 binary64) l) Om) < 1e133

                        1. Initial program 99.5%

                          \[\sqrt{\frac{1}{2} \cdot \left(1 + \frac{1}{\sqrt{1 + {\left(\frac{2 \cdot \ell}{Om}\right)}^{2} \cdot \left({\sin kx}^{2} + {\sin ky}^{2}\right)}}\right)} \]
                        2. Add Preprocessing
                        3. Taylor expanded in kx around 0

                          \[\leadsto \sqrt{\color{blue}{\frac{1}{2} \cdot \left(1 + \sqrt{\frac{1}{1 + 4 \cdot \frac{{\ell}^{2} \cdot {\sin ky}^{2}}{{Om}^{2}}}}\right)}} \]
                        4. Step-by-step derivation
                          1. +-commutativeN/A

                            \[\leadsto \sqrt{\frac{1}{2} \cdot \color{blue}{\left(\sqrt{\frac{1}{1 + 4 \cdot \frac{{\ell}^{2} \cdot {\sin ky}^{2}}{{Om}^{2}}}} + 1\right)}} \]
                          2. distribute-rgt-inN/A

                            \[\leadsto \sqrt{\color{blue}{\sqrt{\frac{1}{1 + 4 \cdot \frac{{\ell}^{2} \cdot {\sin ky}^{2}}{{Om}^{2}}}} \cdot \frac{1}{2} + 1 \cdot \frac{1}{2}}} \]
                          3. metadata-evalN/A

                            \[\leadsto \sqrt{\sqrt{\frac{1}{1 + 4 \cdot \frac{{\ell}^{2} \cdot {\sin ky}^{2}}{{Om}^{2}}}} \cdot \frac{1}{2} + \color{blue}{\frac{1}{2}}} \]
                          4. lower-fma.f64N/A

                            \[\leadsto \sqrt{\color{blue}{\mathsf{fma}\left(\sqrt{\frac{1}{1 + 4 \cdot \frac{{\ell}^{2} \cdot {\sin ky}^{2}}{{Om}^{2}}}}, \frac{1}{2}, \frac{1}{2}\right)}} \]
                        5. Applied rewrites84.8%

                          \[\leadsto \sqrt{\color{blue}{\mathsf{fma}\left(\sqrt{\frac{1}{\mathsf{fma}\left(\frac{{\sin ky}^{2}}{Om} \cdot \frac{\ell \cdot \ell}{Om}, 4, 1\right)}}, 0.5, 0.5\right)}} \]
                        6. Step-by-step derivation
                          1. Applied rewrites93.0%

                            \[\leadsto \sqrt{\frac{0.5}{\sqrt{\mathsf{fma}\left(4, {\sin ky}^{2} \cdot {\left(\frac{Om}{\ell}\right)}^{-2}, 1\right)}} + \color{blue}{0.5}} \]

                          if 1e133 < (/.f64 (*.f64 #s(literal 2 binary64) l) Om)

                          1. Initial program 91.3%

                            \[\sqrt{\frac{1}{2} \cdot \left(1 + \frac{1}{\sqrt{1 + {\left(\frac{2 \cdot \ell}{Om}\right)}^{2} \cdot \left({\sin kx}^{2} + {\sin ky}^{2}\right)}}\right)} \]
                          2. Add Preprocessing
                          3. Taylor expanded in Om around 0

                            \[\leadsto \sqrt{\color{blue}{\frac{1}{2}}} \]
                          4. Step-by-step derivation
                            1. Applied rewrites98.3%

                              \[\leadsto \sqrt{\color{blue}{0.5}} \]
                          5. Recombined 2 regimes into one program.
                          6. Final simplification93.9%

                            \[\leadsto \begin{array}{l} \mathbf{if}\;\frac{\ell \cdot 2}{Om} \leq 10^{+133}:\\ \;\;\;\;\sqrt{\frac{0.5}{\sqrt{\mathsf{fma}\left(4, {\left(\frac{Om}{\ell}\right)}^{-2} \cdot {\sin ky}^{2}, 1\right)}} + 0.5}\\ \mathbf{else}:\\ \;\;\;\;\sqrt{0.5}\\ \end{array} \]
                          7. Add Preprocessing

                          Alternative 7: 55.2% accurate, 52.8× speedup?

                          \[\begin{array}{l} Om_m = \left|Om\right| \\ l_m = \left|\ell\right| \\ \sqrt{0.5} \end{array} \]
                          Om_m = (fabs.f64 Om)
                          l_m = (fabs.f64 l)
                          (FPCore (l_m Om_m kx ky) :precision binary64 (sqrt 0.5))
                          Om_m = fabs(Om);
                          l_m = fabs(l);
                          double code(double l_m, double Om_m, double kx, double ky) {
                          	return sqrt(0.5);
                          }
                          
                          Om_m = abs(om)
                          l_m = abs(l)
                          real(8) function code(l_m, om_m, kx, ky)
                              real(8), intent (in) :: l_m
                              real(8), intent (in) :: om_m
                              real(8), intent (in) :: kx
                              real(8), intent (in) :: ky
                              code = sqrt(0.5d0)
                          end function
                          
                          Om_m = Math.abs(Om);
                          l_m = Math.abs(l);
                          public static double code(double l_m, double Om_m, double kx, double ky) {
                          	return Math.sqrt(0.5);
                          }
                          
                          Om_m = math.fabs(Om)
                          l_m = math.fabs(l)
                          def code(l_m, Om_m, kx, ky):
                          	return math.sqrt(0.5)
                          
                          Om_m = abs(Om)
                          l_m = abs(l)
                          function code(l_m, Om_m, kx, ky)
                          	return sqrt(0.5)
                          end
                          
                          Om_m = abs(Om);
                          l_m = abs(l);
                          function tmp = code(l_m, Om_m, kx, ky)
                          	tmp = sqrt(0.5);
                          end
                          
                          Om_m = N[Abs[Om], $MachinePrecision]
                          l_m = N[Abs[l], $MachinePrecision]
                          code[l$95$m_, Om$95$m_, kx_, ky_] := N[Sqrt[0.5], $MachinePrecision]
                          
                          \begin{array}{l}
                          Om_m = \left|Om\right|
                          \\
                          l_m = \left|\ell\right|
                          
                          \\
                          \sqrt{0.5}
                          \end{array}
                          
                          Derivation
                          1. Initial program 98.0%

                            \[\sqrt{\frac{1}{2} \cdot \left(1 + \frac{1}{\sqrt{1 + {\left(\frac{2 \cdot \ell}{Om}\right)}^{2} \cdot \left({\sin kx}^{2} + {\sin ky}^{2}\right)}}\right)} \]
                          2. Add Preprocessing
                          3. Taylor expanded in Om around 0

                            \[\leadsto \sqrt{\color{blue}{\frac{1}{2}}} \]
                          4. Step-by-step derivation
                            1. Applied rewrites56.8%

                              \[\leadsto \sqrt{\color{blue}{0.5}} \]
                            2. Add Preprocessing

                            Reproduce

                            ?
                            herbie shell --seed 2024243 
                            (FPCore (l Om kx ky)
                              :name "Toniolo and Linder, Equation (3a)"
                              :precision binary64
                              (sqrt (* (/ 1.0 2.0) (+ 1.0 (/ 1.0 (sqrt (+ 1.0 (* (pow (/ (* 2.0 l) Om) 2.0) (+ (pow (sin kx) 2.0) (pow (sin ky) 2.0))))))))))