1 /* Copyright (C) 2007-2023 Free Software Foundation, Inc.
2
3 This file is part of GCC.
4
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 3, or (at your option) any later
8 version.
9
10 GCC is distributed in the hope that it will be useful, but WITHOUT ANY
11 WARRANTY; without even the implied warranty of MERCHANTABILITY or
12 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
13 for more details.
14
15 Under Section 7 of GPL version 3, you are granted additional
16 permissions described in the GCC Runtime Library Exception, version
17 3.1, as published by the Free Software Foundation.
18
19 You should have received a copy of the GNU General Public License and
20 a copy of the GCC Runtime Library Exception along with this program;
21 see the files COPYING3 and COPYING.RUNTIME respectively. If not, see
22 <http://www.gnu.org/licenses/>. */
23
24 #define BID_128RES
25 #include "bid_div_macros.h"
26
27
28 BID128_FUNCTION_ARG2_NORND_CUSTOMRESTYPE (UINT128, bid128_rem, x, y)
29
30 UINT256 P256;
31 UINT128 CX, CY, CX2, CQ, CR, T, CXS, P128, res;
32 UINT64 sign_x, sign_y, valid_y;
33 SINT64 D;
34 int_float f64, fx;
35 int exponent_x, exponent_y, diff_expon, bin_expon_cx, scale,
36 scale0;
37
38 // unpack arguments, check for NaN or Infinity
39
40 valid_y = unpack_BID128_value (&sign_y, &exponent_y, &CY, y);
41
42 if (!unpack_BID128_value (&sign_x, &exponent_x, &CX, x)) {
43 #ifdef SET_STATUS_FLAGS
44 if ((y.w[1] & SNAN_MASK64) == SNAN_MASK64) // y is sNaN
45 __set_status_flags (pfpsf, INVALID_EXCEPTION);
46 #endif
47 // test if x is NaN
48 if ((x.w[1] & 0x7c00000000000000ull) == 0x7c00000000000000ull) {
49 #ifdef SET_STATUS_FLAGS
50 if ((x.w[1] & SNAN_MASK64) == SNAN_MASK64) // y is sNaN
51 __set_status_flags (pfpsf, INVALID_EXCEPTION);
52 #endif
53 res.w[1] = CX.w[1] & QUIET_MASK64;
54 res.w[0] = CX.w[0];
55 BID_RETURN (res);
56 }
57 // x is Infinity?
58 if ((x.w[1] & 0x7800000000000000ull) == 0x7800000000000000ull) {
59 // check if y is Inf.
60 if (((y.w[1] & 0x7c00000000000000ull) != 0x7c00000000000000ull))
61 // return NaN
62 {
63 #ifdef SET_STATUS_FLAGS
64 // set status flags
65 __set_status_flags (pfpsf, INVALID_EXCEPTION);
66 #endif
67 res.w[1] = 0x7c00000000000000ull;
68 res.w[0] = 0;
69 BID_RETURN (res);
70 }
71
72 }
73 // x is 0
74 if ((!CY.w[1]) && (!CY.w[0])) {
75 #ifdef SET_STATUS_FLAGS
76 // set status flags
77 __set_status_flags (pfpsf, INVALID_EXCEPTION);
78 #endif
79 // x=y=0, return NaN
80 res.w[1] = 0x7c00000000000000ull;
81 res.w[0] = 0;
82 BID_RETURN (res);
83 }
84 if (valid_y || ((y.w[1] & NAN_MASK64) == INFINITY_MASK64)) {
85 // return 0
86 if ((exponent_x > exponent_y)
87 && ((y.w[1] & NAN_MASK64) != INFINITY_MASK64))
88 exponent_x = exponent_y;
89
90 res.w[1] = sign_x | (((UINT64) exponent_x) << 49);
91 res.w[0] = 0;
92 BID_RETURN (res);
93 }
94 }
95 if (!valid_y) {
96 // y is Inf. or NaN
97
98 // test if y is NaN
99 if ((y.w[1] & 0x7c00000000000000ull) == 0x7c00000000000000ull) {
100 #ifdef SET_STATUS_FLAGS
101 if ((y.w[1] & SNAN_MASK64) == SNAN_MASK64) // y is sNaN
102 __set_status_flags (pfpsf, INVALID_EXCEPTION);
103 #endif
104 res.w[1] = CY.w[1] & QUIET_MASK64;
105 res.w[0] = CY.w[0];
106 BID_RETURN (res);
107 }
108 // y is Infinity?
109 if ((y.w[1] & 0x7800000000000000ull) == 0x7800000000000000ull) {
110 // return x
111 res.w[1] = x.w[1];
112 res.w[0] = x.w[0];
113 BID_RETURN (res);
114 }
115 // y is 0
116 #ifdef SET_STATUS_FLAGS
117 // set status flags
118 __set_status_flags (pfpsf, INVALID_EXCEPTION);
119 #endif
120 res.w[1] = 0x7c00000000000000ull;
121 res.w[0] = 0;
122 BID_RETURN (res);
123 }
124
125 diff_expon = exponent_x - exponent_y;
126
127 if (diff_expon <= 0) {
128 diff_expon = -diff_expon;
129
130 if (diff_expon > 34) {
131 // |x|<|y| in this case
132 res = x;
133 BID_RETURN (res);
134 }
135 // set exponent of y to exponent_x, scale coefficient_y
136 T = power10_table_128[diff_expon];
137 __mul_128x128_to_256 (P256, CY, T);
138
139 if (P256.w[2] || P256.w[3]) {
140 // |x|<|y| in this case
141 res = x;
142 BID_RETURN (res);
143 }
144
145 CX2.w[1] = (CX.w[1] << 1) | (CX.w[0] >> 63);
146 CX2.w[0] = CX.w[0] << 1;
147 if (__unsigned_compare_ge_128 (P256, CX2)) {
148 // |x|<|y| in this case
149 res = x;
150 BID_RETURN (res);
151 }
152
153 P128.w[0] = P256.w[0];
154 P128.w[1] = P256.w[1];
155 __div_128_by_128 (&CQ, &CR, CX, P128);
156
157 CX2.w[1] = (CR.w[1] << 1) | (CR.w[0] >> 63);
158 CX2.w[0] = CR.w[0] << 1;
159 if ((__unsigned_compare_gt_128 (CX2, P256))
160 || (CX2.w[1] == P256.w[1] && CX2.w[0] == P256.w[0]
161 && (CQ.w[0] & 1))) {
162 __sub_128_128 (CR, P256, CR);
163 sign_x ^= 0x8000000000000000ull;
164 }
165
166 get_BID128_very_fast (&res, sign_x, exponent_x, CR);
167 BID_RETURN (res);
168 }
169 // 2^64
170 f64.i = 0x5f800000;
171
172 scale0 = 38;
173 if (!CY.w[1])
174 scale0 = 34;
175
176 while (diff_expon > 0) {
177 // get number of digits in CX and scale=38-digits
178 // fx ~ CX
179 fx.d = (float) CX.w[1] * f64.d + (float) CX.w[0];
180 bin_expon_cx = ((fx.i >> 23) & 0xff) - 0x7f;
181 scale = scale0 - estimate_decimal_digits[bin_expon_cx];
182 // scale = 38-estimate_decimal_digits[bin_expon_cx];
183 D = CX.w[1] - power10_index_binexp_128[bin_expon_cx].w[1];
184 if (D > 0
185 || (!D && CX.w[0] >= power10_index_binexp_128[bin_expon_cx].w[0]))
186 scale--;
187
188 if (diff_expon >= scale)
189 diff_expon -= scale;
190 else {
191 scale = diff_expon;
192 diff_expon = 0;
193 }
194
195 T = power10_table_128[scale];
196 __mul_128x128_low (CXS, CX, T);
197
198 __div_128_by_128 (&CQ, &CX, CXS, CY);
199
200 // check for remainder == 0
201 if (!CX.w[1] && !CX.w[0]) {
202 get_BID128_very_fast (&res, sign_x, exponent_y, CX);
203 BID_RETURN (res);
204 }
205 }
206
207 CX2.w[1] = (CX.w[1] << 1) | (CX.w[0] >> 63);
208 CX2.w[0] = CX.w[0] << 1;
209 if ((__unsigned_compare_gt_128 (CX2, CY))
210 || (CX2.w[1] == CY.w[1] && CX2.w[0] == CY.w[0] && (CQ.w[0] & 1))) {
211 __sub_128_128 (CX, CY, CX);
212 sign_x ^= 0x8000000000000000ull;
213 }
214
215 get_BID128_very_fast (&res, sign_x, exponent_y, CX);
216 BID_RETURN (res);
217 }