Back to index

glibc  2.9
s_rintl.c
Go to the documentation of this file.
00001 /* s_rintl.c -- long double version of s_rint.c.
00002  * Conversion to long double by Ulrich Drepper,
00003  * Cygnus Support, drepper@cygnus.com.
00004  */
00005 
00006 /*
00007  * ====================================================
00008  * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
00009  *
00010  * Developed at SunPro, a Sun Microsystems, Inc. business.
00011  * Permission to use, copy, modify, and distribute this
00012  * software is freely granted, provided that this notice
00013  * is preserved.
00014  * ====================================================
00015  */
00016 
00017 #if defined(LIBM_SCCS) && !defined(lint)
00018 static char rcsid[] = "$NetBSD: $";
00019 #endif
00020 
00021 /*
00022  * rintl(x)
00023  * Return x rounded to integral value according to the prevailing
00024  * rounding mode.
00025  * Method:
00026  *     Using floating addition.
00027  * Exception:
00028  *     Inexact flag raised if x not equal to rintl(x).
00029  */
00030 
00031 #include "math.h"
00032 #include "math_private.h"
00033 
00034 #ifdef __STDC__
00035 static const long double
00036 #else
00037 static long double
00038 #endif
00039 TWO63[2]={
00040   9.223372036854775808000000e+18, /* 0x403E, 0x00000000, 0x00000000 */
00041  -9.223372036854775808000000e+18  /* 0xC03E, 0x00000000, 0x00000000 */
00042 };
00043 
00044 #ifdef __STDC__
00045        long double __rintl(long double x)
00046 #else
00047        long double __rintl(x)
00048        long double x;
00049 #endif
00050 {
00051        int32_t se,j0,sx;
00052        u_int32_t i,i0,i1;
00053        long double w,t;
00054        GET_LDOUBLE_WORDS(se,i0,i1,x);
00055        sx = (se>>15)&1;
00056        j0 = (se&0x7fff)-0x3fff;
00057        if(j0<31) {
00058            if(j0<0) {
00059               if(((se&0x7fff)|i0|i1)==0) return x;
00060               i1 |= i0;
00061               i0 &= 0xe0000000;
00062               i0 |= (i1|-i1)&0x80000000;
00063               SET_LDOUBLE_MSW(x,i0);
00064                w = TWO63[sx]+x;
00065                t = w-TWO63[sx];
00066               GET_LDOUBLE_EXP(i0,t);
00067               SET_LDOUBLE_EXP(t,(i0&0x7fff)|(sx<<15));
00068                return t;
00069            } else {
00070               i = (0x7fffffff)>>j0;
00071               if(((i0&i)|i1)==0) return x; /* x is integral */
00072               i>>=1;
00073               if(((i0&i)|i1)!=0) {
00074                   if(j0==30) i1 = 0x40000000; else
00075                   i0 = (i0&(~i))|((0x20000000)>>j0);
00076               }
00077            }
00078        } else if (j0>62) {
00079            if(j0==0x4000) return x+x;     /* inf or NaN */
00080            else return x;          /* x is integral */
00081        } else {
00082            i = ((u_int32_t)(0xffffffff))>>(j0-31);
00083            if((i1&i)==0) return x; /* x is integral */
00084            i>>=1;
00085            if((i1&i)!=0) i1 = (i1&(~i))|((0x40000000)>>(j0-31));
00086        }
00087        SET_LDOUBLE_WORDS(x,se,i0,i1);
00088        w = TWO63[sx]+x;
00089        return w-TWO63[sx];
00090 }
00091 weak_alias (__rintl, rintl)