1 /* Copyright (C) 2007 Free Software Foundation, Inc.
3 This file is part of GCC.
5 GCC is free software; you can redistribute it and/or modify it under
6 the terms of the GNU General Public License as published by the Free
7 Software Foundation; either version 2, or (at your option) any later
10 In addition to the permissions in the GNU General Public License, the
11 Free Software Foundation gives you unlimited permission to link the
12 compiled version of this file into combinations with other programs,
13 and to distribute those combinations without any restriction coming
14 from the use of this file. (The General Public License restrictions
15 do apply in other respects; for example, they cover modification of
16 the file, and distribution when not linked into a combine
19 GCC is distributed in the hope that it will be useful, but WITHOUT ANY
20 WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
24 You should have received a copy of the GNU General Public License
25 along with GCC; see the file COPYING. If not, write to the Free
26 Software Foundation, 51 Franklin Street, Fifth Floor, Boston, MA
30 #include "bid_internal.h"
32 BID128_FUNCTION_ARG2(__bid128_quantize
, x
, y
)
35 UINT128 CX
, CY
, T
, CX2
, CR
, Stemp
, res
, REM_H
, C2N
;
36 UINT64 sign_x
, sign_y
, remainder_h
, carry
, CY64
;
38 int exponent_x
= 0, exponent_y
= 0, digits_x
, extra_digits
, amount
;
39 int expon_diff
, total_digits
, bin_expon_cx
, rmode
, status
;
41 // unpack arguments, check for NaN or Infinity
42 if (!unpack_BID128_value (&sign_y
, &exponent_y
, &CY
, y
)) {
44 #ifdef SET_STATUS_FLAGS
45 if ((x
.w
[1] & SNAN_MASK64
) == SNAN_MASK64
) // y is sNaN
46 __set_status_flags (pfpsf
, INVALID_EXCEPTION
);
50 if ((y
.w
[1] & 0x7c00000000000000ull
) == 0x7c00000000000000ull
) {
51 #ifdef SET_STATUS_FLAGS
52 if ((y
.w
[1] & 0x7e00000000000000ull
) == 0x7e00000000000000ull
) {
54 __set_status_flags (pfpsf
, INVALID_EXCEPTION
);
57 res
.w
[1] = y
.w
[1] & QUIET_MASK64
;
62 if ((y
.w
[1] & 0x7800000000000000ull
) == 0x7800000000000000ull
) {
63 // check if x is not Inf.
64 if (((x
.w
[1] & 0x7c00000000000000ull
) < 0x7800000000000000ull
)) {
66 #ifdef SET_STATUS_FLAGS
68 __set_status_flags (pfpsf
, INVALID_EXCEPTION
);
70 res
.w
[1] = 0x7c00000000000000ull
;
73 } else if (((x
.w
[1] & 0x7c00000000000000ull
) <= 0x7800000000000000ull
)) {
81 if (!unpack_BID128_value (&sign_x
, &exponent_x
, &CX
, x
)) {
82 // test if x is NaN or Inf
83 if ((x
.w
[1] & 0x7c00000000000000ull
) == 0x7800000000000000ull
) {
84 #ifdef SET_STATUS_FLAGS
86 __set_status_flags (pfpsf
, INVALID_EXCEPTION
);
88 res
.w
[1] = 0x7c00000000000000ull
;
91 } else if ((x
.w
[1] & 0x7c00000000000000ull
) ==
92 0x7c00000000000000ull
) {
93 if ((x
.w
[1] & 0x7e00000000000000ull
) == 0x7e00000000000000ull
) {
94 #ifdef SET_STATUS_FLAGS
96 __set_status_flags (pfpsf
, INVALID_EXCEPTION
);
99 res
.w
[1] = x
.w
[1] & QUIET_MASK64
;
103 if (!CX
.w
[1] && !CX
.w
[0]) {
104 get_BID128_very_fast (&res
, sign_x
, exponent_y
, CX
);
108 // get number of decimal digits in coefficient_x
110 tempx
.d
= (float) CX
.w
[1];
111 bin_expon_cx
= ((tempx
.i
>> 23) & 0xff) - 0x7f + 64;
113 tempx
.d
= (float) CX
.w
[0];
114 bin_expon_cx
= ((tempx
.i
>> 23) & 0xff) - 0x7f;
116 digits_x
= __bid_estimate_decimal_digits
[bin_expon_cx
];
117 if (CX
.w
[1] > __bid_power10_table_128
[digits_x
].w
[1]
118 || (CX
.w
[1] == __bid_power10_table_128
[digits_x
].w
[1]
119 && CX
.w
[0] >= __bid_power10_table_128
[digits_x
].w
[0]))
122 expon_diff
= exponent_x
- exponent_y
;
123 total_digits
= digits_x
+ expon_diff
;
125 if ((UINT32
) total_digits
<= 34) {
126 if (expon_diff
>= 0) {
127 T
= __bid_power10_table_128
[expon_diff
];
128 __mul_128x128_low (CX2
, T
, CX
);
129 get_BID128_very_fast (&res
, sign_x
, exponent_y
, CX2
);
132 #ifndef IEEE_ROUND_NEAREST_TIES_AWAY
133 #ifndef IEEE_ROUND_NEAREST
135 if (sign_x
&& (unsigned) (rmode
- 1) < 2)
143 // must round off -expon_diff digits
144 extra_digits
= -expon_diff
;
145 __add_128_128 (CX
, CX
, __bid_round_const_table_128
[rmode
][extra_digits
]);
147 // get P*(2^M[extra_digits])/10^extra_digits
148 __mul_128x128_to_256 (CT
, CX
, __bid_reciprocals10_128
[extra_digits
]);
150 // now get P/10^extra_digits: shift C64 right by M[extra_digits]-128
151 amount
= __bid_recip_scale
[extra_digits
];
156 CR
.w
[0] = CX2
.w
[1] >> (amount
- 64);
158 __shr_128 (CR
, CX2
, amount
);
161 #ifndef IEEE_ROUND_NEAREST_TIES_AWAY
162 #ifndef IEEE_ROUND_NEAREST
166 // check whether fractional part of initial_P/10^extra_digits is
167 // exactly .5 this is the same as fractional part of
168 // (initial_P + 0.5*10^extra_digits)/10^extra_digits is exactly zero
172 remainder_h
= CX2
.w
[0] | (CX2
.w
[1] << (128 - amount
));
174 remainder_h
= CX2
.w
[0] << (64 - amount
);
176 // test whether fractional part is 0
178 && (CT
.w
[1] < __bid_reciprocals10_128
[extra_digits
].w
[1]
179 || (CT
.w
[1] == __bid_reciprocals10_128
[extra_digits
].w
[1]
180 && CT
.w
[0] < __bid_reciprocals10_128
[extra_digits
].w
[0]))) {
186 #ifdef SET_STATUS_FLAGS
187 status
= INEXACT_EXCEPTION
;
191 REM_H
.w
[1] = (CX2
.w
[1] << (128 - amount
));
192 REM_H
.w
[0] = CX2
.w
[0];
194 REM_H
.w
[1] = CX2
.w
[0] << (64 - amount
);
199 case ROUNDING_TO_NEAREST
:
200 case ROUNDING_TIES_AWAY
:
201 // test whether fractional part is 0
202 if (REM_H
.w
[1] == 0x8000000000000000ull
&& !REM_H
.w
[0]
203 && (CT
.w
[1] < __bid_reciprocals10_128
[extra_digits
].w
[1]
204 || (CT
.w
[1] == __bid_reciprocals10_128
[extra_digits
].w
[1]
205 && CT
.w
[0] < __bid_reciprocals10_128
[extra_digits
].w
[0])))
206 status
= EXACT_STATUS
;
209 case ROUNDING_TO_ZERO
:
210 if (!(REM_H
.w
[1] | REM_H
.w
[0])
211 && (CT
.w
[1] < __bid_reciprocals10_128
[extra_digits
].w
[1]
212 || (CT
.w
[1] == __bid_reciprocals10_128
[extra_digits
].w
[1]
213 && CT
.w
[0] < __bid_reciprocals10_128
[extra_digits
].w
[0])))
214 status
= EXACT_STATUS
;
218 __add_carry_out (Stemp
.w
[0], CY64
, CT
.w
[0],
219 __bid_reciprocals10_128
[extra_digits
].w
[0]);
220 __add_carry_in_out (Stemp
.w
[1], carry
, CT
.w
[1],
221 __bid_reciprocals10_128
[extra_digits
].w
[1], CY64
);
224 C2N
.w
[0] = ((UINT64
) 1) << amount
;
225 REM_H
.w
[0] = REM_H
.w
[1] >> (64 - amount
);
228 C2N
.w
[1] = ((UINT64
) 1) << (amount
- 64);
230 REM_H
.w
[1] >>= (128 - amount
);
233 if (REM_H
.w
[0] < carry
)
235 if (__unsigned_compare_ge_128 (REM_H
, C2N
))
236 status
= EXACT_STATUS
;
239 __set_status_flags (pfpsf
, status
);
242 get_BID128_very_fast (&res
, sign_x
, exponent_y
, CR
);
245 if (total_digits
< 0) {
246 CR
.w
[1] = CR
.w
[0] = 0;
247 #ifndef IEEE_ROUND_NEAREST_TIES_AWAY
248 #ifndef IEEE_ROUND_NEAREST
250 if (sign_x
&& (unsigned) (rmode
- 1) < 2)
252 if (rmode
== ROUNDING_UP
)
256 #ifdef SET_STATUS_FLAGS
257 __set_status_flags (pfpsf
, INEXACT_EXCEPTION
);
259 get_BID128_very_fast (&res
, sign_x
, exponent_y
, CR
);
262 // else more than 34 digits in coefficient
263 #ifdef SET_STATUS_FLAGS
264 __set_status_flags (pfpsf
, INVALID_EXCEPTION
);
266 res
.w
[1] = 0x7c00000000000000ull
;