R \cdot \left(2 \cdot \tan^{-1}_* \frac{\sqrt{{\left(\sin \left(\frac{\phi_1 - \phi_2}{2}\right)\right)}^{2} + \left(\left(\cos \phi_1 \cdot \cos \phi_2\right) \cdot \sin \left(\frac{\lambda_1 - \lambda_2}{2}\right)\right) \cdot \sin \left(\frac{\lambda_1 - \lambda_2}{2}\right)}}{\sqrt{1 - \left({\left(\sin \left(\frac{\phi_1 - \phi_2}{2}\right)\right)}^{2} + \left(\left(\cos \phi_1 \cdot \cos \phi_2\right) \cdot \sin \left(\frac{\lambda_1 - \lambda_2}{2}\right)\right) \cdot \sin \left(\frac{\lambda_1 - \lambda_2}{2}\right)\right)}}\right)R \cdot \left(2 \cdot \tan^{-1}_* \frac{\sqrt{{\left(\sin \left(\frac{\phi_1}{2}\right) \cdot \cos \left(\frac{\phi_2}{2}\right) - \cos \left(\frac{\phi_1}{2}\right) \cdot \sin \left(\frac{\phi_2}{2}\right)\right)}^{2} + \left(\left(\cos \phi_1 \cdot \cos \phi_2\right) \cdot \mathsf{log1p}\left(\mathsf{expm1}\left(\sin \left(\frac{\lambda_1 - \lambda_2}{2}\right)\right)\right)\right) \cdot \sin \left(\frac{\lambda_1 - \lambda_2}{2}\right)}}{\sqrt{1 - \left({\left(\sin \left(\frac{\phi_1}{2}\right) \cdot \cos \left(\frac{\phi_2}{2}\right) - \cos \left(\frac{\phi_1}{2}\right) \cdot \sin \left(\frac{\phi_2}{2}\right)\right)}^{2} + \left(\left(\cos \phi_1 \cdot \cos \phi_2\right) \cdot \sqrt[3]{{\left(\sin \left(\frac{\lambda_1 - \lambda_2}{2}\right)\right)}^{3}}\right) \cdot \sin \left(\frac{\lambda_1 - \lambda_2}{2}\right)\right)}}\right)double f(double R, double lambda1, double lambda2, double phi1, double phi2) {
double r74802 = R;
double r74803 = 2.0;
double r74804 = phi1;
double r74805 = phi2;
double r74806 = r74804 - r74805;
double r74807 = r74806 / r74803;
double r74808 = sin(r74807);
double r74809 = pow(r74808, r74803);
double r74810 = cos(r74804);
double r74811 = cos(r74805);
double r74812 = r74810 * r74811;
double r74813 = lambda1;
double r74814 = lambda2;
double r74815 = r74813 - r74814;
double r74816 = r74815 / r74803;
double r74817 = sin(r74816);
double r74818 = r74812 * r74817;
double r74819 = r74818 * r74817;
double r74820 = r74809 + r74819;
double r74821 = sqrt(r74820);
double r74822 = 1.0;
double r74823 = r74822 - r74820;
double r74824 = sqrt(r74823);
double r74825 = atan2(r74821, r74824);
double r74826 = r74803 * r74825;
double r74827 = r74802 * r74826;
return r74827;
}
double f(double R, double lambda1, double lambda2, double phi1, double phi2) {
double r74828 = R;
double r74829 = 2.0;
double r74830 = phi1;
double r74831 = r74830 / r74829;
double r74832 = sin(r74831);
double r74833 = phi2;
double r74834 = r74833 / r74829;
double r74835 = cos(r74834);
double r74836 = r74832 * r74835;
double r74837 = cos(r74831);
double r74838 = sin(r74834);
double r74839 = r74837 * r74838;
double r74840 = r74836 - r74839;
double r74841 = pow(r74840, r74829);
double r74842 = cos(r74830);
double r74843 = cos(r74833);
double r74844 = r74842 * r74843;
double r74845 = lambda1;
double r74846 = lambda2;
double r74847 = r74845 - r74846;
double r74848 = r74847 / r74829;
double r74849 = sin(r74848);
double r74850 = expm1(r74849);
double r74851 = log1p(r74850);
double r74852 = r74844 * r74851;
double r74853 = r74852 * r74849;
double r74854 = r74841 + r74853;
double r74855 = sqrt(r74854);
double r74856 = 1.0;
double r74857 = 3.0;
double r74858 = pow(r74849, r74857);
double r74859 = cbrt(r74858);
double r74860 = r74844 * r74859;
double r74861 = r74860 * r74849;
double r74862 = r74841 + r74861;
double r74863 = r74856 - r74862;
double r74864 = sqrt(r74863);
double r74865 = atan2(r74855, r74864);
double r74866 = r74829 * r74865;
double r74867 = r74828 * r74866;
return r74867;
}



Bits error versus R



Bits error versus lambda1



Bits error versus lambda2



Bits error versus phi1



Bits error versus phi2
Results
Initial program 24.4
rmApplied div-sub24.4
Applied sin-diff23.8
rmApplied div-sub23.8
Applied sin-diff13.7
rmApplied log1p-expm1-u13.7
rmApplied add-cbrt-cube13.7
Simplified13.7
Final simplification13.7
herbie shell --seed 2020046 +o rules:numerics
(FPCore (R lambda1 lambda2 phi1 phi2)
:name "Distance on a great circle"
:precision binary64
(* R (* 2 (atan2 (sqrt (+ (pow (sin (/ (- phi1 phi2) 2)) 2) (* (* (* (cos phi1) (cos phi2)) (sin (/ (- lambda1 lambda2) 2))) (sin (/ (- lambda1 lambda2) 2))))) (sqrt (- 1 (+ (pow (sin (/ (- phi1 phi2) 2)) 2) (* (* (* (cos phi1) (cos phi2)) (sin (/ (- lambda1 lambda2) 2))) (sin (/ (- lambda1 lambda2) 2))))))))))