lib/libm/src/e_sinhf.c

   1 /* e_sinhf.c -- float version of e_sinh.c.
   2  * Conversion to float by Ian Lance Taylor, Cygnus Support, ian@cygnus.com.
   3  */
   4
   5 /*
   6  * ====================================================
   7  * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
   8  *
   9  * Developed at SunPro, a Sun Microsystems, Inc. business.
  10  * Permission to use, copy, modify, and distribute this
  11  * software is freely granted, provided that this notice
  12  * is preserved.
  13  * ====================================================
  14  *
  15  * $FreeBSD: head/lib/msun/src/e_sinhf.c 226598 2011-10-21 06:28:47Z das $
  16  */
  17
  18 #include "math.h"
  19 #include "math_private.h"
  20
  21 static const float one = 1.0, shuge = 1.0e37;
  22
  23 float
  24 __ieee754_sinhf(float x)
  25 {
  26         float t,h;
  27         int32_t ix,jx;
  28
  29         GET_FLOAT_WORD(jx,x);
  30         ix = jx&0x7fffffff;
  31
  32     /* x is INF or NaN */
  33         if(ix>=0x7f800000) return x+x;
  34
  35         h = 0.5;
  36         if (jx<0) h = -h;
  37     /* |x| in [0,9], return sign(x)*0.5*(E+E/(E+1))) */
  38         if (ix < 0x41100000) {          /* |x|<9 */
  39             if (ix<0x39800000)          /* |x|<2**-12 */
  40                 if(shuge+x>one) return x;/* sinh(tiny) = tiny with inexact */
  41             t = expm1f(fabsf(x));
  42             if(ix<0x3f800000) return h*((float)2.0*t-t*t/(t+one));
  43             return h*(t+t/(t+one));
  44         }
  45
  46     /* |x| in [9, logf(maxfloat)] return 0.5*exp(|x|) */
  47         if (ix < 0x42b17217)  return h*__ieee754_expf(fabsf(x));
  48
  49     /* |x| in [logf(maxfloat), overflowthresold] */
  50         if (ix<=0x42b2d4fc)
  51             return h*2.0F*__ldexp_expf(fabsf(x), -1);
  52
  53     /* |x| > overflowthresold, sinh(x) overflow */
  54         return x*shuge;
  55 }