double hfsq = 0.5 * f * f;
if (hu == 0) { /* |f| < 2**-20 */
if (f == zero) {
return k == 0 ? zero : k * ln2_hi + (c + k * ln2_lo) ;
}
final double R = hfsq * (1.0 - 0.66666666666666666 * f);
return k == 0 ? f - R : k * ln2_hi - ((R - (k * ln2_lo + c)) - f);
}
final double s = f / (2.0 + f);
final double z = s * s;
final double R =
s*(hfsq +z*(Lp1+z*(Lp2+z*(Lp3+z*(Lp4+z*(Lp5+z*(Lp6+z*Lp7)))))));
return k == 0 ? f - (hfsq - R) :
k * ln2_hi - ((hfsq - (R + (k * ln2_lo + c))) - f);