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 - \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{expm1}\left(\mathsf{log1p}\left(\sin \left(\frac{\lambda_1 - \lambda_2}{2}\right)\right)\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 \sqrt[3]{{\left(\sin \left(\frac{\lambda_1 - \lambda_2}{2}\right)\right)}^{3}}\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 r70184 = R;
double r70185 = 2.0;
double r70186 = phi1;
double r70187 = phi2;
double r70188 = r70186 - r70187;
double r70189 = r70188 / r70185;
double r70190 = sin(r70189);
double r70191 = pow(r70190, r70185);
double r70192 = cos(r70186);
double r70193 = cos(r70187);
double r70194 = r70192 * r70193;
double r70195 = lambda1;
double r70196 = lambda2;
double r70197 = r70195 - r70196;
double r70198 = r70197 / r70185;
double r70199 = sin(r70198);
double r70200 = r70194 * r70199;
double r70201 = r70200 * r70199;
double r70202 = r70191 + r70201;
double r70203 = sqrt(r70202);
double r70204 = 1.0;
double r70205 = r70204 - r70202;
double r70206 = sqrt(r70205);
double r70207 = atan2(r70203, r70206);
double r70208 = r70185 * r70207;
double r70209 = r70184 * r70208;
return r70209;
}
double f(double R, double lambda1, double lambda2, double phi1, double phi2) {
double r70210 = R;
double r70211 = 2.0;
double r70212 = phi1;
double r70213 = phi2;
double r70214 = r70212 - r70213;
double r70215 = r70214 / r70211;
double r70216 = sin(r70215);
double r70217 = pow(r70216, r70211);
double r70218 = cos(r70212);
double r70219 = cos(r70213);
double r70220 = r70218 * r70219;
double r70221 = lambda1;
double r70222 = lambda2;
double r70223 = r70221 - r70222;
double r70224 = r70223 / r70211;
double r70225 = sin(r70224);
double r70226 = r70220 * r70225;
double r70227 = log1p(r70225);
double r70228 = expm1(r70227);
double r70229 = r70226 * r70228;
double r70230 = r70217 + r70229;
double r70231 = sqrt(r70230);
double r70232 = 1.0;
double r70233 = 3.0;
double r70234 = pow(r70225, r70233);
double r70235 = cbrt(r70234);
double r70236 = r70220 * r70235;
double r70237 = expm1(r70225);
double r70238 = log1p(r70237);
double r70239 = r70236 * r70238;
double r70240 = r70217 + r70239;
double r70241 = r70232 - r70240;
double r70242 = sqrt(r70241);
double r70243 = atan2(r70231, r70242);
double r70244 = r70211 * r70243;
double r70245 = r70210 * r70244;
return r70245;
}



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 expm1-log1p-u24.2
rmApplied log1p-expm1-u24.2
rmApplied add-cbrt-cube24.2
Simplified24.2
Final simplification24.2
herbie shell --seed 2019235 +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))))))))))