gas/
[binutils-gdb.git] / libdecnumber / bid / bid2dpd_dpd2bid.c
blob01e83bcd3067ce5a10d53a93f574a31186eb47e8
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
8 version.
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
17 executable.)
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
22 for more details.
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
27 02110-1301, USA. */
29 #undef IN_LIBGCC2
30 #include "bid-dpd.h"
32 /* get full 64x64bit product */
33 #define __mul_64x64_to_128(P, CX, CY) \
34 { \
35 UINT64 CXH, CXL, CYH,CYL,PL,PH,PM,PM2; \
36 CXH = (CX) >> 32; \
37 CXL = (UINT32)(CX); \
38 CYH = (CY) >> 32; \
39 CYL = (UINT32)(CY); \
41 PM = CXH*CYL; \
42 PH = CXH*CYH; \
43 PL = CXL*CYL; \
44 PM2 = CXL*CYH; \
45 PH += (PM>>32); \
46 PM = (UINT64)((UINT32)PM)+PM2+(PL>>32); \
48 (P).w[1] = PH + (PM>>32); \
49 (P).w[0] = (PM<<32)+(UINT32)PL; \
52 /* add 64-bit value to 128-bit */
53 #define __add_128_64(R128, A128, B64) \
54 { \
55 UINT64 R64H; \
56 R64H = (A128).w[1]; \
57 (R128).w[0] = (B64) + (A128).w[0]; \
58 if((R128).w[0] < (B64)) R64H ++; \
59 (R128).w[1] = R64H; \
62 /* add 128-bit value to 128-bit (assume no carry-out) */
63 #define __add_128_128(R128, A128, B128) \
64 { \
65 UINT128 Q128; \
66 Q128.w[1] = (A128).w[1]+(B128).w[1]; \
67 Q128.w[0] = (B128).w[0] + (A128).w[0]; \
68 if(Q128.w[0] < (B128).w[0]) Q128.w[1] ++; \
69 (R128).w[1] = Q128.w[1]; \
70 (R128).w[0] = Q128.w[0]; \
73 #define __mul_128x128_high(Q, A, B) \
74 { \
75 UINT128 ALBL, ALBH, AHBL, AHBH, QM, QM2; \
77 __mul_64x64_to_128(ALBH, (A).w[0], (B).w[1]); \
78 __mul_64x64_to_128(AHBL, (B).w[0], (A).w[1]); \
79 __mul_64x64_to_128(ALBL, (A).w[0], (B).w[0]); \
80 __mul_64x64_to_128(AHBH, (A).w[1],(B).w[1]); \
82 __add_128_128(QM, ALBH, AHBL); \
83 __add_128_64(QM2, QM, ALBL.w[1]); \
84 __add_128_64((Q), AHBH, QM2.w[1]); \
87 #include "bid2dpd_dpd2bid.h"
89 static const unsigned int dm103[] =
90 { 0, 1000, 2000, 3000, 4000, 5000, 6000, 7000, 8000, 9000, 10000, 11000 };
92 void _bid_to_dpd32 (_Decimal32 *, _Decimal32 *);
94 void
95 _bid_to_dpd32 (_Decimal32 *pres, _Decimal32 *px) {
96 unsigned int sign, coefficient_x, exp, dcoeff;
97 unsigned int b2, b1, b0, b01, res;
98 _Decimal32 x = *px;
100 sign = (x & 0x80000000);
101 if ((x & 0x60000000ul) == 0x60000000ul) {
102 /* special encodings */
103 if ((x & 0x78000000ul) == 0x78000000ul) {
104 *pres = x; /* NaN or Infinity */
105 return;
107 /* coefficient */
108 coefficient_x = (x & 0x001ffffful) | 0x00800000ul;
109 if (coefficient_x >= 10000000) coefficient_x = 0;
110 /* get exponent */
111 exp = (x >> 21) & 0xff;
112 } else {
113 exp = (x >> 23) & 0xff;
114 coefficient_x = (x & 0x007ffffful);
116 b01 = coefficient_x / 1000;
117 b2 = coefficient_x - 1000 * b01;
118 b0 = b01 / 1000;
119 b1 = b01 - 1000 * b0;
120 dcoeff = b2d[b2] | b2d2[b1];
121 if (b0 >= 8) { /* is b0 8 or 9? */
122 res = sign | ((0x600 | ((exp >> 6) << 7) |
123 ((b0 & 1) << 6) | (exp & 0x3f)) << 20) | dcoeff;
124 } else { /* else b0 is 0..7 */
125 res = sign | ((((exp >> 6) << 9) | (b0 << 6) |
126 (exp & 0x3f)) << 20) | dcoeff;
128 *pres = res;
131 void _dpd_to_bid32 (_Decimal32 *, _Decimal32 *);
133 void
134 _dpd_to_bid32 (_Decimal32 *pres, _Decimal32 *px) {
135 unsigned int r;
136 unsigned int sign, exp, bcoeff;
137 UINT64 trailing;
138 unsigned int d0, d1, d2;
139 _Decimal32 x = *px;
141 sign = (x & 0x80000000);
142 trailing = (x & 0x000fffff);
143 if ((x & 0x78000000) == 0x78000000) {
144 *pres = x;
145 return;
146 } else { /* normal number */
147 if ((x & 0x60000000) == 0x60000000) { /* G0..G1 = 11 -> d0 = 8 + G4 */
148 d0 = d2b3[((x >> 26) & 1) | 8]; /* d0 = (comb & 0x0100 ? 9 : 8); */
149 exp = (x >> 27) & 3; /* exp leading bits are G2..G3 */
150 } else {
151 d0 = d2b3[(x >> 26) & 0x7];
152 exp = (x >> 29) & 3; /* exp loading bits are G0..G1 */
154 d1 = d2b2[(trailing >> 10) & 0x3ff];
155 d2 = d2b[(trailing) & 0x3ff];
156 bcoeff = d2 + d1 + d0;
157 exp = (exp << 6) + ((x >> 20) & 0x3f);
158 if (bcoeff < (1 << 23)) {
159 r = exp;
160 r <<= 23;
161 r |= (bcoeff | sign);
162 } else {
163 r = exp;
164 r <<= 21;
165 r |= (sign | 0x60000000ul);
166 /* add coeff, without leading bits */
167 r |= (((unsigned int) bcoeff) & 0x1fffff);
170 *pres = r;
173 void _bid_to_dpd64 (_Decimal64 *, _Decimal64 *);
175 void
176 _bid_to_dpd64 (_Decimal64 *pres, _Decimal64 *px) {
177 UINT64 res;
178 UINT64 sign, comb, exp, B34, B01;
179 UINT64 d103, D61;
180 UINT64 b0, b2, b3, b5;
181 unsigned int b1, b4;
182 UINT64 bcoeff;
183 UINT64 dcoeff;
184 unsigned int yhi, ylo;
185 _Decimal64 x = *px;
187 sign = (x & 0x8000000000000000ull);
188 comb = (x & 0x7ffc000000000000ull) >> 51;
189 if ((comb & 0xf00) == 0xf00) {
190 *pres = x;
191 return;
192 } else { /* Normal number */
193 if ((comb & 0xc00) == 0xc00) { /* G0..G1 = 11 -> exp is G2..G11 */
194 exp = (comb) & 0x3ff;
195 bcoeff = (x & 0x0007ffffffffffffull) | 0x0020000000000000ull;
196 } else {
197 exp = (comb >> 2) & 0x3ff;
198 bcoeff = (x & 0x001fffffffffffffull);
200 D61 = 2305843009ull; /* Floor(2^61 / 10^9) */
201 /* Multiply the binary coefficient by ceil(2^64 / 1000), and take the upper
202 64-bits in order to compute a division by 1000. */
203 yhi = (D61 * (UINT64)(bcoeff >> (UINT64)27)) >> (UINT64)34;
204 ylo = bcoeff - 1000000000ull * yhi;
205 if (ylo >= 1000000000) {
206 ylo = ylo - 1000000000;
207 yhi = yhi + 1;
209 d103 = 0x4189374c;
210 B34 = ((UINT64) ylo * d103) >> (32 + 8);
211 B01 = ((UINT64) yhi * d103) >> (32 + 8);
212 b5 = ylo - B34 * 1000;
213 b2 = yhi - B01 * 1000;
214 b3 = ((UINT64) B34 * d103) >> (32 + 8);
215 b0 = ((UINT64) B01 * d103) >> (32 + 8);
216 b4 = (unsigned int) B34 - (unsigned int) b3 *1000;
217 b1 = (unsigned int) B01 - (unsigned int) dm103[b0];
218 dcoeff = b2d[b5] | b2d2[b4] | b2d3[b3] | b2d4[b2] | b2d5[b1];
219 if (b0 >= 8) /* is b0 8 or 9? */
220 res = sign | ((0x1800 | ((exp >> 8) << 9) | ((b0 & 1) << 8) |
221 (exp & 0xff)) << 50) | dcoeff;
222 else /* else b0 is 0..7 */
223 res = sign | ((((exp >> 8) << 11) | (b0 << 8) |
224 (exp & 0xff)) << 50) | dcoeff;
226 *pres = res;
229 void _dpd_to_bid64 (_Decimal64 *, _Decimal64 *);
231 void
232 _dpd_to_bid64 (_Decimal64 *pres, _Decimal64 *px) {
233 UINT64 res;
234 UINT64 sign, comb, exp;
235 UINT64 trailing;
236 UINT64 d0, d1, d2;
237 unsigned int d3, d4, d5;
238 UINT64 bcoeff, mask;
239 _Decimal64 x = *px;
241 sign = (x & 0x8000000000000000ull);
242 comb = (x & 0x7ffc000000000000ull) >> 50;
243 trailing = (x & 0x0003ffffffffffffull);
244 if ((comb & 0x1e00) == 0x1e00) {
245 if ((comb & 0x1f00) == 0x1f00) { /* G0..G4 = 11111 -> NaN */
246 if (comb & 0x0100) { /* G5 = 1 -> sNaN */
247 *pres = x;
248 } else { /* G5 = 0 -> qNaN */
249 *pres = x;
251 } else { /*if ((comb & 0x1e00) == 0x1e00); G0..G4 = 11110 -> INF */
252 *pres = x;
254 return;
255 } else { /* normal number */
256 if ((comb & 0x1800) == 0x1800) { /* G0..G1 = 11 -> d0 = 8 + G4 */
257 d0 = d2b6[((comb >> 8) & 1) | 8]; /* d0 = (comb & 0x0100 ? 9 : 8); */
258 exp = (comb & 0x600) >> 1; /* exp = (comb & 0x0400 ? 1 : 0) * 0x200 +
259 (comb & 0x0200 ? 1 : 0) * 0x100; exp leading bits are G2..G3 */
260 } else {
261 d0 = d2b6[(comb >> 8) & 0x7];
262 exp = (comb & 0x1800) >> 3; /* exp = (comb & 0x1000 ? 1 : 0) * 0x200 +
263 (comb & 0x0800 ? 1 : 0) * 0x100; exp loading bits are G0..G1 */
265 d1 = d2b5[(trailing >> 40) & 0x3ff];
266 d2 = d2b4[(trailing >> 30) & 0x3ff];
267 d3 = d2b3[(trailing >> 20) & 0x3ff];
268 d4 = d2b2[(trailing >> 10) & 0x3ff];
269 d5 = d2b[(trailing) & 0x3ff];
270 bcoeff = (d5 + d4 + d3) + d2 + d1 + d0;
271 exp += (comb & 0xff);
272 mask = 1;
273 mask <<= 53;
274 if (bcoeff < mask) { /* check whether coefficient fits in 10*5+3 bits */
275 res = exp;
276 res <<= 53;
277 res |= (bcoeff | sign);
278 *pres = res;
279 return;
281 /* special format */
282 res = (exp << 51) | (sign | 0x6000000000000000ull);
283 /* add coeff, without leading bits */
284 mask = (mask >> 2) - 1;
285 bcoeff &= mask;
286 res |= bcoeff;
288 *pres = res;
291 void _bid_to_dpd128 (_Decimal128 *, _Decimal128 *);
293 void
294 _bid_to_dpd128 (_Decimal128 *pres, _Decimal128 *px) {
295 UINT128 res;
296 UINT128 sign;
297 unsigned int comb;
298 UINT128 bcoeff;
299 UINT128 dcoeff;
300 UINT128 BH, d1018, BT2, BT1;
301 UINT64 exp, BL, d109;
302 UINT64 d106, d103;
303 UINT64 k1, k2, k4, k5, k7, k8, k10, k11;
304 unsigned int BHH32, BLL32, BHL32, BLH32, k0, k3, k6, k9, amount;
305 _Decimal128 x = *px;
307 sign.w[1] = (x.w[1] & 0x8000000000000000ull);
308 sign.w[0] = 0;
309 comb = (x.w[1] /*& 0x7fffc00000000000ull */ ) >> 46;
310 exp = 0;
311 if ((comb & 0x1e000) == 0x1e000) {
312 if ((comb & 0x1f000) == 0x1f000) { /* G0..G4 = 11111 -> NaN */
313 if (comb & 0x01000) { /* G5 = 1 -> sNaN */
314 res = x;
315 } else { /* G5 = 0 -> qNaN */
316 res = x;
318 } else { /* G0..G4 = 11110 -> INF */
319 res = x;
321 } else { /* normal number */
322 exp = ((x.w[1] & 0x7fff000000000000ull) >> 49) & 0x3fff;
323 bcoeff.w[1] = (x.w[1] & 0x0001ffffffffffffull);
324 bcoeff.w[0] = x.w[0];
325 d1018 = reciprocals10_128[18];
326 __mul_128x128_high (BH, bcoeff, d1018);
327 amount = recip_scale[18];
328 BH.w[0] = (BH.w[0] >> amount) | (BH.w[1] << (64 - amount));
329 BL = bcoeff.w[0] - BH.w[0] * 1000000000000000000ull;
330 d109 = reciprocals10_64[9];
331 __mul_64x64_to_128 (BT1, BH.w[0], d109);
332 BHH32 = (unsigned int) (BT1.w[1] >> short_recip_scale[9]);
333 BHL32 = (unsigned int) BH.w[0] - BHH32 * 1000000000;
334 __mul_64x64_to_128 (BT2, BL, d109);
335 BLH32 = (unsigned int) (BT2.w[1] >> short_recip_scale[9]);
336 BLL32 = (unsigned int) BL - BLH32 * 1000000000;
337 d106 = 0x431BDE83;
338 d103 = 0x4189374c;
339 k0 = ((UINT64) BHH32 * d106) >> (32 + 18);
340 BHH32 -= (unsigned int) k0 *1000000;
341 k1 = ((UINT64) BHH32 * d103) >> (32 + 8);
342 k2 = BHH32 - (unsigned int) k1 *1000;
343 k3 = ((UINT64) BHL32 * d106) >> (32 + 18);
344 BHL32 -= (unsigned int) k3 *1000000;
345 k4 = ((UINT64) BHL32 * d103) >> (32 + 8);
346 k5 = BHL32 - (unsigned int) k4 *1000;
347 k6 = ((UINT64) BLH32 * d106) >> (32 + 18);
348 BLH32 -= (unsigned int) k6 *1000000;
349 k7 = ((UINT64) BLH32 * d103) >> (32 + 8);
350 k8 = BLH32 - (unsigned int) k7 *1000;
351 k9 = ((UINT64) BLL32 * d106) >> (32 + 18);
352 BLL32 -= (unsigned int) k9 *1000000;
353 k10 = ((UINT64) BLL32 * d103) >> (32 + 8);
354 k11 = BLL32 - (unsigned int) k10 *1000;
355 dcoeff.w[1] = (b2d[k5] >> 4) | (b2d[k4] << 6) | (b2d[k3] << 16) |
356 (b2d[k2] << 26) | (b2d[k1] << 36);
357 dcoeff.w[0] = b2d[k11] | (b2d[k10] << 10) | (b2d[k9] << 20) |
358 (b2d[k8] << 30) | (b2d[k7] << 40) | (b2d[k6] << 50) | (b2d[k5] << 60);
359 res.w[0] = dcoeff.w[0];
360 if (k0 >= 8) {
361 res.w[1] = sign.w[1] | ((0x18000 | ((exp >> 12) << 13) |
362 ((k0 & 1) << 12) | (exp & 0xfff)) << 46) | dcoeff.w[1];
363 } else {
364 res.w[1] = sign.w[1] | ((((exp >> 12) << 15) | (k0 << 12) |
365 (exp & 0xfff)) << 46) | dcoeff.w[1];
368 *pres = res;
371 void _dpd_to_bid128 (_Decimal128 *, _Decimal128 *);
373 void
374 _dpd_to_bid128 (_Decimal128 *pres, _Decimal128 *px) {
375 UINT128 res;
376 UINT128 sign;
377 UINT64 exp, comb;
378 UINT128 trailing;
379 UINT64 d0, d1, d2, d3, d4, d5, d6, d7, d8, d9, d10, d11;
380 UINT128 bcoeff;
381 UINT64 tl, th;
382 _Decimal128 x = *px;
384 sign.w[1] = (x.w[1] & 0x8000000000000000ull);
385 sign.w[0] = 0;
386 comb = (x.w[1] & 0x7fffc00000000000ull) >> 46;
387 trailing.w[1] = x.w[1];
388 trailing.w[0] = x.w[0];
389 if ((comb & 0x1e000) == 0x1e000) {
390 if ((comb & 0x1f000) == 0x1f000) { /* G0..G4 = 11111 -> NaN */
391 if (comb & 0x01000) { /* G5 = 1 -> sNaN */
392 *pres = x;
393 } else { /* G5 = 0 -> qNaN */
394 *pres = x;
396 } else { /* G0..G4 = 11110 -> INF */
397 *pres = x;
399 return;
400 } else { /* Normal number */
401 if ((comb & 0x18000) == 0x18000) { /* G0..G1 = 11 -> d0 = 8 + G4 */
402 d0 = d2b6[8 + ((comb & 0x01000) >> 12)];
403 exp = (comb & 0x06000) >> 1; /* exp leading bits are G2..G3 */
404 } else {
405 d0 = d2b6[((comb & 0x07000) >> 12)];
406 exp = (comb & 0x18000) >> 3; /* exp loading bits are G0..G1 */
408 d11 = d2b[(trailing.w[0]) & 0x3ff];
409 d10 = d2b2[(trailing.w[0] >> 10) & 0x3ff];
410 d9 = d2b3[(trailing.w[0] >> 20) & 0x3ff];
411 d8 = d2b4[(trailing.w[0] >> 30) & 0x3ff];
412 d7 = d2b5[(trailing.w[0] >> 40) & 0x3ff];
413 d6 = d2b6[(trailing.w[0] >> 50) & 0x3ff];
414 d5 = d2b[(trailing.w[0] >> 60) | ((trailing.w[1] & 0x3f) << 4)];
415 d4 = d2b2[(trailing.w[1] >> 6) & 0x3ff];
416 d3 = d2b3[(trailing.w[1] >> 16) & 0x3ff];
417 d2 = d2b4[(trailing.w[1] >> 26) & 0x3ff];
418 d1 = d2b5[(trailing.w[1] >> 36) & 0x3ff];
419 tl = d11 + d10 + d9 + d8 + d7 + d6;
420 th = d5 + d4 + d3 + d2 + d1 + d0;
421 __mul_64x64_to_128 (bcoeff, th, 1000000000000000000ull);
422 __add_128_64 (bcoeff, bcoeff, tl);
423 exp += (comb & 0xfff);
424 res.w[0] = bcoeff.w[0];
425 res.w[1] = (exp << 49) | sign.w[1] | bcoeff.w[1];
427 *pres = res;