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 \mathsf{log1p}\left(\mathsf{expm1}\left(\sin \left(\frac{\lambda_1 - \lambda_2}{2}\right)\right)\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 \log \left(e^{\sin \left(\frac{\lambda_1 - \lambda_2}{2}\right)}\right)\right)}}\right)double f(double R, double lambda1, double lambda2, double phi1, double phi2) {
double r79823 = R;
double r79824 = 2.0;
double r79825 = phi1;
double r79826 = phi2;
double r79827 = r79825 - r79826;
double r79828 = r79827 / r79824;
double r79829 = sin(r79828);
double r79830 = pow(r79829, r79824);
double r79831 = cos(r79825);
double r79832 = cos(r79826);
double r79833 = r79831 * r79832;
double r79834 = lambda1;
double r79835 = lambda2;
double r79836 = r79834 - r79835;
double r79837 = r79836 / r79824;
double r79838 = sin(r79837);
double r79839 = r79833 * r79838;
double r79840 = r79839 * r79838;
double r79841 = r79830 + r79840;
double r79842 = sqrt(r79841);
double r79843 = 1.0;
double r79844 = r79843 - r79841;
double r79845 = sqrt(r79844);
double r79846 = atan2(r79842, r79845);
double r79847 = r79824 * r79846;
double r79848 = r79823 * r79847;
return r79848;
}
double f(double R, double lambda1, double lambda2, double phi1, double phi2) {
double r79849 = R;
double r79850 = 2.0;
double r79851 = phi1;
double r79852 = r79851 / r79850;
double r79853 = sin(r79852);
double r79854 = phi2;
double r79855 = r79854 / r79850;
double r79856 = cos(r79855);
double r79857 = r79853 * r79856;
double r79858 = cos(r79852);
double r79859 = sin(r79855);
double r79860 = r79858 * r79859;
double r79861 = r79857 - r79860;
double r79862 = pow(r79861, r79850);
double r79863 = cos(r79851);
double r79864 = cos(r79854);
double r79865 = r79863 * r79864;
double r79866 = lambda1;
double r79867 = lambda2;
double r79868 = r79866 - r79867;
double r79869 = r79868 / r79850;
double r79870 = sin(r79869);
double r79871 = r79865 * r79870;
double r79872 = expm1(r79870);
double r79873 = log1p(r79872);
double r79874 = r79871 * r79873;
double r79875 = r79862 + r79874;
double r79876 = sqrt(r79875);
double r79877 = 1.0;
double r79878 = exp(r79870);
double r79879 = log(r79878);
double r79880 = r79871 * r79879;
double r79881 = r79862 + r79880;
double r79882 = r79877 - r79881;
double r79883 = sqrt(r79882);
double r79884 = atan2(r79876, r79883);
double r79885 = r79850 * r79884;
double r79886 = r79849 * r79885;
return r79886;
}



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-diff13.8
rmApplied add-log-exp13.8
rmApplied log1p-expm1-u13.8
Final simplification13.8
herbie shell --seed 2020065 +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))))))))))