sysdeps/ieee754/ldbl-96/s_rintl.c

   1 /* s_rintl.c -- long double version of s_rint.c.
   2  * Conversion to long double by Ulrich Drepper,
   3  * Cygnus Support, drepper@cygnus.com.
   4  */
   5
   6 /*
   7  * ====================================================
   8  * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
   9  *
  10  * Developed at SunPro, a Sun Microsystems, Inc. business.
  11  * Permission to use, copy, modify, and distribute this
  12  * software is freely granted, provided that this notice
  13  * is preserved.
  14  * ====================================================
  15  */
  16
  17 #if defined(LIBM_SCCS) && !defined(lint)
  18 static char rcsid[] = "$NetBSD: $";
  19 #endif
  20
  21 /*
  22  * rintl(x)
  23  * Return x rounded to integral value according to the prevailing
  24  * rounding mode.
  25  * Method:
  26  *      Using floating addition.
  27  * Exception:
  28  *      Inexact flag raised if x not equal to rintl(x).
  29  */
  30
  31 #include "math.h"
  32 #include "math_private.h"
  33
  34 #ifdef __STDC__
  35 static const long double
  36 #else
  37 static long double
  38 #endif
  39 TWO63[2]={
  40   9.223372036854775808000000e+18, /* 0x403E, 0x00000000, 0x00000000 */
  41  -9.223372036854775808000000e+18  /* 0xC03E, 0x00000000, 0x00000000 */
  42 };
  43
  44 #ifdef __STDC__
  45         long double __rintl(long double x)
  46 #else
  47         long double __rintl(x)
  48         long double x;
  49 #endif
  50 {
  51         int32_t se,j0,sx;
  52         u_int32_t i,i0,i1;
  53         long double w,t;
  54         GET_LDOUBLE_WORDS(se,i0,i1,x);
  55         sx = (se>>15)&1;
  56         j0 = (se&0x7fff)-0x3fff;
  57         if(j0<31) {
  58             if(j0<0) {
  59                 if(((se&0x7fff)|i0|i1)==0) return x;
  60                 i1 |= i0;
  61                 i0 &= 0xe0000000;
  62                 i0 |= (i1|-i1)&0x80000000;
  63                 SET_LDOUBLE_MSW(x,i0);
  64                 w = TWO63[sx]+x;
  65                 t = w-TWO63[sx];
  66                 GET_LDOUBLE_EXP(i0,t);
  67                 SET_LDOUBLE_EXP(t,(i0&0x7fff)|(sx<<15));
  68                 return t;
  69             } else {
  70                 i = (0x7fffffff)>>j0;
  71                 if(((i0&i)|i1)==0) return x; /* x is integral */
  72                 i>>=1;
  73                 if(((i0&i)|i1)!=0) {
  74                     if(j0==30) i1 = 0x40000000; else
  75                     i0 = (i0&(~i))|((0x20000000)>>j0);
  76                 }
  77             }
  78         } else if (j0>62) {
  79             if(j0==0x4000) return x+x;  /* inf or NaN */
  80             else return x;              /* x is integral */
  81         } else {
  82             i = ((u_int32_t)(0xffffffff))>>(j0-31);
  83             if((i1&i)==0) return x;     /* x is integral */
  84             i>>=1;
  85             if((i1&i)!=0) i1 = (i1&(~i))|((0x40000000)>>(j0-31));
  86         }
  87         SET_LDOUBLE_WORDS(x,se,i0,i1);
  88         w = TWO63[sx]+x;
  89         return w-TWO63[sx];
  90 }
  91 weak_alias (__rintl, rintl)