sysdeps/ieee754/dbl-64/wordsize-64/s_nearbyint.c

   1 /* Adapted for use as nearbyint by Ulrich Drepper <drepper@cygnus.com>.  */
   2 /*
   3  * ====================================================
   4  * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
   5  *
   6  * Developed at SunPro, a Sun Microsystems, Inc. business.
   7  * Permission to use, copy, modify, and distribute this
   8  * software is freely granted, provided that this notice
   9  * is preserved.
  10  * ====================================================
  11  */
  12
  13 /*
  14  * rint(x)
  15  * Return x rounded to integral value according to the prevailing
  16  * rounding mode.
  17  * Method:
  18  *      Using floating addition.
  19  * Exception:
  20  *      Inexact flag raised if x not equal to rint(x).
  21  */
  22
  23 #include <fenv.h>
  24 #include <math.h>
  25 #include <math_private.h>
  26
  27 static const double
  28 TWO52[2]={
  29   4.50359962737049600000e+15, /* 0x43300000, 0x00000000 */
  30  -4.50359962737049600000e+15, /* 0xC3300000, 0x00000000 */
  31 };
  32
  33 double
  34 __nearbyint(double x)
  35 {
  36         fenv_t env;
  37         int64_t i0,sx;
  38         int32_t j0;
  39         EXTRACT_WORDS64(i0,x);
  40         sx = (i0>>63)&1;
  41         j0 = ((i0>>52)&0x7ff)-0x3ff;
  42         if(__builtin_expect(j0<52, 1)) {
  43             if(j0<0) {
  44                 libc_feholdexcept (&env);
  45                 double w = TWO52[sx]+x;
  46                 double t =  w-TWO52[sx];
  47                 math_opt_barrier(t);
  48                 libc_fesetenv (&env);
  49                 return copysign(t, x);
  50             }
  51         } else {
  52             if(j0==0x400) return x+x;   /* inf or NaN */
  53             else return x;              /* x is integral */
  54         }
  55         libc_feholdexcept (&env);
  56         double w = TWO52[sx]+x;
  57         double t = w-TWO52[sx];
  58         math_opt_barrier (t);
  59         libc_fesetenv (&env);
  60         return t;
  61 }
  62 weak_alias (__nearbyint, nearbyint)
  63 #ifdef NO_LONG_DOUBLE
  64 strong_alias (__nearbyint, __nearbyintl)
  65 weak_alias (__nearbyint, nearbyintl)
  66 #endif