sysdeps/ieee754/dbl-64/s_nearbyint.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  * rint(x)
  14  * Return x rounded to integral value according to the prevailing
  15  * rounding mode.
  16  * Method:
  17  *      Using floating addition.
  18  * Exception:
  19  *      Inexact flag raised if x not equal to rint(x).
  20  */
  21
  22 #include <fenv.h>
  23 #include <math.h>
  24 #include <math-barriers.h>
  25 #include <math_private.h>
  26 #include <fenv_private.h>
  27 #include <libm-alias-double.h>
  28 #include <math-use-builtins.h>
  29
  30 double
  31 __nearbyint (double x)
  32 {
  33 #if USE_NEARBYINT_BUILTIN
  34   return __builtin_nearbyint (x);
  35 #else
  36   /* Use generic implementation.  */
  37   static const double
  38     TWO52[2] = {
  39                 4.50359962737049600000e+15, /* 0x43300000, 0x00000000 */
  40                 -4.50359962737049600000e+15, /* 0xC3300000, 0x00000000 */
  41   };
  42   fenv_t env;
  43   int64_t i0, sx;
  44   int32_t j0;
  45   EXTRACT_WORDS64 (i0, x);
  46   sx = (i0 >> 63) & 1;
  47   j0 = ((i0 >> 52) & 0x7ff) - 0x3ff;
  48   if (__glibc_likely (j0 < 52))
  49     {
  50       if (j0 < 0)
  51         {
  52           libc_feholdexcept (&env);
  53           double w = TWO52[sx] + math_opt_barrier (x);
  54           double t =  w - TWO52[sx];
  55           math_force_eval (t);
  56           libc_fesetenv (&env);
  57           return copysign (t, x);
  58         }
  59     }
  60   else
  61     {
  62       if (j0 == 0x400)
  63         return x + x;                   /* inf or NaN  */
  64       else
  65         return x;                       /* x is integral  */
  66     }
  67   libc_feholdexcept (&env);
  68   double w = TWO52[sx] + math_opt_barrier (x);
  69   double t = w - TWO52[sx];
  70   math_force_eval (t);
  71   libc_fesetenv (&env);
  72   return t;
  73 #endif /* ! USE_NEARBYINT_BUILTIN  */
  74 }
  75 libm_alias_double (__nearbyint, nearbyint)