\sqrt[3]{\frac{1}{2 \cdot a} \cdot \left(\left(-g\right) + \sqrt{g \cdot g - h \cdot h}\right)} + \sqrt[3]{\frac{1}{2 \cdot a} \cdot \left(\left(-g\right) - \sqrt{g \cdot g - h \cdot h}\right)}\begin{array}{l}
\mathbf{if}\;g \le 1.7601370127167322 \cdot 10^{-167}:\\
\;\;\;\;\frac{\sqrt[3]{1 \cdot \left(\left(-g\right) + -1 \cdot g\right)}}{\sqrt[3]{2 \cdot a}} + \sqrt[3]{\frac{1}{2 \cdot a} \cdot \left(\left(-g\right) - \sqrt{g \cdot g - h \cdot h}\right)}\\
\mathbf{else}:\\
\;\;\;\;\frac{\sqrt[3]{1 \cdot \left(\left(-g\right) + \sqrt{g \cdot g - h \cdot h}\right)}}{\sqrt[3]{2 \cdot a}} + \sqrt[3]{\frac{1}{2 \cdot a}} \cdot \sqrt[3]{\left(-g\right) - \left(\sqrt[3]{\sqrt{g \cdot g - h \cdot h}} \cdot \sqrt[3]{\sqrt{g \cdot g - h \cdot h}}\right) \cdot \sqrt[3]{\sqrt{g \cdot g - h \cdot h}}}\\
\end{array}double f(double g, double h, double a) {
double r435 = 1.0;
double r436 = 2.0;
double r437 = a;
double r438 = r436 * r437;
double r439 = r435 / r438;
double r440 = g;
double r441 = -r440;
double r442 = r440 * r440;
double r443 = h;
double r444 = r443 * r443;
double r445 = r442 - r444;
double r446 = sqrt(r445);
double r447 = r441 + r446;
double r448 = r439 * r447;
double r449 = cbrt(r448);
double r450 = r441 - r446;
double r451 = r439 * r450;
double r452 = cbrt(r451);
double r453 = r449 + r452;
return r453;
}
double f(double g, double h, double a) {
double r454 = g;
double r455 = 1.760137012716732e-167;
bool r456 = r454 <= r455;
double r457 = 1.0;
double r458 = -r454;
double r459 = -1.0;
double r460 = r459 * r454;
double r461 = r458 + r460;
double r462 = r457 * r461;
double r463 = cbrt(r462);
double r464 = 2.0;
double r465 = a;
double r466 = r464 * r465;
double r467 = cbrt(r466);
double r468 = r463 / r467;
double r469 = r457 / r466;
double r470 = r454 * r454;
double r471 = h;
double r472 = r471 * r471;
double r473 = r470 - r472;
double r474 = sqrt(r473);
double r475 = r458 - r474;
double r476 = r469 * r475;
double r477 = cbrt(r476);
double r478 = r468 + r477;
double r479 = r458 + r474;
double r480 = r457 * r479;
double r481 = cbrt(r480);
double r482 = r481 / r467;
double r483 = cbrt(r469);
double r484 = cbrt(r474);
double r485 = r484 * r484;
double r486 = r485 * r484;
double r487 = r458 - r486;
double r488 = cbrt(r487);
double r489 = r483 * r488;
double r490 = r482 + r489;
double r491 = r456 ? r478 : r490;
return r491;
}



Bits error versus g



Bits error versus h



Bits error versus a
Results
if g < 1.760137012716732e-167Initial program 37.6
rmApplied associate-*l/37.6
Applied cbrt-div34.2
Taylor expanded around -inf 32.9
if 1.760137012716732e-167 < g Initial program 35.5
rmApplied associate-*l/35.5
Applied cbrt-div35.4
rmApplied cbrt-prod31.4
rmApplied add-cube-cbrt31.4
Final simplification32.2
herbie shell --seed 2020025 +o rules:numerics
(FPCore (g h a)
:name "2-ancestry mixing, positive discriminant"
:precision binary64
(+ (cbrt (* (/ 1 (* 2 a)) (+ (- g) (sqrt (- (* g g) (* h h)))))) (cbrt (* (/ 1 (* 2 a)) (- (- g) (sqrt (- (* g g) (* h h))))))))