/ base / math / logb.c
logb.c
  1  #if defined(__ppc__)
  2  /*******************************************************************************
  3  *                                                                              *
  4  *      File logb.c,                                                            *
  5  *      Functions logb.                                                         *
  6  *      Implementation of logb for the PowerPC.                                 *
  7  *                                                                              *
  8  *      Copyright � 1991 Apple Computer, Inc.  All rights reserved.             *
  9  *                                                                              *
 10  *      Written by Ali Sazegari, started on June 1991,                          *
 11  *                                                                              *
 12  *      August   26 1991: removed CFront Version 1.1d17 warnings.               *
 13  *      August   27 1991: no errors reported by the test suite.                 *
 14  *      November 11 1991: changed CLASSEXTENDED to the macro CLASSIFY and       *
 15  *                        + or - infinity to constants.                         *
 16  *      November 18 1991: changed the macro CLASSIFY to CLASSEXTENDEDint to     *
 17  *                        improve performance.                                  *
 18  *      February 07 1992: changed bit operations to macros (  object size is    *
 19  *                        unchanged  ).                                         *
 20  *      September24 1992: took the "#include support.h" out.                    *
 21  *      December 03 1992: first rs/6000 port.                                   *
 22  *      August   30 1992: set the divide by zero for the zero argument case.    *
 23  *      October  05 1993: corrected the environment.                            *
 24  *      October  17 1994: replaced all environmental functions with __setflm.   *
 25  *      May      28 1997: made speed improvements.                              *
 26  *      April    30 2001: forst mac os x port using gcc.                        *
 27  *                                                                              *
 28  ********************************************************************************
 29  *     The C math library offers a similar function called "frexp".  It is      *
 30  *     different in details from logb, but similar in spirit.  This current     *
 31  *     implementation of logb follows the recommendation in IEEE Standard 854   *
 32  *     which is different in its handling of denormalized numbers from the IEEE *
 33  *     Standard 754.                                                            *
 34  *******************************************************************************/
 35  
 36  typedef union           
 37        { 
 38        struct {
 39  #if defined(__BIG_ENDIAN__)
 40          unsigned long int hi;
 41          unsigned long int lo;
 42  #else
 43          unsigned long int lo;
 44          unsigned long int hi;
 45  #endif
 46        } words;
 47        double dbl;
 48        } DblInHex;
 49  
 50  static const double twoTo52 = 4.50359962737049600e15;              // 0x1p52
 51  static const double klTod = 4503601774854144.0;                    // 0x1.000008p52
 52  static const unsigned long int signMask = 0x80000000ul;
 53  static const DblInHex minusInf  = {{ 0xFFF00000, 0x00000000 }};
 54  
 55  
 56  /*******************************************************************************
 57  ********************************************************************************
 58  *                                    L  O  G  B                                *
 59  ********************************************************************************
 60  *******************************************************************************/
 61  
 62  double logb (  double x  )
 63        {
 64        DblInHex xInHex;
 65        long int shiftedExp;
 66        
 67        xInHex.dbl = x;
 68        shiftedExp = ( xInHex.words.hi & 0x7ff00000UL ) >> 20;
 69        
 70        if ( shiftedExp == 2047 ) 
 71              {                                            // NaN or INF
 72              if ( ( ( xInHex.words.hi & signMask ) == 0 ) || ( x != x ) )
 73                    return x;                              // NaN or +INF return x
 74              else
 75                    return -x;                             // -INF returns +INF
 76              }
 77        
 78        if ( shiftedExp != 0 )                             // normal number
 79              shiftedExp -= 1023;                          // unbias exponent
 80        
 81        else if ( x == 0.0 ) 
 82              {                                            // zero
 83              xInHex.words.hi = 0x0UL;                      // return -infinity
 84              return (  minusInf.dbl  );
 85              }
 86        
 87        else 
 88              {                                            // subnormal number
 89              xInHex.dbl *= twoTo52;                       // scale up
 90              shiftedExp = ( xInHex.words.hi & 0x7ff00000UL ) >> 20;
 91              shiftedExp -= 1075;                          // unbias exponent
 92              }
 93        
 94        if ( shiftedExp == 0 )                             // zero result
 95              return ( 0.0 );
 96        
 97        else 
 98              {                                            // nonzero result
 99              xInHex.dbl = klTod;
100              xInHex.words.lo += shiftedExp;
101              return ( xInHex.dbl - klTod );
102              }
103        }
104  #endif /* __ppc__ */