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 FUNCTION
13 <<round>>, <<roundf>>---round to integer, to nearest
14 INDEX
15 round
16 INDEX
17 roundf
18
19 SYNOPSIS
20 #include <math.h>
21 double round(double <[x]>);
22 float roundf(float <[x]>);
23
24 DESCRIPTION
25 The <<round>> functions round their argument to the nearest integer
26 value in floating-point format, rounding halfway cases away from zero,
27 regardless of the current rounding direction. (While the "inexact"
28 floating-point exception behavior is unspecified by the C standard, the
29 <<round>> functions are written so that "inexact" is not raised if the
30 result does not equal the argument, which behavior is as recommended by
31 IEEE 754 for its related functions.)
32
33 RETURNS
34 <[x]> rounded to an integral value.
35
36 PORTABILITY
37 ANSI C, POSIX
38
39 SEEALSO
40 <<nearbyint>>, <<rint>>
41
42 */
43
44 #include "fdlibm.h"
45
46 #ifdef _NEED_FLOAT64
47
48 __float64
round64(__float64 x)49 round64(__float64 x)
50 {
51 /* Most significant word, least significant word. */
52 __int32_t msw, exponent_less_1023;
53 __uint32_t lsw;
54
55 EXTRACT_WORDS(msw, lsw, x);
56
57 /* Extract exponent field. */
58 exponent_less_1023 = ((msw & 0x7ff00000) >> 20) - 1023;
59
60 if (exponent_less_1023 < 20)
61 {
62 if (exponent_less_1023 < 0)
63 {
64 msw &= 0x80000000;
65 if (exponent_less_1023 == -1)
66 /* Result is +1.0 or -1.0. */
67 msw |= ((__int32_t)1023 << 20);
68 lsw = 0;
69 }
70 else
71 {
72 __uint32_t exponent_mask = 0x000fffff >> exponent_less_1023;
73 if ((msw & exponent_mask) == 0 && lsw == 0)
74 /* x in an integral value. */
75 return x;
76
77 msw += 0x00080000 >> exponent_less_1023;
78 msw &= ~exponent_mask;
79 lsw = 0;
80 }
81 }
82 else if (exponent_less_1023 > 51)
83 {
84 if (exponent_less_1023 == 1024)
85 /* x is NaN or infinite. */
86 return x + x;
87 else
88 return x;
89 }
90 else
91 {
92 __uint32_t exponent_mask = 0xffffffff >> (exponent_less_1023 - 20);
93 __uint32_t tmp;
94
95 if ((lsw & exponent_mask) == 0)
96 /* x is an integral value. */
97 return x;
98
99 tmp = lsw + (1 << (51 - exponent_less_1023));
100 if (tmp < lsw)
101 msw += 1;
102 lsw = tmp;
103
104 lsw &= ~exponent_mask;
105 }
106 INSERT_WORDS(x, msw, lsw);
107
108 return x;
109 }
110
111 _MATH_ALIAS_d_d(round)
112
113 #endif /* _NEED_FLOAT64 */
114