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 \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}{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 \sin \left(\frac{\lambda_1 - \lambda_2}{2}\right)\right) \cdot \mathsf{log1p}\left(\mathsf{expm1}\left(\sqrt[3]{{\left(\sin \left(\frac{\lambda_1 - \lambda_2}{2}\right)\right)}^{3}}\right)\right)\right)}}\right)double f(double R, double lambda1, double lambda2, double phi1, double phi2) {
double r119770 = R;
double r119771 = 2.0;
double r119772 = phi1;
double r119773 = phi2;
double r119774 = r119772 - r119773;
double r119775 = r119774 / r119771;
double r119776 = sin(r119775);
double r119777 = pow(r119776, r119771);
double r119778 = cos(r119772);
double r119779 = cos(r119773);
double r119780 = r119778 * r119779;
double r119781 = lambda1;
double r119782 = lambda2;
double r119783 = r119781 - r119782;
double r119784 = r119783 / r119771;
double r119785 = sin(r119784);
double r119786 = r119780 * r119785;
double r119787 = r119786 * r119785;
double r119788 = r119777 + r119787;
double r119789 = sqrt(r119788);
double r119790 = 1.0;
double r119791 = r119790 - r119788;
double r119792 = sqrt(r119791);
double r119793 = atan2(r119789, r119792);
double r119794 = r119771 * r119793;
double r119795 = r119770 * r119794;
return r119795;
}
double f(double R, double lambda1, double lambda2, double phi1, double phi2) {
double r119796 = R;
double r119797 = 2.0;
double r119798 = phi1;
double r119799 = r119798 / r119797;
double r119800 = sin(r119799);
double r119801 = phi2;
double r119802 = r119801 / r119797;
double r119803 = cos(r119802);
double r119804 = r119800 * r119803;
double r119805 = cos(r119799);
double r119806 = sin(r119802);
double r119807 = r119805 * r119806;
double r119808 = r119804 - r119807;
double r119809 = pow(r119808, r119797);
double r119810 = cos(r119798);
double r119811 = cos(r119801);
double r119812 = r119810 * r119811;
double r119813 = lambda1;
double r119814 = lambda2;
double r119815 = r119813 - r119814;
double r119816 = r119815 / r119797;
double r119817 = sin(r119816);
double r119818 = r119812 * r119817;
double r119819 = r119818 * r119817;
double r119820 = r119809 + r119819;
double r119821 = sqrt(r119820);
double r119822 = 1.0;
double r119823 = 3.0;
double r119824 = pow(r119817, r119823);
double r119825 = cbrt(r119824);
double r119826 = expm1(r119825);
double r119827 = log1p(r119826);
double r119828 = r119818 * r119827;
double r119829 = r119809 + r119828;
double r119830 = r119822 - r119829;
double r119831 = sqrt(r119830);
double r119832 = atan2(r119821, r119831);
double r119833 = r119797 * r119832;
double r119834 = r119796 * r119833;
return r119834;
}



Bits error versus R



Bits error versus lambda1



Bits error versus lambda2



Bits error versus phi1



Bits error versus phi2
Results
Initial program 24.3
rmApplied div-sub24.3
Applied sin-diff23.7
rmApplied div-sub23.7
Applied sin-diff14.2
rmApplied log1p-expm1-u14.2
rmApplied add-cbrt-cube14.2
Simplified14.2
Final simplification14.2
herbie shell --seed 2020027 +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))))))))))