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 <float.h> 23 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 rint(double x) 35 { 36 int32_t i0,j0,sx; 37 u_int32_t i,i1; 38 double w,t; 39 EXTRACT_WORDS(i0,i1,x); 40 sx = (i0>>31)&1; 41 j0 = ((i0>>20)&0x7ff)-0x3ff; 42 if(j0<20) { 43 if(j0<0) { 44 if(((i0&0x7fffffff)|i1)==0) return x; 45 i1 |= (i0&0x0fffff); 46 i0 &= 0xfffe0000; 47 i0 |= ((i1|-i1)>>12)&0x80000; 48 SET_HIGH_WORD(x,i0); 49 STRICT_ASSIGN(double,w,TWO52[sx]+x); 50 t = w-TWO52[sx]; 51 GET_HIGH_WORD(i0,t); 52 SET_HIGH_WORD(t,(i0&0x7fffffff)|(sx<<31)); 53 return t; 54 } else { 55 i = (0x000fffff)>>j0; 56 if(((i0&i)|i1)==0) return x; /* x is integral */ 57 i>>=1; 58 if(((i0&i)|i1)!=0) { 59 /* 60 * Some bit is set after the 0.5 bit. To avoid the 61 * possibility of errors from double rounding in 62 * w = TWO52[sx]+x, adjust the 0.25 bit to a lower 63 * guard bit. We do this for all j0<=51. The 64 * adjustment is trickiest for j0==18 and j0==19 65 * since then it spans the word boundary. 66 */ 67 if(j0==19) i1 = 0x40000000; else 68 if(j0==18) i1 = 0x80000000; else 69 i0 = (i0&(~i))|((0x20000)>>j0); 70 } 71 } 72 } else if (j0>51) { 73 if(j0==0x400) return x+x; /* inf or NaN */ 74 else return x; /* x is integral */ 75 } else { 76 i = ((u_int32_t)(0xffffffff))>>(j0-20); 77 if((i1&i)==0) return x; /* x is integral */ 78 i>>=1; 79 if((i1&i)!=0) i1 = (i1&(~i))|((0x40000000)>>(j0-20)); 80 } 81 INSERT_WORDS(x,i0,i1); 82 STRICT_ASSIGN(double,w,TWO52[sx]+x); 83 return w-TWO52[sx]; 84 } 85 86 #if (LDBL_MANT_DIG == 53) 87 __weak_reference(rint, rintl); 88 #endif 89