github.com/afumu/libc@v0.0.6/musl/src/math/sinh.c (about)

     1  #include "libm.h"
     2  
     3  /* sinh(x) = (exp(x) - 1/exp(x))/2
     4   *         = (exp(x)-1 + (exp(x)-1)/exp(x))/2
     5   *         = x + x^3/6 + o(x^5)
     6   */
     7  double sinh(double x)
     8  {
     9  	union {double f; uint64_t i;} u = {.f = x};
    10  	uint32_t w;
    11  	double t, h, absx;
    12  
    13  	h = 0.5;
    14  	if (u.i >> 63)
    15  		h = -h;
    16  	/* |x| */
    17  	u.i &= (uint64_t)-1/2;
    18  	absx = u.f;
    19  	w = u.i >> 32;
    20  
    21  	/* |x| < log(DBL_MAX) */
    22  	if (w < 0x40862e42) {
    23  		t = expm1(absx);
    24  		if (w < 0x3ff00000) {
    25  			if (w < 0x3ff00000 - (26<<20))
    26  				/* note: inexact and underflow are raised by expm1 */
    27  				/* note: this branch avoids spurious underflow */
    28  				return x;
    29  			return h*(2*t - t*t/(t+1));
    30  		}
    31  		/* note: |x|>log(0x1p26)+eps could be just h*exp(x) */
    32  		return h*(t + t/(t+1));
    33  	}
    34  
    35  	/* |x| > log(DBL_MAX) or nan */
    36  	/* note: the result is stored to handle overflow */
    37  	t = __expo2(absx, 2*h);
    38  	return t;
    39  }