1  /* s_scalbnl.c -- long double version of s_scalbn.c.
       2   */
       3  
       4  /*
       5   * ====================================================
       6   * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
       7   *
       8   * Developed at SunPro, a Sun Microsystems, Inc. business.
       9   * Permission to use, copy, modify, and distribute this
      10   * software is freely granted, provided that this notice
      11   * is preserved.
      12   * ====================================================
      13   */
      14  
      15  /*
      16   * scalbnl (long double x, int n)
      17   * scalbnl(x,n) returns x* 2**n  computed by  exponent
      18   * manipulation rather than by actually performing an
      19   * exponentiation or a multiplication.
      20   */
      21  
      22  #include <math.h>
      23  #include <math_private.h>
      24  
      25  static const long double
      26  two63   =  0x1p63L,
      27  twom64  =  0x1p-64L,
      28  huge   = 1.0e+4900L,
      29  tiny   = 1.0e-4900L;
      30  
      31  long double
      32  __scalblnl (long double x, long int n)
      33  {
      34  	int32_t k,es,hx,lx;
      35  	GET_LDOUBLE_WORDS(es,hx,lx,x);
      36  	k = es&0x7fff;				/* extract exponent */
      37  	if (__builtin_expect(k==0, 0)) {	/* 0 or subnormal x */
      38  	    if ((lx|(hx&0x7fffffff))==0) return x; /* +-0 */
      39  	    x *= two63;
      40  	    GET_LDOUBLE_EXP(es,x);
      41  	    k = (es&0x7fff) - 63;
      42  	    }
      43  	if (__builtin_expect(k==0x7fff, 0)) return x+x;	/* NaN or Inf */
      44  	if (__builtin_expect(n< -50000, 0))
      45  	  return tiny*copysignl(tiny,x);
      46  	if (__builtin_expect(n> 50000 || k+n > 0x7ffe, 0))
      47  	  return huge*copysignl(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  	    {SET_LDOUBLE_EXP(x,(es&0x8000)|k); return x;}
      53  	if (k <= -64)
      54  	    return tiny*copysignl(tiny,x); 	/*underflow*/
      55  	k += 64;				/* subnormal result */
      56  	SET_LDOUBLE_EXP(x,(es&0x8000)|k);
      57  	return x*twom64;
      58  }