
(FPCore (lambda1 lambda2 phi1 phi2) :precision binary64 (atan2 (* (sin (- lambda1 lambda2)) (cos phi2)) (- (* (cos phi1) (sin phi2)) (* (* (sin phi1) (cos phi2)) (cos (- lambda1 lambda2))))))
double code(double lambda1, double lambda2, double phi1, double phi2) {
return atan2((sin((lambda1 - lambda2)) * cos(phi2)), ((cos(phi1) * sin(phi2)) - ((sin(phi1) * cos(phi2)) * cos((lambda1 - lambda2)))));
}
real(8) function code(lambda1, lambda2, phi1, phi2)
real(8), intent (in) :: lambda1
real(8), intent (in) :: lambda2
real(8), intent (in) :: phi1
real(8), intent (in) :: phi2
code = atan2((sin((lambda1 - lambda2)) * cos(phi2)), ((cos(phi1) * sin(phi2)) - ((sin(phi1) * cos(phi2)) * cos((lambda1 - lambda2)))))
end function
public static double code(double lambda1, double lambda2, double phi1, double phi2) {
return Math.atan2((Math.sin((lambda1 - lambda2)) * Math.cos(phi2)), ((Math.cos(phi1) * Math.sin(phi2)) - ((Math.sin(phi1) * Math.cos(phi2)) * Math.cos((lambda1 - lambda2)))));
}
def code(lambda1, lambda2, phi1, phi2): return math.atan2((math.sin((lambda1 - lambda2)) * math.cos(phi2)), ((math.cos(phi1) * math.sin(phi2)) - ((math.sin(phi1) * math.cos(phi2)) * math.cos((lambda1 - lambda2)))))
function code(lambda1, lambda2, phi1, phi2) return atan(Float64(sin(Float64(lambda1 - lambda2)) * cos(phi2)), Float64(Float64(cos(phi1) * sin(phi2)) - Float64(Float64(sin(phi1) * cos(phi2)) * cos(Float64(lambda1 - lambda2))))) end
function tmp = code(lambda1, lambda2, phi1, phi2) tmp = atan2((sin((lambda1 - lambda2)) * cos(phi2)), ((cos(phi1) * sin(phi2)) - ((sin(phi1) * cos(phi2)) * cos((lambda1 - lambda2))))); end
code[lambda1_, lambda2_, phi1_, phi2_] := N[ArcTan[N[(N[Sin[N[(lambda1 - lambda2), $MachinePrecision]], $MachinePrecision] * N[Cos[phi2], $MachinePrecision]), $MachinePrecision] / N[(N[(N[Cos[phi1], $MachinePrecision] * N[Sin[phi2], $MachinePrecision]), $MachinePrecision] - N[(N[(N[Sin[phi1], $MachinePrecision] * N[Cos[phi2], $MachinePrecision]), $MachinePrecision] * N[Cos[N[(lambda1 - lambda2), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]
\begin{array}{l}
\\
\tan^{-1}_* \frac{\sin \left(\lambda_1 - \lambda_2\right) \cdot \cos \phi_2}{\cos \phi_1 \cdot \sin \phi_2 - \left(\sin \phi_1 \cdot \cos \phi_2\right) \cdot \cos \left(\lambda_1 - \lambda_2\right)}
\end{array}
Sampling outcomes in binary64 precision:
Herbie found 35 alternatives:
| Alternative | Accuracy | Speedup |
|---|
(FPCore (lambda1 lambda2 phi1 phi2) :precision binary64 (atan2 (* (sin (- lambda1 lambda2)) (cos phi2)) (- (* (cos phi1) (sin phi2)) (* (* (sin phi1) (cos phi2)) (cos (- lambda1 lambda2))))))
double code(double lambda1, double lambda2, double phi1, double phi2) {
return atan2((sin((lambda1 - lambda2)) * cos(phi2)), ((cos(phi1) * sin(phi2)) - ((sin(phi1) * cos(phi2)) * cos((lambda1 - lambda2)))));
}
real(8) function code(lambda1, lambda2, phi1, phi2)
real(8), intent (in) :: lambda1
real(8), intent (in) :: lambda2
real(8), intent (in) :: phi1
real(8), intent (in) :: phi2
code = atan2((sin((lambda1 - lambda2)) * cos(phi2)), ((cos(phi1) * sin(phi2)) - ((sin(phi1) * cos(phi2)) * cos((lambda1 - lambda2)))))
end function
public static double code(double lambda1, double lambda2, double phi1, double phi2) {
return Math.atan2((Math.sin((lambda1 - lambda2)) * Math.cos(phi2)), ((Math.cos(phi1) * Math.sin(phi2)) - ((Math.sin(phi1) * Math.cos(phi2)) * Math.cos((lambda1 - lambda2)))));
}
def code(lambda1, lambda2, phi1, phi2): return math.atan2((math.sin((lambda1 - lambda2)) * math.cos(phi2)), ((math.cos(phi1) * math.sin(phi2)) - ((math.sin(phi1) * math.cos(phi2)) * math.cos((lambda1 - lambda2)))))
function code(lambda1, lambda2, phi1, phi2) return atan(Float64(sin(Float64(lambda1 - lambda2)) * cos(phi2)), Float64(Float64(cos(phi1) * sin(phi2)) - Float64(Float64(sin(phi1) * cos(phi2)) * cos(Float64(lambda1 - lambda2))))) end
function tmp = code(lambda1, lambda2, phi1, phi2) tmp = atan2((sin((lambda1 - lambda2)) * cos(phi2)), ((cos(phi1) * sin(phi2)) - ((sin(phi1) * cos(phi2)) * cos((lambda1 - lambda2))))); end
code[lambda1_, lambda2_, phi1_, phi2_] := N[ArcTan[N[(N[Sin[N[(lambda1 - lambda2), $MachinePrecision]], $MachinePrecision] * N[Cos[phi2], $MachinePrecision]), $MachinePrecision] / N[(N[(N[Cos[phi1], $MachinePrecision] * N[Sin[phi2], $MachinePrecision]), $MachinePrecision] - N[(N[(N[Sin[phi1], $MachinePrecision] * N[Cos[phi2], $MachinePrecision]), $MachinePrecision] * N[Cos[N[(lambda1 - lambda2), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]
\begin{array}{l}
\\
\tan^{-1}_* \frac{\sin \left(\lambda_1 - \lambda_2\right) \cdot \cos \phi_2}{\cos \phi_1 \cdot \sin \phi_2 - \left(\sin \phi_1 \cdot \cos \phi_2\right) \cdot \cos \left(\lambda_1 - \lambda_2\right)}
\end{array}
(FPCore (lambda1 lambda2 phi1 phi2)
:precision binary64
(atan2
(*
(fma (sin lambda1) (cos lambda2) (- (* (cos lambda1) (sin lambda2))))
(cos phi2))
(-
(* (cos phi1) (sin phi2))
(*
(* (cos phi2) (sin phi1))
(fma (sin lambda2) (sin lambda1) (* (cos lambda2) (cos lambda1)))))))
double code(double lambda1, double lambda2, double phi1, double phi2) {
return atan2((fma(sin(lambda1), cos(lambda2), -(cos(lambda1) * sin(lambda2))) * cos(phi2)), ((cos(phi1) * sin(phi2)) - ((cos(phi2) * sin(phi1)) * fma(sin(lambda2), sin(lambda1), (cos(lambda2) * cos(lambda1))))));
}
function code(lambda1, lambda2, phi1, phi2) return atan(Float64(fma(sin(lambda1), cos(lambda2), Float64(-Float64(cos(lambda1) * sin(lambda2)))) * cos(phi2)), Float64(Float64(cos(phi1) * sin(phi2)) - Float64(Float64(cos(phi2) * sin(phi1)) * fma(sin(lambda2), sin(lambda1), Float64(cos(lambda2) * cos(lambda1)))))) end
code[lambda1_, lambda2_, phi1_, phi2_] := N[ArcTan[N[(N[(N[Sin[lambda1], $MachinePrecision] * N[Cos[lambda2], $MachinePrecision] + (-N[(N[Cos[lambda1], $MachinePrecision] * N[Sin[lambda2], $MachinePrecision]), $MachinePrecision])), $MachinePrecision] * N[Cos[phi2], $MachinePrecision]), $MachinePrecision] / N[(N[(N[Cos[phi1], $MachinePrecision] * N[Sin[phi2], $MachinePrecision]), $MachinePrecision] - N[(N[(N[Cos[phi2], $MachinePrecision] * N[Sin[phi1], $MachinePrecision]), $MachinePrecision] * N[(N[Sin[lambda2], $MachinePrecision] * N[Sin[lambda1], $MachinePrecision] + N[(N[Cos[lambda2], $MachinePrecision] * N[Cos[lambda1], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]
\begin{array}{l}
\\
\tan^{-1}_* \frac{\mathsf{fma}\left(\sin \lambda_1, \cos \lambda_2, -\cos \lambda_1 \cdot \sin \lambda_2\right) \cdot \cos \phi_2}{\cos \phi_1 \cdot \sin \phi_2 - \left(\cos \phi_2 \cdot \sin \phi_1\right) \cdot \mathsf{fma}\left(\sin \lambda_2, \sin \lambda_1, \cos \lambda_2 \cdot \cos \lambda_1\right)}
\end{array}
Initial program 79.2%
lift-sin.f64N/A
lift--.f64N/A
sub-negN/A
+-commutativeN/A
sin-sumN/A
cos-negN/A
*-commutativeN/A
lower-fma.f64N/A
lower-sin.f64N/A
lower-neg.f64N/A
lower-cos.f64N/A
lower-*.f64N/A
lower-sin.f64N/A
lower-cos.f6490.0
Applied rewrites90.0%
lift-cos.f64N/A
lift--.f64N/A
cos-diffN/A
+-commutativeN/A
lift-sin.f64N/A
*-commutativeN/A
lower-fma.f64N/A
lower-sin.f64N/A
lift-cos.f64N/A
lift-cos.f64N/A
lower-*.f6499.8
Applied rewrites99.8%
lift-fma.f64N/A
+-commutativeN/A
*-commutativeN/A
lift-sin.f64N/A
lift-neg.f64N/A
sin-negN/A
lift-sin.f64N/A
distribute-rgt-neg-outN/A
lift-cos.f64N/A
lift-sin.f64N/A
sub-negN/A
lift-cos.f64N/A
lift-sin.f64N/A
cancel-sign-sub-invN/A
lift-*.f64N/A
lower-fma.f64N/A
lower-*.f64N/A
lower-neg.f6499.8
Applied rewrites99.8%
Final simplification99.8%
(FPCore (lambda1 lambda2 phi1 phi2)
:precision binary64
(let* ((t_0 (* (cos phi1) (sin phi2)))
(t_1 (* (sin lambda1) (cos lambda2)))
(t_2 (* (cos phi2) (sin phi1)))
(t_3 (sin (- lambda2)))
(t_4 (* (cos phi2) (fma t_3 (cos lambda1) t_1))))
(if (<= phi2 -750.0)
(atan2
(fma t_1 (cos phi2) (* (cos phi2) (* (cos lambda1) t_3)))
(- t_0 (* t_2 (cos (- lambda1 lambda2)))))
(if (<= phi2 9000000000.0)
(atan2
t_4
(-
t_0
(*
(sin phi1)
(fma (cos lambda1) (cos lambda2) (* (sin lambda1) (sin lambda2))))))
(atan2
t_4
(-
t_0
(*
t_2
(cos
(*
(+ lambda1 lambda2)
(/ (- lambda1 lambda2) (+ lambda1 lambda2)))))))))))
double code(double lambda1, double lambda2, double phi1, double phi2) {
double t_0 = cos(phi1) * sin(phi2);
double t_1 = sin(lambda1) * cos(lambda2);
double t_2 = cos(phi2) * sin(phi1);
double t_3 = sin(-lambda2);
double t_4 = cos(phi2) * fma(t_3, cos(lambda1), t_1);
double tmp;
if (phi2 <= -750.0) {
tmp = atan2(fma(t_1, cos(phi2), (cos(phi2) * (cos(lambda1) * t_3))), (t_0 - (t_2 * cos((lambda1 - lambda2)))));
} else if (phi2 <= 9000000000.0) {
tmp = atan2(t_4, (t_0 - (sin(phi1) * fma(cos(lambda1), cos(lambda2), (sin(lambda1) * sin(lambda2))))));
} else {
tmp = atan2(t_4, (t_0 - (t_2 * cos(((lambda1 + lambda2) * ((lambda1 - lambda2) / (lambda1 + lambda2)))))));
}
return tmp;
}
function code(lambda1, lambda2, phi1, phi2) t_0 = Float64(cos(phi1) * sin(phi2)) t_1 = Float64(sin(lambda1) * cos(lambda2)) t_2 = Float64(cos(phi2) * sin(phi1)) t_3 = sin(Float64(-lambda2)) t_4 = Float64(cos(phi2) * fma(t_3, cos(lambda1), t_1)) tmp = 0.0 if (phi2 <= -750.0) tmp = atan(fma(t_1, cos(phi2), Float64(cos(phi2) * Float64(cos(lambda1) * t_3))), Float64(t_0 - Float64(t_2 * cos(Float64(lambda1 - lambda2))))); elseif (phi2 <= 9000000000.0) tmp = atan(t_4, Float64(t_0 - Float64(sin(phi1) * fma(cos(lambda1), cos(lambda2), Float64(sin(lambda1) * sin(lambda2)))))); else tmp = atan(t_4, Float64(t_0 - Float64(t_2 * cos(Float64(Float64(lambda1 + lambda2) * Float64(Float64(lambda1 - lambda2) / Float64(lambda1 + lambda2))))))); end return tmp end
code[lambda1_, lambda2_, phi1_, phi2_] := Block[{t$95$0 = N[(N[Cos[phi1], $MachinePrecision] * N[Sin[phi2], $MachinePrecision]), $MachinePrecision]}, Block[{t$95$1 = N[(N[Sin[lambda1], $MachinePrecision] * N[Cos[lambda2], $MachinePrecision]), $MachinePrecision]}, Block[{t$95$2 = N[(N[Cos[phi2], $MachinePrecision] * N[Sin[phi1], $MachinePrecision]), $MachinePrecision]}, Block[{t$95$3 = N[Sin[(-lambda2)], $MachinePrecision]}, Block[{t$95$4 = N[(N[Cos[phi2], $MachinePrecision] * N[(t$95$3 * N[Cos[lambda1], $MachinePrecision] + t$95$1), $MachinePrecision]), $MachinePrecision]}, If[LessEqual[phi2, -750.0], N[ArcTan[N[(t$95$1 * N[Cos[phi2], $MachinePrecision] + N[(N[Cos[phi2], $MachinePrecision] * N[(N[Cos[lambda1], $MachinePrecision] * t$95$3), $MachinePrecision]), $MachinePrecision]), $MachinePrecision] / N[(t$95$0 - N[(t$95$2 * N[Cos[N[(lambda1 - lambda2), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision], If[LessEqual[phi2, 9000000000.0], N[ArcTan[t$95$4 / N[(t$95$0 - N[(N[Sin[phi1], $MachinePrecision] * N[(N[Cos[lambda1], $MachinePrecision] * N[Cos[lambda2], $MachinePrecision] + N[(N[Sin[lambda1], $MachinePrecision] * N[Sin[lambda2], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision], N[ArcTan[t$95$4 / N[(t$95$0 - N[(t$95$2 * N[Cos[N[(N[(lambda1 + lambda2), $MachinePrecision] * N[(N[(lambda1 - lambda2), $MachinePrecision] / N[(lambda1 + lambda2), $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]), $MachinePrecision]), $MachinePrecision]], $MachinePrecision]]]]]]]]
\begin{array}{l}
\\
\begin{array}{l}
t_0 := \cos \phi_1 \cdot \sin \phi_2\\
t_1 := \sin \lambda_1 \cdot \cos \lambda_2\\
t_2 := \cos \phi_2 \cdot \sin \phi_1\\
t_3 := \sin \left(-\lambda_2\right)\\
t_4 := \cos \phi_2 \cdot \mathsf{fma}\left(t\_3, \cos \lambda_1, t\_1\right)\\
\mathbf{if}\;\phi_2 \leq -750:\\
\;\;\;\;\tan^{-1}_* \frac{\mathsf{fma}\left(t\_1, \cos \phi_2, \cos \phi_2 \cdot \left(\cos \lambda_1 \cdot t\_3\right)\right)}{t\_0 - t\_2 \cdot \cos \left(\lambda_1 - \lambda_2\right)}\\
\mathbf{elif}\;\phi_2 \leq 9000000000:\\
\;\;\;\;\tan^{-1}_* \frac{t\_4}{t\_0 - \sin \phi_1 \cdot \mathsf{fma}\left(\cos \lambda_1, \cos \lambda_2, \sin \lambda_1 \cdot \sin \lambda_2\right)}\\
\mathbf{else}:\\
\;\;\;\;\tan^{-1}_* \frac{t\_4}{t\_0 - t\_2 \cdot \cos \left(\left(\lambda_1 + \lambda_2\right) \cdot \frac{\lambda_1 - \lambda_2}{\lambda_1 + \lambda_2}\right)}\\
\end{array}
\end{array}
if phi2 < -750Initial program 78.3%
lift-*.f64N/A
*-commutativeN/A
lift-sin.f64N/A
lift--.f64N/A
sin-diffN/A
sub-negN/A
distribute-rgt-inN/A
lower-fma.f64N/A
lower-*.f64N/A
lower-sin.f64N/A
lower-cos.f64N/A
lower-*.f64N/A
distribute-rgt-neg-inN/A
sin-negN/A
lower-*.f64N/A
lower-cos.f64N/A
lower-sin.f64N/A
lower-neg.f6490.3
Applied rewrites90.3%
if -750 < phi2 < 9e9Initial program 82.3%
lift-sin.f64N/A
lift--.f64N/A
sub-negN/A
+-commutativeN/A
sin-sumN/A
cos-negN/A
*-commutativeN/A
lower-fma.f64N/A
lower-sin.f64N/A
lower-neg.f64N/A
lower-cos.f64N/A
lower-*.f64N/A
lower-sin.f64N/A
lower-cos.f6489.3
Applied rewrites89.3%
lift-cos.f64N/A
lift--.f64N/A
cos-diffN/A
+-commutativeN/A
lift-sin.f64N/A
*-commutativeN/A
lower-fma.f64N/A
lower-sin.f64N/A
lift-cos.f64N/A
lift-cos.f64N/A
lower-*.f6499.9
Applied rewrites99.9%
Taylor expanded in phi2 around 0
lower-*.f64N/A
lower-sin.f64N/A
lower-fma.f64N/A
lower-cos.f64N/A
lower-cos.f64N/A
lower-*.f64N/A
lower-sin.f64N/A
lower-sin.f6498.6
Applied rewrites98.6%
if 9e9 < phi2 Initial program 76.3%
lift-sin.f64N/A
lift--.f64N/A
sub-negN/A
+-commutativeN/A
sin-sumN/A
cos-negN/A
*-commutativeN/A
lower-fma.f64N/A
lower-sin.f64N/A
lower-neg.f64N/A
lower-cos.f64N/A
lower-*.f64N/A
lower-sin.f64N/A
lower-cos.f6489.7
Applied rewrites89.7%
lift--.f64N/A
flip--N/A
flip3-+N/A
associate-/r/N/A
lower-*.f64N/A
lower-/.f64N/A
difference-of-squaresN/A
lift-+.f64N/A
lift--.f64N/A
lower-*.f64N/A
lift-+.f64N/A
+-commutativeN/A
lower-+.f64N/A
cube-multN/A
lower-fma.f64N/A
lower-*.f64N/A
cube-multN/A
lower-*.f64N/A
lower-*.f64N/A
+-commutativeN/A
distribute-rgt-out--N/A
Applied rewrites38.8%
lift-*.f64N/A
lift-/.f64N/A
lift-fma.f64N/A
+-commutativeN/A
lift-*.f64N/A
lift--.f64N/A
distribute-rgt-out--N/A
associate-/r/N/A
lift-*.f64N/A
lift-fma.f64N/A
lift-*.f64N/A
cube-unmultN/A
lift-*.f64N/A
lift-*.f64N/A
cube-unmultN/A
flip3-+N/A
+-commutativeN/A
lift-+.f64N/A
Applied rewrites89.7%
Final simplification94.3%
herbie shell --seed 2024227
(FPCore (lambda1 lambda2 phi1 phi2)
:name "Bearing on a great circle"
:precision binary64
(atan2 (* (sin (- lambda1 lambda2)) (cos phi2)) (- (* (cos phi1) (sin phi2)) (* (* (sin phi1) (cos phi2)) (cos (- lambda1 lambda2))))))