libdecnumber/bid/bid2dpd_dpd2bid.c

   1 /* Copyright (C) 2007-2017 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   }
 142   /* normal number */
 143   if ((x & 0x60000000) == 0x60000000) { /* G0..G1 = 11 -> d0 = 8 + G4 */
 144     d0 = d2b3[((x >> 26) & 1) | 8]; /* d0 = (comb & 0x0100 ? 9 : 8); */
 145     exp = (x >> 27) & 3; /* exp leading bits are G2..G3 */
 146   } else {
 147     d0 = d2b3[(x >> 26) & 0x7];
 148     exp = (x >> 29) & 3; /* exp loading bits are G0..G1 */
 149   }
 150   d1 = d2b2[(trailing >> 10) & 0x3ff];
 151   d2 = d2b[(trailing) & 0x3ff];
 152   bcoeff = d2 + d1 + d0;
 153   exp = (exp << 6) + ((x >> 20) & 0x3f);
 154   if (bcoeff < (1 << 23)) {
 155     r = exp;
 156     r <<= 23;
 157     r |= (bcoeff | sign);
 158   } else {
 159     r = exp;
 160     r <<= 21;
 161     r |= (sign | 0x60000000ul);
 162     /* add coeff, without leading bits */
 163     r |= (((unsigned int) bcoeff) & 0x1fffff);
 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   }
 188   /* Normal number */
 189   if ((comb & 0xc00) == 0xc00) { /* G0..G1 = 11 -> exp is G2..G11 */
 190     exp = (comb) & 0x3ff;
 191     bcoeff = (x & 0x0007ffffffffffffull) | 0x0020000000000000ull;
 192   } else {
 193     exp = (comb >> 2) & 0x3ff;
 194     bcoeff = (x & 0x001fffffffffffffull);
 195   }
 196   D61 = 2305843009ull; /* Floor(2^61 / 10^9) */
 197   /* Multiply the binary coefficient by ceil(2^64 / 1000), and take the upper
 198      64-bits in order to compute a division by 1000. */
 199   yhi = (D61 * (UINT64)(bcoeff >> (UINT64)27)) >> (UINT64)34;
 200   ylo = bcoeff - 1000000000ull * yhi;
 201   if (ylo >= 1000000000) {
 202     ylo = ylo - 1000000000;
 203     yhi = yhi + 1;
 204   }
 205   d103 = 0x4189374c;
 206   B34 = ((UINT64) ylo * d103) >> (32 + 8);
 207   B01 = ((UINT64) yhi * d103) >> (32 + 8);
 208   b5 = ylo - B34 * 1000;
 209   b2 = yhi - B01 * 1000;
 210   b3 = ((UINT64) B34 * d103) >> (32 + 8);
 211   b0 = ((UINT64) B01 * d103) >> (32 + 8);
 212   b4 = (unsigned int) B34 - (unsigned int) b3 *1000;
 213   b1 = (unsigned int) B01 - (unsigned int) dm103[b0];
 214   dcoeff = b2d[b5] | b2d2[b4] | b2d3[b3] | b2d4[b2] | b2d5[b1];
 215   if (b0 >= 8) /* is b0 8 or 9? */
 216     res = sign | ((0x1800 | ((exp >> 8) << 9) | ((b0 & 1) << 8) |
 217                    (exp & 0xff)) << 50) | dcoeff;
 218   else /* else b0 is 0..7 */
 219     res = sign | ((((exp >> 8) << 11) | (b0 << 8) |
 220                      (exp & 0xff)) << 50) | dcoeff;
 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     *pres = x;
 241     return;
 242   }
 243   /* normal number */
 244   if ((comb & 0x1800) == 0x1800) { /* G0..G1 = 11 -> d0 = 8 + G4 */
 245     d0 = d2b6[((comb >> 8) & 1) | 8]; /* d0 = (comb & 0x0100 ? 9 : 8); */
 246     exp = (comb & 0x600) >> 1; /* exp = (comb & 0x0400 ? 1 : 0) * 0x200 +
 247         (comb & 0x0200 ? 1 : 0) * 0x100; exp leading bits are G2..G3 */
 248   } else {
 249     d0 = d2b6[(comb >> 8) & 0x7];
 250     exp = (comb & 0x1800) >> 3; /* exp = (comb & 0x1000 ? 1 : 0) * 0x200 +
 251         (comb & 0x0800 ? 1 : 0) * 0x100; exp loading bits are G0..G1 */
 252   }
 253   d1 = d2b5[(trailing >> 40) & 0x3ff];
 254   d2 = d2b4[(trailing >> 30) & 0x3ff];
 255   d3 = d2b3[(trailing >> 20) & 0x3ff];
 256   d4 = d2b2[(trailing >> 10) & 0x3ff];
 257   d5 = d2b[(trailing) & 0x3ff];
 258   bcoeff = (d5 + d4 + d3) + d2 + d1 + d0;
 259   exp += (comb & 0xff);
 260   mask = 1;
 261   mask <<= 53;
 262   if (bcoeff < mask) { /* check whether coefficient fits in 10*5+3 bits */
 263     res = exp;
 264     res <<= 53;
 265     res |= (bcoeff | sign);
 266     *pres = res;
 267     return;
 268   }
 269   /* special format */
 270   res = (exp << 51) | (sign | 0x6000000000000000ull);
 271   /* add coeff, without leading bits */
 272   mask = (mask >> 2) - 1;
 273   bcoeff &= mask;
 274   res |= bcoeff;
 275   *pres = res;
 276 }
 277
 278 void _bid_to_dpd128 (_Decimal128 *, _Decimal128 *);
 279
 280 void
 281 _bid_to_dpd128 (_Decimal128 *pres, _Decimal128 *px) {
 282   UINT128 res;
 283   UINT128 sign;
 284   unsigned int comb;
 285   UINT128 bcoeff;
 286   UINT128 dcoeff;
 287   UINT128 BH, d1018, BT2, BT1;
 288   UINT64 exp, BL, d109;
 289   UINT64 d106, d103;
 290   UINT64 k1, k2, k4, k5, k7, k8, k10, k11;
 291   unsigned int BHH32, BLL32, BHL32, BLH32, k0, k3, k6, k9, amount;
 292   _Decimal128 x = *px;
 293
 294   sign.w[1] = (x.w[1] & 0x8000000000000000ull);
 295   sign.w[0] = 0;
 296   comb = (x.w[1] /*& 0x7fffc00000000000ull */ ) >> 46;
 297   exp = 0;
 298   if ((comb & 0x1e000) == 0x1e000) {
 299     res = x;
 300   } else { /* normal number */
 301     exp = ((x.w[1] & 0x7fff000000000000ull) >> 49) & 0x3fff;
 302     bcoeff.w[1] = (x.w[1] & 0x0001ffffffffffffull);
 303     bcoeff.w[0] = x.w[0];
 304     d1018 = reciprocals10_128[18];
 305     __mul_128x128_high (BH, bcoeff, d1018);
 306     amount = recip_scale[18];
 307     BH.w[0] = (BH.w[0] >> amount) | (BH.w[1] << (64 - amount));
 308     BL = bcoeff.w[0] - BH.w[0] * 1000000000000000000ull;
 309     d109 = reciprocals10_64[9];
 310     __mul_64x64_to_128 (BT1, BH.w[0], d109);
 311     BHH32 = (unsigned int) (BT1.w[1] >> short_recip_scale[9]);
 312     BHL32 = (unsigned int) BH.w[0] - BHH32 * 1000000000;
 313     __mul_64x64_to_128 (BT2, BL, d109);
 314     BLH32 = (unsigned int) (BT2.w[1] >> short_recip_scale[9]);
 315     BLL32 = (unsigned int) BL - BLH32 * 1000000000;
 316     d106 = 0x431BDE83;
 317     d103 = 0x4189374c;
 318     k0 = ((UINT64) BHH32 * d106) >> (32 + 18);
 319     BHH32 -= (unsigned int) k0 *1000000;
 320     k1 = ((UINT64) BHH32 * d103) >> (32 + 8);
 321     k2 = BHH32 - (unsigned int) k1 *1000;
 322     k3 = ((UINT64) BHL32 * d106) >> (32 + 18);
 323     BHL32 -= (unsigned int) k3 *1000000;
 324     k4 = ((UINT64) BHL32 * d103) >> (32 + 8);
 325     k5 = BHL32 - (unsigned int) k4 *1000;
 326     k6 = ((UINT64) BLH32 * d106) >> (32 + 18);
 327     BLH32 -= (unsigned int) k6 *1000000;
 328     k7 = ((UINT64) BLH32 * d103) >> (32 + 8);
 329     k8 = BLH32 - (unsigned int) k7 *1000;
 330     k9 = ((UINT64) BLL32 * d106) >> (32 + 18);
 331     BLL32 -= (unsigned int) k9 *1000000;
 332     k10 = ((UINT64) BLL32 * d103) >> (32 + 8);
 333     k11 = BLL32 - (unsigned int) k10 *1000;
 334     dcoeff.w[1] = (b2d[k5] >> 4) | (b2d[k4] << 6) | (b2d[k3] << 16) |
 335       (b2d[k2] << 26) | (b2d[k1] << 36);
 336     dcoeff.w[0] = b2d[k11] | (b2d[k10] << 10) | (b2d[k9] << 20) |
 337       (b2d[k8] << 30) | (b2d[k7] << 40) | (b2d[k6] << 50) | (b2d[k5] << 60);
 338     res.w[0] = dcoeff.w[0];
 339     if (k0 >= 8) {
 340       res.w[1] = sign.w[1] | ((0x18000 | ((exp >> 12) << 13) |
 341           ((k0 & 1) << 12) | (exp & 0xfff)) << 46) | dcoeff.w[1];
 342     } else {
 343       res.w[1] = sign.w[1] | ((((exp >> 12) << 15) | (k0 << 12) |
 344           (exp & 0xfff)) << 46) | dcoeff.w[1];
 345     }
 346   }
 347   *pres = res;
 348 }
 349
 350 void _dpd_to_bid128 (_Decimal128 *, _Decimal128 *);
 351
 352 void
 353 _dpd_to_bid128 (_Decimal128 *pres, _Decimal128 *px) {
 354   UINT128 res;
 355   UINT128 sign;
 356   UINT64 exp, comb;
 357   UINT128 trailing;
 358   UINT64 d0, d1, d2, d3, d4, d5, d6, d7, d8, d9, d10, d11;
 359   UINT128 bcoeff;
 360   UINT64 tl, th;
 361   _Decimal128 x = *px;
 362
 363   sign.w[1] = (x.w[1] & 0x8000000000000000ull);
 364   sign.w[0] = 0;
 365   comb = (x.w[1] & 0x7fffc00000000000ull) >> 46;
 366   trailing.w[1] = x.w[1];
 367   trailing.w[0] = x.w[0];
 368   if ((comb & 0x1e000) == 0x1e000) {
 369       *pres = x;
 370       return;
 371   }
 372   if ((comb & 0x18000) == 0x18000) { /* G0..G1 = 11 -> d0 = 8 + G4 */
 373     d0 = d2b6[8 + ((comb & 0x01000) >> 12)];
 374     exp = (comb & 0x06000) >> 1;  /* exp leading bits are G2..G3 */
 375   } else {
 376     d0 = d2b6[((comb & 0x07000) >> 12)];
 377     exp = (comb & 0x18000) >> 3;  /* exp loading bits are G0..G1 */
 378   }
 379   d11 = d2b[(trailing.w[0]) & 0x3ff];
 380   d10 = d2b2[(trailing.w[0] >> 10) & 0x3ff];
 381   d9 = d2b3[(trailing.w[0] >> 20) & 0x3ff];
 382   d8 = d2b4[(trailing.w[0] >> 30) & 0x3ff];
 383   d7 = d2b5[(trailing.w[0] >> 40) & 0x3ff];
 384   d6 = d2b6[(trailing.w[0] >> 50) & 0x3ff];
 385   d5 = d2b[(trailing.w[0] >> 60) | ((trailing.w[1] & 0x3f) << 4)];
 386   d4 = d2b2[(trailing.w[1] >> 6) & 0x3ff];
 387   d3 = d2b3[(trailing.w[1] >> 16) & 0x3ff];
 388   d2 = d2b4[(trailing.w[1] >> 26) & 0x3ff];
 389   d1 = d2b5[(trailing.w[1] >> 36) & 0x3ff];
 390   tl = d11 + d10 + d9 + d8 + d7 + d6;
 391   th = d5 + d4 + d3 + d2 + d1 + d0;
 392   __mul_64x64_to_128 (bcoeff, th, 1000000000000000000ull);
 393   __add_128_64 (bcoeff, bcoeff, tl);
 394   exp += (comb & 0xfff);
 395   res.w[0] = bcoeff.w[0];
 396   res.w[1] = (exp << 49) | sign.w[1] | bcoeff.w[1];
 397   *pres = res;
 398 }