1/*
2 * Copyright (c) 1998, 2001, Oracle and/or its affiliates. All rights reserved.
3 * DO NOT ALTER OR REMOVE COPYRIGHT NOTICES OR THIS FILE HEADER.
4 *
5 * This code is free software; you can redistribute it and/or modify it
6 * under the terms of the GNU General Public License version 2 only, as
7 * published by the Free Software Foundation. Oracle designates this
8 * particular file as subject to the "Classpath" exception as provided
9 * by Oracle in the LICENSE file that accompanied this code.
10 *
11 * This code is distributed in the hope that it will be useful, but WITHOUT
12 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
13 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
14 * version 2 for more details (a copy is included in the LICENSE file that
15 * accompanied this code).
16 *
17 * You should have received a copy of the GNU General Public License version
18 * 2 along with this work; if not, write to the Free Software Foundation,
19 * Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA.
20 *
21 * Please contact Oracle, 500 Oracle Parkway, Redwood Shores, CA 94065 USA
22 * or visit www.oracle.com if you need additional information or have any
23 * questions.
24 */
25
26/*
27 * rint(x)
28 * Return x rounded to integral value according to the prevailing
29 * rounding mode.
30 * Method:
31 * Using floating addition.
32 * Exception:
33 * Inexact flag raised if x not equal to rint(x).
34 */
35
36#include "fdlibm.h"
37
38#ifdef __STDC__
39static const double
40#else
41static double
42#endif
43TWO52[2]={
44 4.50359962737049600000e+15, /* 0x43300000, 0x00000000 */
45 -4.50359962737049600000e+15, /* 0xC3300000, 0x00000000 */
46};
47
48#ifdef __STDC__
49 double rint(double x)
50#else
51 double rint(x)
52 double x;
53#endif
54{
55 int i0,j0,sx;
56 unsigned i,i1;
57 double w,t;
58 i0 = __HI(x);
59 sx = (i0>>31)&1;
60 i1 = __LO(x);
61 j0 = ((i0>>20)&0x7ff)-0x3ff;
62 if(j0<20) {
63 if(j0<0) {
64 if(((i0&0x7fffffff)|i1)==0) return x;
65 i1 |= (i0&0x0fffff);
66 i0 &= 0xfffe0000;
67 i0 |= ((i1|-i1)>>12)&0x80000;
68 __HI(x)=i0;
69 w = TWO52[sx]+x;
70 t = w-TWO52[sx];
71 i0 = __HI(t);
72 __HI(t) = (i0&0x7fffffff)|(sx<<31);
73 return t;
74 } else {
75 i = (0x000fffff)>>j0;
76 if(((i0&i)|i1)==0) return x; /* x is integral */
77 i>>=1;
78 if(((i0&i)|i1)!=0) {
79 if(j0==19) i1 = 0x40000000; else
80 i0 = (i0&(~i))|((0x20000)>>j0);
81 }
82 }
83 } else if (j0>51) {
84 if(j0==0x400) return x+x; /* inf or NaN */
85 else return x; /* x is integral */
86 } else {
87 i = ((unsigned)(0xffffffff))>>(j0-20);
88 if((i1&i)==0) return x; /* x is integral */
89 i>>=1;
90 if((i1&i)!=0) i1 = (i1&(~i))|((0x40000000)>>(j0-20));
91 }
92 __HI(x) = i0;
93 __LO(x) = i1;
94 w = TWO52[sx]+x;
95 return w-TWO52[sx];
96}
97