toolchain/gcc/newlib/libm/common/s_round.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 #include "fdlibm.h"
  13
  14 #ifndef _DOUBLE_IS_32BITS
  15
  16 #ifdef __STDC__
  17         double round(double x)
  18 #else
  19         double round(x)
  20         double x;
  21 #endif
  22 {
  23   /* Most significant word, least significant word. */
  24   __int32_t msw, exponent_less_1023;
  25   __uint32_t lsw;
  26
  27   EXTRACT_WORDS(msw, lsw, x);
  28
  29   /* Extract exponent field. */
  30   exponent_less_1023 = ((msw & 0x7ff00000) >> 20) - 1023;
  31
  32   if (exponent_less_1023 < 20)
  33     {
  34       if (exponent_less_1023 < 0)
  35         {
  36           msw &= 0x80000000;
  37           if (exponent_less_1023 == -1)
  38             /* Result is +1.0 or -1.0. */
  39             msw |= (1023 << 20);
  40           lsw = 0;
  41         }
  42       else
  43         {
  44           __uint32_t exponent_mask = 0x000fffff >> exponent_less_1023;
  45           if ((msw & exponent_mask) == 0 && lsw == 0)
  46             /* x in an integral value. */
  47             return x;
  48
  49           msw += 0x00080000 >> exponent_less_1023;
  50           msw &= ~exponent_mask;
  51           lsw = 0;
  52         }
  53     }
  54   else if (exponent_less_1023 > 51)
  55     {
  56       if (exponent_less_1023 == 1024)
  57         /* x is NaN or infinite. */
  58         return x + x;
  59       else
  60         return x;
  61     }
  62   else
  63     {
  64       __uint32_t exponent_mask = 0xffffffff >> (exponent_less_1023 - 20);
  65       __uint32_t tmp;
  66
  67       if ((lsw & exponent_mask) == 0)
  68         /* x is an integral value. */
  69         return x;
  70
  71       tmp = lsw + (1 << (51 - exponent_less_1023));
  72       if (tmp < lsw)
  73         msw += 1;
  74       lsw = tmp;
  75
  76       lsw &= ~exponent_mask;
  77     }
  78   INSERT_WORDS(x, msw, lsw);
  79
  80   return x;
  81 }
  82
  83 #endif /* _DOUBLE_IS_32BITS */