lib/ack/math/ldexp.c

   1 /*
   2  * (c) copyright 1987 by the Vrije Universiteit, Amsterdam, The Netherlands.
   3  * See the copyright notice in the ACK home directory, in the file "Copyright".
   4  */
   5 /* $Header$ */
   6
   7 #include        <math.h>
   8 #include        <float.h>
   9 #include        <errno.h>
  10
  11 double
  12 ldexp(double fl, int exp)
  13 {
  14         int sign = 1;
  15         int currexp;
  16
  17         if (__IsNan(fl)) {
  18                 errno = EDOM;
  19                 return fl;
  20         }
  21         if (fl == 0.0) return 0.0;
  22         if (fl<0) {
  23                 fl = -fl;
  24                 sign = -1;
  25         }
  26         if (fl > DBL_MAX) {             /* for infinity */
  27                 errno = ERANGE;
  28                 return sign * fl;
  29         }
  30         fl = frexp(fl,&currexp);
  31         exp += currexp;
  32         if (exp > 0) {
  33                 if (exp > DBL_MAX_EXP) {
  34                         errno = ERANGE;
  35                         return sign * HUGE_VAL;
  36                 }
  37                 while (exp>30) {
  38                         fl *= (double) (1L << 30);
  39                         exp -= 30;
  40                 }
  41                 fl *= (double) (1L << exp);
  42         }
  43         else    {
  44                 /* number need not be normalized */
  45                 if (exp < DBL_MIN_EXP - DBL_MANT_DIG) {
  46                         return 0.0;
  47                 }
  48                 while (exp<-30) {
  49                         fl /= (double) (1L << 30);
  50                         exp += 30;
  51                 }
  52                 fl /= (double) (1L << -exp);
  53         }
  54         return sign * fl;
  55 }