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 \mathsf{expm1}\left(\mathsf{log1p}\left(\sin \left(\frac{\lambda_1 - \lambda_2}{2}\right)\right)\right)\right) \cdot \mathsf{log1p}\left(\mathsf{expm1}\left(\sin \left(\frac{\lambda_1 - \lambda_2}{2}\right)\right)\right)\right)}}\right)double f(double R, double lambda1, double lambda2, double phi1, double phi2) {
double r9007177 = R;
double r9007178 = 2.0;
double r9007179 = phi1;
double r9007180 = phi2;
double r9007181 = r9007179 - r9007180;
double r9007182 = r9007181 / r9007178;
double r9007183 = sin(r9007182);
double r9007184 = pow(r9007183, r9007178);
double r9007185 = cos(r9007179);
double r9007186 = cos(r9007180);
double r9007187 = r9007185 * r9007186;
double r9007188 = lambda1;
double r9007189 = lambda2;
double r9007190 = r9007188 - r9007189;
double r9007191 = r9007190 / r9007178;
double r9007192 = sin(r9007191);
double r9007193 = r9007187 * r9007192;
double r9007194 = r9007193 * r9007192;
double r9007195 = r9007184 + r9007194;
double r9007196 = sqrt(r9007195);
double r9007197 = 1.0;
double r9007198 = r9007197 - r9007195;
double r9007199 = sqrt(r9007198);
double r9007200 = atan2(r9007196, r9007199);
double r9007201 = r9007178 * r9007200;
double r9007202 = r9007177 * r9007201;
return r9007202;
}
double f(double R, double lambda1, double lambda2, double phi1, double phi2) {
double r9007203 = R;
double r9007204 = 2.0;
double r9007205 = phi1;
double r9007206 = r9007205 / r9007204;
double r9007207 = sin(r9007206);
double r9007208 = phi2;
double r9007209 = r9007208 / r9007204;
double r9007210 = cos(r9007209);
double r9007211 = r9007207 * r9007210;
double r9007212 = cos(r9007206);
double r9007213 = sin(r9007209);
double r9007214 = r9007212 * r9007213;
double r9007215 = r9007211 - r9007214;
double r9007216 = pow(r9007215, r9007204);
double r9007217 = cos(r9007205);
double r9007218 = cos(r9007208);
double r9007219 = r9007217 * r9007218;
double r9007220 = lambda1;
double r9007221 = lambda2;
double r9007222 = r9007220 - r9007221;
double r9007223 = r9007222 / r9007204;
double r9007224 = sin(r9007223);
double r9007225 = r9007219 * r9007224;
double r9007226 = r9007225 * r9007224;
double r9007227 = r9007216 + r9007226;
double r9007228 = sqrt(r9007227);
double r9007229 = 1.0;
double r9007230 = log1p(r9007224);
double r9007231 = expm1(r9007230);
double r9007232 = r9007219 * r9007231;
double r9007233 = expm1(r9007224);
double r9007234 = log1p(r9007233);
double r9007235 = r9007232 * r9007234;
double r9007236 = r9007216 + r9007235;
double r9007237 = r9007229 - r9007236;
double r9007238 = sqrt(r9007237);
double r9007239 = atan2(r9007228, r9007238);
double r9007240 = r9007204 * r9007239;
double r9007241 = r9007203 * r9007240;
return r9007241;
}



Bits error versus R



Bits error versus lambda1



Bits error versus lambda2



Bits error versus phi1



Bits error versus phi2
Results
Initial program 24.2
rmApplied div-sub24.2
Applied sin-diff23.6
rmApplied div-sub23.6
Applied sin-diff13.7
rmApplied log1p-expm1-u13.7
rmApplied expm1-log1p-u13.7
Final simplification13.7
herbie shell --seed 2019173 +o rules:numerics
(FPCore (R lambda1 lambda2 phi1 phi2)
:name "Distance on a great circle"
(* R (* 2.0 (atan2 (sqrt (+ (pow (sin (/ (- phi1 phi2) 2.0)) 2.0) (* (* (* (cos phi1) (cos phi2)) (sin (/ (- lambda1 lambda2) 2.0))) (sin (/ (- lambda1 lambda2) 2.0))))) (sqrt (- 1.0 (+ (pow (sin (/ (- phi1 phi2) 2.0)) 2.0) (* (* (* (cos phi1) (cos phi2)) (sin (/ (- lambda1 lambda2) 2.0))) (sin (/ (- lambda1 lambda2) 2.0))))))))))