libdecnumber/bid/bid2dpd_dpd2bid.c

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