
(FPCore (lambda1 phi1 phi2 delta theta)
:precision binary64
(+
lambda1
(atan2
(* (* (sin theta) (sin delta)) (cos phi1))
(-
(cos delta)
(*
(sin phi1)
(sin
(asin
(+
(* (sin phi1) (cos delta))
(* (* (cos phi1) (sin delta)) (cos theta))))))))))
double code(double lambda1, double phi1, double phi2, double delta, double theta) {
return lambda1 + atan2(((sin(theta) * sin(delta)) * cos(phi1)), (cos(delta) - (sin(phi1) * sin(asin(((sin(phi1) * cos(delta)) + ((cos(phi1) * sin(delta)) * cos(theta))))))));
}
real(8) function code(lambda1, phi1, phi2, delta, theta)
real(8), intent (in) :: lambda1
real(8), intent (in) :: phi1
real(8), intent (in) :: phi2
real(8), intent (in) :: delta
real(8), intent (in) :: theta
code = lambda1 + atan2(((sin(theta) * sin(delta)) * cos(phi1)), (cos(delta) - (sin(phi1) * sin(asin(((sin(phi1) * cos(delta)) + ((cos(phi1) * sin(delta)) * cos(theta))))))))
end function
public static double code(double lambda1, double phi1, double phi2, double delta, double theta) {
return lambda1 + Math.atan2(((Math.sin(theta) * Math.sin(delta)) * Math.cos(phi1)), (Math.cos(delta) - (Math.sin(phi1) * Math.sin(Math.asin(((Math.sin(phi1) * Math.cos(delta)) + ((Math.cos(phi1) * Math.sin(delta)) * Math.cos(theta))))))));
}
def code(lambda1, phi1, phi2, delta, theta): return lambda1 + math.atan2(((math.sin(theta) * math.sin(delta)) * math.cos(phi1)), (math.cos(delta) - (math.sin(phi1) * math.sin(math.asin(((math.sin(phi1) * math.cos(delta)) + ((math.cos(phi1) * math.sin(delta)) * math.cos(theta))))))))
function code(lambda1, phi1, phi2, delta, theta) return Float64(lambda1 + atan(Float64(Float64(sin(theta) * sin(delta)) * cos(phi1)), Float64(cos(delta) - Float64(sin(phi1) * sin(asin(Float64(Float64(sin(phi1) * cos(delta)) + Float64(Float64(cos(phi1) * sin(delta)) * cos(theta))))))))) end
function tmp = code(lambda1, phi1, phi2, delta, theta) tmp = lambda1 + atan2(((sin(theta) * sin(delta)) * cos(phi1)), (cos(delta) - (sin(phi1) * sin(asin(((sin(phi1) * cos(delta)) + ((cos(phi1) * sin(delta)) * cos(theta)))))))); end
code[lambda1_, phi1_, phi2_, delta_, theta_] := N[(lambda1 + N[ArcTan[N[(N[(N[Sin[theta], $MachinePrecision] * N[Sin[delta], $MachinePrecision]), $MachinePrecision] * N[Cos[phi1], $MachinePrecision]), $MachinePrecision] / N[(N[Cos[delta], $MachinePrecision] - N[(N[Sin[phi1], $MachinePrecision] * N[Sin[N[ArcSin[N[(N[(N[Sin[phi1], $MachinePrecision] * N[Cos[delta], $MachinePrecision]), $MachinePrecision] + N[(N[(N[Cos[phi1], $MachinePrecision] * N[Sin[delta], $MachinePrecision]), $MachinePrecision] * N[Cos[theta], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
\lambda_1 + \tan^{-1}_* \frac{\left(\sin theta \cdot \sin delta\right) \cdot \cos \phi_1}{\cos delta - \sin \phi_1 \cdot \sin \sin^{-1} \left(\sin \phi_1 \cdot \cos delta + \left(\cos \phi_1 \cdot \sin delta\right) \cdot \cos theta\right)}
\end{array}
Sampling outcomes in binary64 precision:
Herbie found 12 alternatives:
| Alternative | Accuracy | Speedup |
|---|
(FPCore (lambda1 phi1 phi2 delta theta)
:precision binary64
(+
lambda1
(atan2
(* (* (sin theta) (sin delta)) (cos phi1))
(-
(cos delta)
(*
(sin phi1)
(sin
(asin
(+
(* (sin phi1) (cos delta))
(* (* (cos phi1) (sin delta)) (cos theta))))))))))
double code(double lambda1, double phi1, double phi2, double delta, double theta) {
return lambda1 + atan2(((sin(theta) * sin(delta)) * cos(phi1)), (cos(delta) - (sin(phi1) * sin(asin(((sin(phi1) * cos(delta)) + ((cos(phi1) * sin(delta)) * cos(theta))))))));
}
real(8) function code(lambda1, phi1, phi2, delta, theta)
real(8), intent (in) :: lambda1
real(8), intent (in) :: phi1
real(8), intent (in) :: phi2
real(8), intent (in) :: delta
real(8), intent (in) :: theta
code = lambda1 + atan2(((sin(theta) * sin(delta)) * cos(phi1)), (cos(delta) - (sin(phi1) * sin(asin(((sin(phi1) * cos(delta)) + ((cos(phi1) * sin(delta)) * cos(theta))))))))
end function
public static double code(double lambda1, double phi1, double phi2, double delta, double theta) {
return lambda1 + Math.atan2(((Math.sin(theta) * Math.sin(delta)) * Math.cos(phi1)), (Math.cos(delta) - (Math.sin(phi1) * Math.sin(Math.asin(((Math.sin(phi1) * Math.cos(delta)) + ((Math.cos(phi1) * Math.sin(delta)) * Math.cos(theta))))))));
}
def code(lambda1, phi1, phi2, delta, theta): return lambda1 + math.atan2(((math.sin(theta) * math.sin(delta)) * math.cos(phi1)), (math.cos(delta) - (math.sin(phi1) * math.sin(math.asin(((math.sin(phi1) * math.cos(delta)) + ((math.cos(phi1) * math.sin(delta)) * math.cos(theta))))))))
function code(lambda1, phi1, phi2, delta, theta) return Float64(lambda1 + atan(Float64(Float64(sin(theta) * sin(delta)) * cos(phi1)), Float64(cos(delta) - Float64(sin(phi1) * sin(asin(Float64(Float64(sin(phi1) * cos(delta)) + Float64(Float64(cos(phi1) * sin(delta)) * cos(theta))))))))) end
function tmp = code(lambda1, phi1, phi2, delta, theta) tmp = lambda1 + atan2(((sin(theta) * sin(delta)) * cos(phi1)), (cos(delta) - (sin(phi1) * sin(asin(((sin(phi1) * cos(delta)) + ((cos(phi1) * sin(delta)) * cos(theta)))))))); end
code[lambda1_, phi1_, phi2_, delta_, theta_] := N[(lambda1 + N[ArcTan[N[(N[(N[Sin[theta], $MachinePrecision] * N[Sin[delta], $MachinePrecision]), $MachinePrecision] * N[Cos[phi1], $MachinePrecision]), $MachinePrecision] / N[(N[Cos[delta], $MachinePrecision] - N[(N[Sin[phi1], $MachinePrecision] * N[Sin[N[ArcSin[N[(N[(N[Sin[phi1], $MachinePrecision] * N[Cos[delta], $MachinePrecision]), $MachinePrecision] + N[(N[(N[Cos[phi1], $MachinePrecision] * N[Sin[delta], $MachinePrecision]), $MachinePrecision] * N[Cos[theta], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
\lambda_1 + \tan^{-1}_* \frac{\left(\sin theta \cdot \sin delta\right) \cdot \cos \phi_1}{\cos delta - \sin \phi_1 \cdot \sin \sin^{-1} \left(\sin \phi_1 \cdot \cos delta + \left(\cos \phi_1 \cdot \sin delta\right) \cdot \cos theta\right)}
\end{array}
(FPCore (lambda1 phi1 phi2 delta theta)
:precision binary64
(+
lambda1
(atan2
(* (sin theta) (* (cos phi1) (sin delta)))
(fma
(pow (cos phi1) 2.0)
(cos delta)
(* (* (- (sin phi1)) (* (cos theta) (cos phi1))) (sin delta))))))
double code(double lambda1, double phi1, double phi2, double delta, double theta) {
return lambda1 + atan2((sin(theta) * (cos(phi1) * sin(delta))), fma(pow(cos(phi1), 2.0), cos(delta), ((-sin(phi1) * (cos(theta) * cos(phi1))) * sin(delta))));
}
function code(lambda1, phi1, phi2, delta, theta) return Float64(lambda1 + atan(Float64(sin(theta) * Float64(cos(phi1) * sin(delta))), fma((cos(phi1) ^ 2.0), cos(delta), Float64(Float64(Float64(-sin(phi1)) * Float64(cos(theta) * cos(phi1))) * sin(delta))))) end
code[lambda1_, phi1_, phi2_, delta_, theta_] := N[(lambda1 + N[ArcTan[N[(N[Sin[theta], $MachinePrecision] * N[(N[Cos[phi1], $MachinePrecision] * N[Sin[delta], $MachinePrecision]), $MachinePrecision]), $MachinePrecision] / N[(N[Power[N[Cos[phi1], $MachinePrecision], 2.0], $MachinePrecision] * N[Cos[delta], $MachinePrecision] + N[(N[((-N[Sin[phi1], $MachinePrecision]) * N[(N[Cos[theta], $MachinePrecision] * N[Cos[phi1], $MachinePrecision]), $MachinePrecision]), $MachinePrecision] * N[Sin[delta], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]
\begin{array}{l}
\\
\lambda_1 + \tan^{-1}_* \frac{\sin theta \cdot \left(\cos \phi_1 \cdot \sin delta\right)}{\mathsf{fma}\left({\cos \phi_1}^{2}, \cos delta, \left(\left(-\sin \phi_1\right) \cdot \left(\cos theta \cdot \cos \phi_1\right)\right) \cdot \sin delta\right)}
\end{array}
Initial program 99.7%
Applied rewrites99.6%
lift-+.f64N/A
lift-neg.f64N/A
lift-*.f64N/A
distribute-lft-neg-inN/A
distribute-lft1-inN/A
lower-*.f64N/A
Applied rewrites99.6%
Applied rewrites99.7%
lift-fma.f64N/A
+-commutativeN/A
lift-cos.f64N/A
lift-+.f64N/A
count-2N/A
sqr-cos-aN/A
lift-cos.f64N/A
lift-cos.f64N/A
pow2N/A
lower-pow.f6499.8
Applied rewrites99.8%
Final simplification99.8%
(FPCore (lambda1 phi1 phi2 delta theta)
:precision binary64
(let* ((t_1 (* (cos phi1) (sin delta)))
(t_2 (* (sin theta) t_1))
(t_3
(+
(atan2
(* (* (sin theta) (sin delta)) (cos phi1))
(-
(cos delta)
(*
(sin (asin (+ (* (cos theta) t_1) (* (sin phi1) (cos delta)))))
(sin phi1))))
lambda1)))
(if (<= t_3 5e-15)
(+
(atan2
t_2
(fma
(fma 0.5 (cos (+ phi1 phi1)) 0.5)
(cos delta)
(* (* (- (sin phi1)) (cos phi1)) (sin delta))))
lambda1)
(if (<= t_3 5.0)
(atan2
t_2
(fma
(fma (cos (* 2.0 phi1)) 0.5 0.5)
(cos delta)
(* (* (* (sin phi1) (sin delta)) (cos phi1)) (- (cos theta)))))
(+ (atan2 (* theta (sin delta)) (cos delta)) lambda1)))))
double code(double lambda1, double phi1, double phi2, double delta, double theta) {
double t_1 = cos(phi1) * sin(delta);
double t_2 = sin(theta) * t_1;
double t_3 = atan2(((sin(theta) * sin(delta)) * cos(phi1)), (cos(delta) - (sin(asin(((cos(theta) * t_1) + (sin(phi1) * cos(delta))))) * sin(phi1)))) + lambda1;
double tmp;
if (t_3 <= 5e-15) {
tmp = atan2(t_2, fma(fma(0.5, cos((phi1 + phi1)), 0.5), cos(delta), ((-sin(phi1) * cos(phi1)) * sin(delta)))) + lambda1;
} else if (t_3 <= 5.0) {
tmp = atan2(t_2, fma(fma(cos((2.0 * phi1)), 0.5, 0.5), cos(delta), (((sin(phi1) * sin(delta)) * cos(phi1)) * -cos(theta))));
} else {
tmp = atan2((theta * sin(delta)), cos(delta)) + lambda1;
}
return tmp;
}
function code(lambda1, phi1, phi2, delta, theta) t_1 = Float64(cos(phi1) * sin(delta)) t_2 = Float64(sin(theta) * t_1) t_3 = Float64(atan(Float64(Float64(sin(theta) * sin(delta)) * cos(phi1)), Float64(cos(delta) - Float64(sin(asin(Float64(Float64(cos(theta) * t_1) + Float64(sin(phi1) * cos(delta))))) * sin(phi1)))) + lambda1) tmp = 0.0 if (t_3 <= 5e-15) tmp = Float64(atan(t_2, fma(fma(0.5, cos(Float64(phi1 + phi1)), 0.5), cos(delta), Float64(Float64(Float64(-sin(phi1)) * cos(phi1)) * sin(delta)))) + lambda1); elseif (t_3 <= 5.0) tmp = atan(t_2, fma(fma(cos(Float64(2.0 * phi1)), 0.5, 0.5), cos(delta), Float64(Float64(Float64(sin(phi1) * sin(delta)) * cos(phi1)) * Float64(-cos(theta))))); else tmp = Float64(atan(Float64(theta * sin(delta)), cos(delta)) + lambda1); end return tmp end
code[lambda1_, phi1_, phi2_, delta_, theta_] := Block[{t$95$1 = N[(N[Cos[phi1], $MachinePrecision] * N[Sin[delta], $MachinePrecision]), $MachinePrecision]}, Block[{t$95$2 = N[(N[Sin[theta], $MachinePrecision] * t$95$1), $MachinePrecision]}, Block[{t$95$3 = N[(N[ArcTan[N[(N[(N[Sin[theta], $MachinePrecision] * N[Sin[delta], $MachinePrecision]), $MachinePrecision] * N[Cos[phi1], $MachinePrecision]), $MachinePrecision] / N[(N[Cos[delta], $MachinePrecision] - N[(N[Sin[N[ArcSin[N[(N[(N[Cos[theta], $MachinePrecision] * t$95$1), $MachinePrecision] + N[(N[Sin[phi1], $MachinePrecision] * N[Cos[delta], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]], $MachinePrecision] * N[Sin[phi1], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision] + lambda1), $MachinePrecision]}, If[LessEqual[t$95$3, 5e-15], N[(N[ArcTan[t$95$2 / N[(N[(0.5 * N[Cos[N[(phi1 + phi1), $MachinePrecision]], $MachinePrecision] + 0.5), $MachinePrecision] * N[Cos[delta], $MachinePrecision] + N[(N[((-N[Sin[phi1], $MachinePrecision]) * N[Cos[phi1], $MachinePrecision]), $MachinePrecision] * N[Sin[delta], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision] + lambda1), $MachinePrecision], If[LessEqual[t$95$3, 5.0], N[ArcTan[t$95$2 / N[(N[(N[Cos[N[(2.0 * phi1), $MachinePrecision]], $MachinePrecision] * 0.5 + 0.5), $MachinePrecision] * N[Cos[delta], $MachinePrecision] + N[(N[(N[(N[Sin[phi1], $MachinePrecision] * N[Sin[delta], $MachinePrecision]), $MachinePrecision] * N[Cos[phi1], $MachinePrecision]), $MachinePrecision] * (-N[Cos[theta], $MachinePrecision])), $MachinePrecision]), $MachinePrecision]], $MachinePrecision], N[(N[ArcTan[N[(theta * N[Sin[delta], $MachinePrecision]), $MachinePrecision] / N[Cos[delta], $MachinePrecision]], $MachinePrecision] + lambda1), $MachinePrecision]]]]]]
\begin{array}{l}
\\
\begin{array}{l}
t_1 := \cos \phi_1 \cdot \sin delta\\
t_2 := \sin theta \cdot t\_1\\
t_3 := \tan^{-1}_* \frac{\left(\sin theta \cdot \sin delta\right) \cdot \cos \phi_1}{\cos delta - \sin \sin^{-1} \left(\cos theta \cdot t\_1 + \sin \phi_1 \cdot \cos delta\right) \cdot \sin \phi_1} + \lambda_1\\
\mathbf{if}\;t\_3 \leq 5 \cdot 10^{-15}:\\
\;\;\;\;\tan^{-1}_* \frac{t\_2}{\mathsf{fma}\left(\mathsf{fma}\left(0.5, \cos \left(\phi_1 + \phi_1\right), 0.5\right), \cos delta, \left(\left(-\sin \phi_1\right) \cdot \cos \phi_1\right) \cdot \sin delta\right)} + \lambda_1\\
\mathbf{elif}\;t\_3 \leq 5:\\
\;\;\;\;\tan^{-1}_* \frac{t\_2}{\mathsf{fma}\left(\mathsf{fma}\left(\cos \left(2 \cdot \phi_1\right), 0.5, 0.5\right), \cos delta, \left(\left(\sin \phi_1 \cdot \sin delta\right) \cdot \cos \phi_1\right) \cdot \left(-\cos theta\right)\right)}\\
\mathbf{else}:\\
\;\;\;\;\tan^{-1}_* \frac{theta \cdot \sin delta}{\cos delta} + \lambda_1\\
\end{array}
\end{array}
if (+.f64 lambda1 (atan2.f64 (*.f64 (*.f64 (sin.f64 theta) (sin.f64 delta)) (cos.f64 phi1)) (-.f64 (cos.f64 delta) (*.f64 (sin.f64 phi1) (sin.f64 (asin.f64 (+.f64 (*.f64 (sin.f64 phi1) (cos.f64 delta)) (*.f64 (*.f64 (cos.f64 phi1) (sin.f64 delta)) (cos.f64 theta))))))))) < 4.99999999999999999e-15Initial program 99.7%
Applied rewrites99.7%
lift-+.f64N/A
lift-neg.f64N/A
lift-*.f64N/A
distribute-lft-neg-inN/A
distribute-lft1-inN/A
lower-*.f64N/A
Applied rewrites99.7%
Applied rewrites99.8%
Taylor expanded in theta around 0
associate-*r*N/A
lower-*.f64N/A
neg-mul-1N/A
lower-neg.f64N/A
lower-cos.f64N/A
lower-sin.f6495.9
Applied rewrites95.9%
if 4.99999999999999999e-15 < (+.f64 lambda1 (atan2.f64 (*.f64 (*.f64 (sin.f64 theta) (sin.f64 delta)) (cos.f64 phi1)) (-.f64 (cos.f64 delta) (*.f64 (sin.f64 phi1) (sin.f64 (asin.f64 (+.f64 (*.f64 (sin.f64 phi1) (cos.f64 delta)) (*.f64 (*.f64 (cos.f64 phi1) (sin.f64 delta)) (cos.f64 theta))))))))) < 5Initial program 99.6%
Applied rewrites99.6%
lift-+.f64N/A
lift-neg.f64N/A
lift-*.f64N/A
distribute-lft-neg-inN/A
distribute-lft1-inN/A
lower-*.f64N/A
Applied rewrites99.6%
Taylor expanded in lambda1 around 0
Applied rewrites90.6%
if 5 < (+.f64 lambda1 (atan2.f64 (*.f64 (*.f64 (sin.f64 theta) (sin.f64 delta)) (cos.f64 phi1)) (-.f64 (cos.f64 delta) (*.f64 (sin.f64 phi1) (sin.f64 (asin.f64 (+.f64 (*.f64 (sin.f64 phi1) (cos.f64 delta)) (*.f64 (*.f64 (cos.f64 phi1) (sin.f64 delta)) (cos.f64 theta))))))))) Initial program 100.0%
Taylor expanded in phi1 around 0
lower-cos.f6499.7
Applied rewrites99.7%
Taylor expanded in phi1 around 0
*-commutativeN/A
lower-*.f64N/A
lower-sin.f64N/A
lower-sin.f6499.7
Applied rewrites99.7%
Taylor expanded in theta around 0
Applied rewrites99.0%
Final simplification96.1%
herbie shell --seed 2024234
(FPCore (lambda1 phi1 phi2 delta theta)
:name "Destination given bearing on a great circle"
:precision binary64
(+ lambda1 (atan2 (* (* (sin theta) (sin delta)) (cos phi1)) (- (cos delta) (* (sin phi1) (sin (asin (+ (* (sin phi1) (cos delta)) (* (* (cos phi1) (sin delta)) (cos theta))))))))))