
(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));
}
real(8) function code(a1, a2, th)
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:
Herbie found 11 alternatives:
| Alternative | Accuracy | Speedup |
|---|
(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));
}
real(8) function code(a1, a2, th)
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}
NOTE: a1, a2, and th should be sorted in increasing order before calling this function. (FPCore (a1 a2 th) :precision binary64 (fma (/ (* (cos th) a2) (sqrt 2.0)) a2 (/ (* a1 (* (cos th) a1)) (sqrt 2.0))))
assert(a1 < a2 && a2 < th);
double code(double a1, double a2, double th) {
return fma(((cos(th) * a2) / sqrt(2.0)), a2, ((a1 * (cos(th) * a1)) / sqrt(2.0)));
}
a1, a2, th = sort([a1, a2, th]) function code(a1, a2, th) return fma(Float64(Float64(cos(th) * a2) / sqrt(2.0)), a2, Float64(Float64(a1 * Float64(cos(th) * a1)) / sqrt(2.0))) end
NOTE: a1, a2, and th should be sorted in increasing order before calling this function. code[a1_, a2_, th_] := N[(N[(N[(N[Cos[th], $MachinePrecision] * a2), $MachinePrecision] / N[Sqrt[2.0], $MachinePrecision]), $MachinePrecision] * a2 + N[(N[(a1 * N[(N[Cos[th], $MachinePrecision] * a1), $MachinePrecision]), $MachinePrecision] / N[Sqrt[2.0], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
[a1, a2, th] = \mathsf{sort}([a1, a2, th])\\
\\
\mathsf{fma}\left(\frac{\cos th \cdot a2}{\sqrt{2}}, a2, \frac{a1 \cdot \left(\cos th \cdot a1\right)}{\sqrt{2}}\right)
\end{array}
Initial program 99.2%
lift-+.f64N/A
+-commutativeN/A
lift-*.f64N/A
lift-*.f64N/A
associate-*r*N/A
lower-fma.f64N/A
lift-/.f64N/A
div-invN/A
associate-*l*N/A
lower-*.f64N/A
associate-*l/N/A
*-lft-identityN/A
lower-/.f6499.5
lift-*.f64N/A
lift-/.f64N/A
associate-*l/N/A
associate-/l*N/A
lower-*.f64N/A
Applied rewrites99.5%
lift-*.f64N/A
lift-/.f64N/A
associate-*r/N/A
lower-/.f64N/A
lower-*.f6499.5
lift-*.f64N/A
lift-*.f64N/A
lift-/.f64N/A
associate-*r/N/A
lift-*.f64N/A
associate-*r/N/A
lower-/.f64N/A
lift-*.f64N/A
associate-*r*N/A
lower-*.f64N/A
lower-*.f6499.5
Applied rewrites99.5%
Final simplification99.5%
NOTE: a1, a2, and th should be sorted in increasing order before calling this function.
(FPCore (a1 a2 th)
:precision binary64
(let* ((t_1 (/ (cos th) (sqrt 2.0))))
(if (<= (+ (* t_1 (* a1 a1)) (* (* a2 a2) t_1)) -4e-88)
(* (fma (* th th) -0.5 1.0) (/ (* a2 a2) (sqrt 2.0)))
(/ 1.0 (/ (sqrt 2.0) (fma a1 a1 (* a2 a2)))))))assert(a1 < a2 && a2 < th);
double code(double a1, double a2, double th) {
double t_1 = cos(th) / sqrt(2.0);
double tmp;
if (((t_1 * (a1 * a1)) + ((a2 * a2) * t_1)) <= -4e-88) {
tmp = fma((th * th), -0.5, 1.0) * ((a2 * a2) / sqrt(2.0));
} else {
tmp = 1.0 / (sqrt(2.0) / fma(a1, a1, (a2 * a2)));
}
return tmp;
}
a1, a2, th = sort([a1, a2, th]) function code(a1, a2, th) t_1 = Float64(cos(th) / sqrt(2.0)) tmp = 0.0 if (Float64(Float64(t_1 * Float64(a1 * a1)) + Float64(Float64(a2 * a2) * t_1)) <= -4e-88) tmp = Float64(fma(Float64(th * th), -0.5, 1.0) * Float64(Float64(a2 * a2) / sqrt(2.0))); else tmp = Float64(1.0 / Float64(sqrt(2.0) / fma(a1, a1, Float64(a2 * a2)))); end return tmp end
NOTE: a1, a2, and th should be sorted in increasing order before calling this function.
code[a1_, a2_, th_] := Block[{t$95$1 = N[(N[Cos[th], $MachinePrecision] / N[Sqrt[2.0], $MachinePrecision]), $MachinePrecision]}, If[LessEqual[N[(N[(t$95$1 * N[(a1 * a1), $MachinePrecision]), $MachinePrecision] + N[(N[(a2 * a2), $MachinePrecision] * t$95$1), $MachinePrecision]), $MachinePrecision], -4e-88], N[(N[(N[(th * th), $MachinePrecision] * -0.5 + 1.0), $MachinePrecision] * N[(N[(a2 * a2), $MachinePrecision] / N[Sqrt[2.0], $MachinePrecision]), $MachinePrecision]), $MachinePrecision], N[(1.0 / N[(N[Sqrt[2.0], $MachinePrecision] / N[(a1 * a1 + N[(a2 * a2), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]]]
\begin{array}{l}
[a1, a2, th] = \mathsf{sort}([a1, a2, th])\\
\\
\begin{array}{l}
t_1 := \frac{\cos th}{\sqrt{2}}\\
\mathbf{if}\;t\_1 \cdot \left(a1 \cdot a1\right) + \left(a2 \cdot a2\right) \cdot t\_1 \leq -4 \cdot 10^{-88}:\\
\;\;\;\;\mathsf{fma}\left(th \cdot th, -0.5, 1\right) \cdot \frac{a2 \cdot a2}{\sqrt{2}}\\
\mathbf{else}:\\
\;\;\;\;\frac{1}{\frac{\sqrt{2}}{\mathsf{fma}\left(a1, a1, a2 \cdot a2\right)}}\\
\end{array}
\end{array}
if (+.f64 (*.f64 (/.f64 (cos.f64 th) (sqrt.f64 #s(literal 2 binary64))) (*.f64 a1 a1)) (*.f64 (/.f64 (cos.f64 th) (sqrt.f64 #s(literal 2 binary64))) (*.f64 a2 a2))) < -3.99999999999999974e-88Initial program 99.6%
lift-+.f64N/A
+-commutativeN/A
lift-*.f64N/A
lift-*.f64N/A
associate-*r*N/A
lower-fma.f64N/A
lift-/.f64N/A
div-invN/A
associate-*l*N/A
lower-*.f64N/A
associate-*l/N/A
*-lft-identityN/A
lower-/.f6499.6
lift-*.f64N/A
lift-/.f64N/A
associate-*l/N/A
associate-/l*N/A
lower-*.f64N/A
Applied rewrites99.6%
lift-*.f64N/A
lift-/.f64N/A
associate-*r/N/A
lower-/.f64N/A
lower-*.f6499.6
lift-*.f64N/A
lift-*.f64N/A
lift-/.f64N/A
associate-*r/N/A
lift-*.f64N/A
associate-*r/N/A
lower-/.f64N/A
lift-*.f64N/A
associate-*r*N/A
lower-*.f64N/A
lower-*.f6499.6
Applied rewrites99.6%
Taylor expanded in th around 0
distribute-lft-outN/A
associate-*r*N/A
*-commutativeN/A
distribute-lft1-inN/A
+-commutativeN/A
lower-*.f64N/A
+-commutativeN/A
*-commutativeN/A
lower-fma.f64N/A
unpow2N/A
lower-*.f64N/A
Applied rewrites55.9%
Taylor expanded in a1 around 0
Applied rewrites42.2%
if -3.99999999999999974e-88 < (+.f64 (*.f64 (/.f64 (cos.f64 th) (sqrt.f64 #s(literal 2 binary64))) (*.f64 a1 a1)) (*.f64 (/.f64 (cos.f64 th) (sqrt.f64 #s(literal 2 binary64))) (*.f64 a2 a2))) Initial program 99.6%
lift-+.f64N/A
lift-*.f64N/A
lift-*.f64N/A
distribute-lft-outN/A
lift-/.f64N/A
clear-numN/A
associate-*l/N/A
clear-numN/A
lower-/.f64N/A
*-lft-identityN/A
lower-/.f64N/A
lower-/.f64N/A
lift-*.f64N/A
lower-fma.f6499.3
Applied rewrites99.3%
Taylor expanded in th around 0
lower-sqrt.f6483.1
Applied rewrites83.1%
Final simplification74.7%
herbie shell --seed 2024222
(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))))