\frac{c0}{2 \cdot w} \cdot \left(\frac{c0 \cdot \left(d \cdot d\right)}{\left(w \cdot h\right) \cdot \left(D \cdot D\right)} + \sqrt{\frac{c0 \cdot \left(d \cdot d\right)}{\left(w \cdot h\right) \cdot \left(D \cdot D\right)} \cdot \frac{c0 \cdot \left(d \cdot d\right)}{\left(w \cdot h\right) \cdot \left(D \cdot D\right)} - M \cdot M}\right)\begin{array}{l}
\mathbf{if}\;d \le 1.028014946374977112793459567345674978253 \cdot 10^{-41} \lor \neg \left(d \le 1.879132275263087348510262122271649682884 \cdot 10^{52}\right):\\
\;\;\;\;\mathsf{log1p}\left(\mathsf{expm1}\left(0\right)\right)\\
\mathbf{else}:\\
\;\;\;\;\frac{c0}{2 \cdot w} \cdot \left(2 \cdot \frac{{d}^{2} \cdot c0}{w \cdot \left({D}^{2} \cdot h\right)}\right)\\
\end{array}double f(double c0, double w, double h, double D, double d, double M) {
double r139439 = c0;
double r139440 = 2.0;
double r139441 = w;
double r139442 = r139440 * r139441;
double r139443 = r139439 / r139442;
double r139444 = d;
double r139445 = r139444 * r139444;
double r139446 = r139439 * r139445;
double r139447 = h;
double r139448 = r139441 * r139447;
double r139449 = D;
double r139450 = r139449 * r139449;
double r139451 = r139448 * r139450;
double r139452 = r139446 / r139451;
double r139453 = r139452 * r139452;
double r139454 = M;
double r139455 = r139454 * r139454;
double r139456 = r139453 - r139455;
double r139457 = sqrt(r139456);
double r139458 = r139452 + r139457;
double r139459 = r139443 * r139458;
return r139459;
}
double f(double c0, double w, double h, double D, double d, double __attribute__((unused)) M) {
double r139460 = d;
double r139461 = 1.0280149463749771e-41;
bool r139462 = r139460 <= r139461;
double r139463 = 1.8791322752630873e+52;
bool r139464 = r139460 <= r139463;
double r139465 = !r139464;
bool r139466 = r139462 || r139465;
double r139467 = 0.0;
double r139468 = expm1(r139467);
double r139469 = log1p(r139468);
double r139470 = c0;
double r139471 = 2.0;
double r139472 = w;
double r139473 = r139471 * r139472;
double r139474 = r139470 / r139473;
double r139475 = 2.0;
double r139476 = pow(r139460, r139475);
double r139477 = r139476 * r139470;
double r139478 = D;
double r139479 = pow(r139478, r139475);
double r139480 = h;
double r139481 = r139479 * r139480;
double r139482 = r139472 * r139481;
double r139483 = r139477 / r139482;
double r139484 = r139475 * r139483;
double r139485 = r139474 * r139484;
double r139486 = r139466 ? r139469 : r139485;
return r139486;
}



Bits error versus c0



Bits error versus w



Bits error versus h



Bits error versus D



Bits error versus d



Bits error versus M
Results
if d < 1.0280149463749771e-41 or 1.8791322752630873e+52 < d Initial program 59.9
Taylor expanded around inf 34.6
rmApplied log1p-expm1-u34.6
Simplified32.7
if 1.0280149463749771e-41 < d < 1.8791322752630873e+52Initial program 52.7
rmApplied associate-/l*53.8
Taylor expanded around inf 54.7
Final simplification34.5
herbie shell --seed 2019323 +o rules:numerics
(FPCore (c0 w h D d M)
:name "Henrywood and Agarwal, Equation (13)"
:precision binary64
(* (/ c0 (* 2 w)) (+ (/ (* c0 (* d d)) (* (* w h) (* D D))) (sqrt (- (* (/ (* c0 (* d d)) (* (* w h) (* D D))) (/ (* c0 (* d d)) (* (* w h) (* D D)))) (* M M))))))