sysdeps/ieee754/dbl-64/s_scalbln.c

   1 /*
   2  * ====================================================
   3  * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
   4  *
   5  * Developed at SunPro, a Sun Microsystems, Inc. business.
   6  * Permission to use, copy, modify, and distribute this
   7  * software is freely granted, provided that this notice
   8  * is preserved.
   9  * ====================================================
  10  */
  11
  12 /*
  13  * scalbn (double x, int n)
  14  * scalbn(x,n) returns x* 2**n  computed by  exponent
  15  * manipulation rather than by actually performing an
  16  * exponentiation or a multiplication.
  17  */
  18
  19 #include <math.h>
  20 #include <math_private.h>
  21
  22 static const double
  23   two54 = 1.80143985094819840000e+16,  /* 0x43500000, 0x00000000 */
  24   twom54 = 5.55111512312578270212e-17, /* 0x3C900000, 0x00000000 */
  25   huge = 1.0e+300,
  26   tiny = 1.0e-300;
  27
  28 double
  29 __scalbln (double x, long int n)
  30 {
  31   int32_t k, hx, lx;
  32   EXTRACT_WORDS (hx, lx, x);
  33   k = (hx & 0x7ff00000) >> 20;                  /* extract exponent */
  34   if (__builtin_expect (k == 0, 0))             /* 0 or subnormal x */
  35     {
  36       if ((lx | (hx & 0x7fffffff)) == 0)
  37         return x;                                  /* +-0 */
  38       x *= two54;
  39       GET_HIGH_WORD (hx, x);
  40       k = ((hx & 0x7ff00000) >> 20) - 54;
  41     }
  42   if (__builtin_expect (k == 0x7ff, 0))
  43     return x + x;                                       /* NaN or Inf */
  44   if (__builtin_expect (n < -50000, 0))
  45     return tiny * __copysign (tiny, x);   /*underflow*/
  46   if (__builtin_expect (n > 50000 || k + n > 0x7fe, 0))
  47     return huge * __copysign (huge, x);   /* overflow  */
  48   /* Now k and n are bounded we know that k = k+n does not
  49      overflow.  */
  50   k = k + n;
  51   if (__builtin_expect (k > 0, 1))              /* normal result */
  52     {
  53       SET_HIGH_WORD (x, (hx & 0x800fffff) | (k << 20)); return x;
  54     }
  55   if (k <= -54)
  56     return tiny * __copysign (tiny, x);         /*underflow*/
  57   k += 54;                                      /* subnormal result */
  58   SET_HIGH_WORD (x, (hx & 0x800fffff) | (k << 20));
  59   return x * twom54;
  60 }
  61 weak_alias (__scalbln, scalbln)
  62 #ifdef NO_LONG_DOUBLE
  63 strong_alias (__scalbln, __scalblnl)
  64 weak_alias (__scalbln, scalblnl)
  65 #endif