Subversion Repositories Kolibri OS

Rev

Blame | Last modification | View Log | RSS feed

  1. /* Copyright (C) 2007-2015 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  DECNUMDIGITS 34        // work with up to 34 digits
  25.  
  26. #include "bid_internal.h"
  27. #include "bid_b2d.h"
  28.  
  29. UINT32
  30. bid_to_bid32 (UINT32 ba) {
  31.   UINT32 res;
  32.   UINT32 sign, comb, exp;
  33.   UINT32 trailing;
  34.   UINT32 bcoeff;
  35.  
  36.   sign = (ba & 0x80000000ul);
  37.   comb = (ba & 0x7ff00000ul) >> 20;
  38.   trailing = (ba & 0x000ffffful);
  39.  
  40.   if ((comb & 0x780) == 0x780) {
  41.     ba &= 0xfff0000ul;
  42.     return ba;
  43.   } else {
  44.     if ((comb & 0x600) == 0x600) {      // G0..G1 = 11 -> exp is G2..G11
  45.       exp = (comb >> 1) & 0xff;
  46.       bcoeff = ((8 + (comb & 1)) << 20) | trailing;
  47.     } else {
  48.       exp = (comb >> 3) & 0xff;
  49.       bcoeff = ((comb & 7) << 20) | trailing;
  50.     }
  51.     if (bcoeff >= 10000000)
  52.       bcoeff = 0;
  53.     res = very_fast_get_BID32 (sign, exp, bcoeff);
  54.   }
  55.   return res;
  56. }
  57.  
  58. UINT64
  59. bid_to_bid64 (UINT64 ba) {
  60.   UINT64 res;
  61.   UINT64 sign, comb, exp;
  62.   UINT64 trailing;
  63.   UINT64 bcoeff;
  64.  
  65.   sign = (ba & 0x8000000000000000ull);
  66.   comb = (ba & 0x7ffc000000000000ull) >> 50;
  67.   trailing = (ba & 0x0003ffffffffffffull);
  68.  
  69.   if ((comb & 0x1e00) == 0x1e00) {
  70.     ba &= 0xfff000000000000ULL;
  71.     return ba;
  72.   } else {
  73.     if ((comb & 0x1800) == 0x1800) {    // G0..G1 = 11 -> exp is G2..G11
  74.       exp = (comb >> 1) & 0x3ff;
  75.       bcoeff = ((8 + (comb & 1)) << 50) | trailing;
  76.     } else {
  77.       exp = (comb >> 3) & 0x3ff;
  78.       bcoeff = ((comb & 7) << 50) | trailing;
  79.     }
  80.     if (bcoeff >= 10000000000000000ull)
  81.       bcoeff = 0ull;
  82.     res = very_fast_get_BID64 (sign, exp, bcoeff);
  83.   }
  84.   return res;
  85. }
  86.  
  87. #if DECIMAL_CALL_BY_REFERENCE
  88. void
  89. bid_to_dpd32 (UINT32 * pres, UINT32 * pba) {
  90.   UINT32 ba = *pba;
  91. #else
  92. UINT32
  93. bid_to_dpd32 (UINT32 ba) {
  94. #endif
  95.   UINT32 res;
  96.  
  97.   UINT32 sign, comb, exp, trailing;
  98.   UINT32 b0, b1, b2;
  99.   UINT32 bcoeff, dcoeff;
  100.   UINT32 nanb = 0;
  101.  
  102.   sign = (ba & 0x80000000);
  103.   comb = (ba & 0x7ff00000) >> 20;
  104.   trailing = (ba & 0xfffff);
  105.  
  106.   // Detect infinity, and return canonical infinity
  107.   if ((comb & 0x7c0) == 0x780) {
  108.     res = sign | 0x78000000;
  109.     BID_RETURN (res);
  110.     // Detect NaN, and canonicalize trailing
  111.   } else if ((comb & 0x7c0) == 0x7c0) {
  112.     if (trailing > 999999)
  113.       trailing = 0;
  114.     nanb = ba & 0xfe000000;
  115.     exp = 0;
  116.     bcoeff = trailing;
  117.   } else {      // Normal number
  118.     if ((comb & 0x600) == 0x600) {      // G0..G1 = 11 -> exp is G2..G11
  119.       exp = (comb >> 1) & 0xff;
  120.       bcoeff = ((8 + (comb & 1)) << 20) | trailing;
  121.     } else {
  122.       exp = (comb >> 3) & 0xff;
  123.       bcoeff = ((comb & 7) << 20) | trailing;
  124.     }
  125.     // Zero the coefficient if non-canonical (>= 10^7)
  126.     if (bcoeff >= 10000000)
  127.       bcoeff = 0;
  128.   }
  129.  
  130.   b0 = bcoeff / 1000000;
  131.   b1 = (bcoeff / 1000) % 1000;
  132.   b2 = bcoeff % 1000;
  133.   dcoeff = (b2d[b1] << 10) | b2d[b2];
  134.  
  135.   if (b0 >= 8)  // is b0 8 or 9?
  136.     res =
  137.       sign |
  138.       ((0x600 | ((exp >> 6) << 7) | ((b0 & 1) << 6) | (exp & 0x3f)) <<
  139.        20) | dcoeff;
  140.   else  // else b0 is 0..7
  141.     res =
  142.       sign | ((((exp >> 6) << 9) | (b0 << 6) | (exp & 0x3f)) << 20) |
  143.       dcoeff;
  144.  
  145.   res |= nanb;
  146.  
  147.   BID_RETURN (res);
  148. }
  149.  
  150. #if DECIMAL_CALL_BY_REFERENCE
  151. void
  152. bid_to_dpd64 (UINT64 * pres, UINT64 * pba) {
  153.   UINT64 ba = *pba;
  154. #else
  155. UINT64
  156. bid_to_dpd64 (UINT64 ba) {
  157. #endif
  158.   UINT64 res;
  159.  
  160.   UINT64 sign, comb, exp;
  161.   UINT64 trailing;
  162.   UINT32 b0, b1, b2, b3, b4, b5;
  163.   UINT64 bcoeff;
  164.   UINT64 dcoeff;
  165.   UINT32 yhi, ylo;
  166.   UINT64 nanb = 0;
  167.  
  168. //printf("arg bid "FMT_LLX16" \n", ba);
  169.   sign = (ba & 0x8000000000000000ull);
  170.   comb = (ba & 0x7ffc000000000000ull) >> 50;
  171.   trailing = (ba & 0x0003ffffffffffffull);
  172.  
  173.   // Detect infinity, and return canonical infinity
  174.   if ((comb & 0x1f00) == 0x1e00) {
  175.     res = sign | 0x7800000000000000ull;
  176.     BID_RETURN (res);
  177.     // Detect NaN, and canonicalize trailing
  178.   } else if ((comb & 0x1e00) == 0x1e00) {
  179.     if (trailing > 999999999999999ull)
  180.       trailing = 0;
  181.     nanb = ba & 0xfe00000000000000ull;
  182.     exp = 0;
  183.     bcoeff = trailing;
  184.   } else {      // Normal number
  185.     if ((comb & 0x1800) == 0x1800) {    // G0..G1 = 11 -> exp is G2..G11
  186.       exp = (comb >> 1) & 0x3ff;
  187.       bcoeff = ((8 + (comb & 1)) << 50) | trailing;
  188.     } else {
  189.       exp = (comb >> 3) & 0x3ff;
  190.       bcoeff = ((comb & 7) << 50) | trailing;
  191.     }
  192.  
  193.     // Zero the coefficient if it is non-canonical (>= 10^16)
  194.     if (bcoeff >= 10000000000000000ull)
  195.       bcoeff = 0;
  196.   }
  197.  
  198. // Floor(2^61 / 10^9)
  199. #define D61 (2305843009ull)
  200.  
  201. // Multipy the binary coefficient by ceil(2^64 / 1000), and take the upper
  202. // 64-bits in order to compute a division by 1000.
  203.  
  204. #if 1
  205.   yhi =
  206.     ((UINT64) D61 *
  207.      (UINT64) (UINT32) (bcoeff >> (UINT64) 27)) >> (UINT64) 34;
  208.   ylo = bcoeff - 1000000000ull * yhi;
  209.   if (ylo >= 1000000000) {
  210.     ylo = ylo - 1000000000;
  211.     yhi = yhi + 1;
  212.   }
  213. #else
  214.   yhi = bcoeff / 1000000000ull;
  215.   ylo = bcoeff % 1000000000ull;
  216. #endif
  217.  
  218.   // yhi = ABBBCCC ylo = DDDEEEFFF
  219.   b5 = ylo % 1000;      // b5 = FFF
  220.   b3 = ylo / 1000000;   // b3 = DDD
  221.   b4 = (ylo / 1000) - (1000 * b3);      // b4 = EEE
  222.   b2 = yhi % 1000;      // b2 = CCC
  223.   b0 = yhi / 1000000;   // b0 = A
  224.   b1 = (yhi / 1000) - (1000 * b0);      // b1 = BBB
  225.  
  226.   dcoeff = b2d[b5] | b2d2[b4] | b2d3[b3] | b2d4[b2] | b2d5[b1];
  227.  
  228.   if (b0 >= 8)  // is b0 8 or 9?
  229.     res =
  230.       sign |
  231.       ((0x1800 | ((exp >> 8) << 9) | ((b0 & 1) << 8) | (exp & 0xff)) <<
  232.        50) | dcoeff;
  233.   else  // else b0 is 0..7
  234.     res =
  235.       sign | ((((exp >> 8) << 11) | (b0 << 8) | (exp & 0xff)) << 50) |
  236.       dcoeff;
  237.  
  238.   res |= nanb;
  239.  
  240.   BID_RETURN (res);
  241. }
  242.  
  243. #if DECIMAL_CALL_BY_REFERENCE
  244. void
  245. dpd_to_bid32 (UINT32 * pres, UINT32 * pda) {
  246.   UINT32 da = *pda;
  247. #else
  248. UINT32
  249. dpd_to_bid32 (UINT32 da) {
  250. #endif
  251.   UINT32 in = *(UINT32 *) & da;
  252.   UINT32 res;
  253.  
  254.   UINT32 sign, comb, exp;
  255.   UINT32 trailing;
  256.   UINT32 d0 = 0, d1, d2;
  257.   UINT64 bcoeff;
  258.   UINT32 nanb = 0;
  259.  
  260.   sign = (in & 0x80000000);
  261.   comb = (in & 0x7ff00000) >> 20;
  262.   trailing = (in & 0x000fffff);
  263.  
  264.   if ((comb & 0x7e0) == 0x780) {        // G0..G4 = 1111 -> Inf
  265.     res = in & 0xf8000000;
  266.     BID_RETURN (res);
  267.   } else if ((comb & 0x7c0) == 0x7c0) { // G0..G5 = 11111 -> NaN
  268.     nanb = in & 0xfe000000;
  269.     exp = 0;
  270.   } else {      // Normal number
  271.     if ((comb & 0x600) == 0x600) {      // G0..G1 = 11 -> d0 = 8 + G4
  272.       d0 = ((comb >> 6) & 1) | 8;
  273.       exp = ((comb & 0x180) >> 1) | (comb & 0x3f);
  274.     } else {
  275.       d0 = (comb >> 6) & 0x7;
  276.       exp = ((comb & 0x600) >> 3) | (comb & 0x3f);
  277.     }
  278.   }
  279.   d1 = d2b2[(trailing >> 10) & 0x3ff];
  280.   d2 = d2b[(trailing) & 0x3ff];
  281.  
  282.   bcoeff = d2 + d1 + (1000000 * d0);
  283.   if (bcoeff < 0x800000) {
  284.     res = (exp << 23) | bcoeff | sign;
  285.   } else {
  286.     res = (exp << 21) | sign | 0x60000000 | (bcoeff & 0x1fffff);
  287.   }
  288.  
  289.   res |= nanb;
  290.  
  291.   BID_RETURN (res);
  292. }
  293.  
  294. #if DECIMAL_CALL_BY_REFERENCE
  295. void
  296. dpd_to_bid64 (UINT64 * pres, UINT64 * pda) {
  297.   UINT64 da = *pda;
  298. #else
  299. UINT64
  300. dpd_to_bid64 (UINT64 da) {
  301. #endif
  302.   UINT64 in = *(UINT64 *) & da;
  303.   UINT64 res;
  304.  
  305.   UINT64 sign, comb, exp;
  306.   UINT64 trailing;
  307.   // UINT64 d0, d1, d2, d3, d4, d5;
  308.  
  309.   UINT64 d1, d2;
  310.   UINT32 d0, d3, d4, d5;
  311.   UINT64 bcoeff;
  312.   UINT64 nanb = 0;
  313.  
  314. //printf("arg dpd "FMT_LLX16" \n", in);
  315.   sign = (in & 0x8000000000000000ull);
  316.   comb = (in & 0x7ffc000000000000ull) >> 50;
  317.   trailing = (in & 0x0003ffffffffffffull);
  318.  
  319.   if ((comb & 0x1f00) == 0x1e00) {      // G0..G4 = 1111 -> Inf
  320.     res = in & 0xf800000000000000ull;
  321.     BID_RETURN (res);
  322.   } else if ((comb & 0x1f00) == 0x1f00) {       // G0..G5 = 11111 -> NaN
  323.     nanb = in & 0xfe00000000000000ull;
  324.     exp = 0;
  325.     d0 = 0;
  326.   } else {      // Normal number
  327.     if ((comb & 0x1800) == 0x1800) {    // G0..G1 = 11 -> d0 = 8 + G4
  328.       d0 = ((comb >> 8) & 1) | 8;
  329.       // d0 = (comb & 0x0100 ? 9 : 8);
  330.       exp = (comb & 0x600) >> 1;
  331.       // exp = (comb & 0x0400 ? 1 : 0) * 0x200 + (comb & 0x0200 ? 1 : 0) * 0x100; // exp leading bits are G2..G3
  332.     } else {
  333.       d0 = (comb >> 8) & 0x7;
  334.       exp = (comb & 0x1800) >> 3;
  335.       // exp = (comb & 0x1000 ? 1 : 0) * 0x200 + (comb & 0x0800 ? 1 : 0) * 0x100; // exp loading bits are G0..G1
  336.     }
  337.   }
  338.   d1 = d2b5[(trailing >> 40) & 0x3ff];
  339.   d2 = d2b4[(trailing >> 30) & 0x3ff];
  340.   d3 = d2b3[(trailing >> 20) & 0x3ff];
  341.   d4 = d2b2[(trailing >> 10) & 0x3ff];
  342.   d5 = d2b[(trailing) & 0x3ff];
  343.  
  344.   bcoeff = (d5 + d4 + d3) + d2 + d1 + (1000000000000000ull * d0);
  345.   exp += (comb & 0xff);
  346.   res = very_fast_get_BID64 (sign, exp, bcoeff);
  347.  
  348.   res |= nanb;
  349.  
  350.   BID_RETURN (res);
  351. }
  352.  
  353. #if DECIMAL_CALL_BY_REFERENCE
  354. void
  355. bid_to_dpd128 (UINT128 * pres, UINT128 * pba) {
  356.   UINT128 ba = *pba;
  357. #else
  358. UINT128
  359. bid_to_dpd128 (UINT128 ba) {
  360. #endif
  361.   UINT128 res;
  362.  
  363.   UINT128 sign;
  364.   UINT32 comb, exp;
  365.   UINT128 trailing;
  366.   UINT128 d0, d1, d2, d3, d4, d5, d6, d7, d8, d9, d10, d11;
  367.   UINT128 bcoeff;
  368.   UINT128 dcoeff;
  369.   UINT64 nanb = 0;
  370.  
  371.   sign.w[1] = (ba.w[HIGH_128W] & 0x8000000000000000ull);
  372.   sign.w[0] = 0;
  373.   comb = (ba.w[HIGH_128W] & 0x7fffc00000000000ull) >> 46;
  374.   trailing.w[1] = (ba.w[HIGH_128W] & 0x00003fffffffffffull);
  375.   trailing.w[0] = ba.w[LOW_128W];
  376.   exp = 0;
  377.  
  378.   if ((comb & 0x1f000) == 0x1e000) {    // G0..G4 = 1111 -> Inf
  379.     res.w[HIGH_128W] = ba.w[HIGH_128W] & 0xf800000000000000ull;
  380.     res.w[LOW_128W] = 0;
  381.     BID_RETURN (res);
  382.     // Detect NaN, and canonicalize trailing
  383.   } else if ((comb & 0x1f000) == 0x1f000) {
  384.     if ((trailing.w[1] > 0x0000314dc6448d93ULL) ||      // significand is non-canonical
  385.         ((trailing.w[1] == 0x0000314dc6448d93ULL)
  386.          && (trailing.w[0] >= 0x38c15b0a00000000ULL))
  387.         // significand is non-canonical
  388.       ) {
  389.       trailing.w[1] = trailing.w[0] = 0ull;
  390.     }
  391.     bcoeff.w[1] = trailing.w[1];
  392.     bcoeff.w[0] = trailing.w[0];
  393.     nanb = ba.w[HIGH_128W] & 0xfe00000000000000ull;
  394.     exp = 0;
  395.   } else {      // Normal number
  396.     if ((comb & 0x18000) == 0x18000) {  // G0..G1 = 11 -> exp is G2..G11
  397.       exp = (comb >> 1) & 0x3fff;
  398.       bcoeff.w[1] =
  399.         ((UINT64) (8 + (comb & 1)) << (UINT64) 46) | trailing.w[1];
  400.       bcoeff.w[0] = trailing.w[0];
  401.     } else {
  402.       exp = (comb >> 3) & 0x3fff;
  403.       bcoeff.w[1] =
  404.         ((UINT64) (comb & 7) << (UINT64) 46) | trailing.w[1];
  405.       bcoeff.w[0] = trailing.w[0];
  406.     }
  407.     // Zero the coefficient if non-canonical (>= 10^34)
  408.     if (bcoeff.w[1] > 0x1ed09bead87c0ull ||
  409.         (bcoeff.w[1] == 0x1ed09bead87c0ull
  410.          && bcoeff.w[0] >= 0x378D8E6400000000ull)) {
  411.       bcoeff.w[0] = bcoeff.w[1] = 0;
  412.     }
  413.   }
  414.   // Constant 2^128 / 1000 + 1
  415.   {
  416.     UINT128 t;
  417.     UINT64 t2;
  418.     UINT128 d1000;
  419.     UINT128 b11, b10, b9, b8, b7, b6, b5, b4, b3, b2, b1;
  420.     d1000.w[1] = 0x4189374BC6A7EFull;
  421.     d1000.w[0] = 0x9DB22D0E56041894ull;
  422.     __mul_128x128_high (b11, bcoeff, d1000);
  423.     __mul_128x128_high (b10, b11, d1000);
  424.     __mul_128x128_high (b9, b10, d1000);
  425.     __mul_128x128_high (b8, b9, d1000);
  426.     __mul_128x128_high (b7, b8, d1000);
  427.     __mul_128x128_high (b6, b7, d1000);
  428.     __mul_128x128_high (b5, b6, d1000);
  429.     __mul_128x128_high (b4, b5, d1000);
  430.     __mul_128x128_high (b3, b4, d1000);
  431.     __mul_128x128_high (b2, b3, d1000);
  432.     __mul_128x128_high (b1, b2, d1000);
  433.  
  434.  
  435.     __mul_64x128_full (t2, t, 1000ull, b11);
  436.     __sub_128_128 (d11, bcoeff, t);
  437.     __mul_64x128_full (t2, t, 1000ull, b10);
  438.     __sub_128_128 (d10, b11, t);
  439.     __mul_64x128_full (t2, t, 1000ull, b9);
  440.     __sub_128_128 (d9, b10, t);
  441.     __mul_64x128_full (t2, t, 1000ull, b8);
  442.     __sub_128_128 (d8, b9, t);
  443.     __mul_64x128_full (t2, t, 1000ull, b7);
  444.     __sub_128_128 (d7, b8, t);
  445.     __mul_64x128_full (t2, t, 1000ull, b6);
  446.     __sub_128_128 (d6, b7, t);
  447.     __mul_64x128_full (t2, t, 1000ull, b5);
  448.     __sub_128_128 (d5, b6, t);
  449.     __mul_64x128_full (t2, t, 1000ull, b4);
  450.     __sub_128_128 (d4, b5, t);
  451.     __mul_64x128_full (t2, t, 1000ull, b3);
  452.     __sub_128_128 (d3, b4, t);
  453.     __mul_64x128_full (t2, t, 1000ull, b2);
  454.     __sub_128_128 (d2, b3, t);
  455.     __mul_64x128_full (t2, t, 1000ull, b1);
  456.     __sub_128_128 (d1, b2, t);
  457.     d0 = b1;
  458.  
  459.   }
  460.  
  461.   dcoeff.w[0] = b2d[d11.w[0]] | (b2d[d10.w[0]] << 10) |
  462.     (b2d[d9.w[0]] << 20) | (b2d[d8.w[0]] << 30) | (b2d[d7.w[0]] << 40) |
  463.     (b2d[d6.w[0]] << 50) | (b2d[d5.w[0]] << 60);
  464.   dcoeff.w[1] =
  465.     (b2d[d5.w[0]] >> 4) | (b2d[d4.w[0]] << 6) | (b2d[d3.w[0]] << 16) |
  466.     (b2d[d2.w[0]] << 26) | (b2d[d1.w[0]] << 36);
  467.  
  468.   res.w[0] = dcoeff.w[0];
  469.   if (d0.w[0] >= 8) {
  470.     res.w[1] =
  471.       sign.
  472.       w[1] |
  473.       ((0x18000 | ((exp >> 12) << 13) | ((d0.w[0] & 1) << 12) |
  474.         (exp & 0xfff)) << 46) | dcoeff.w[1];
  475.   } else {
  476.     res.w[1] =
  477.       sign.
  478.       w[1] | ((((exp >> 12) << 15) | (d0.w[0] << 12) | (exp & 0xfff))
  479.               << 46) | dcoeff.w[1];
  480.   }
  481.  
  482.   res.w[1] |= nanb;
  483.  
  484.   BID_SWAP128 (res);
  485.   BID_RETURN (res);
  486. }
  487.  
  488. #if DECIMAL_CALL_BY_REFERENCE
  489. void
  490. dpd_to_bid128 (UINT128 * pres, UINT128 * pda) {
  491.   UINT128 da = *pda;
  492. #else
  493. UINT128
  494. dpd_to_bid128 (UINT128 da) {
  495. #endif
  496.   UINT128 in = *(UINT128 *) & da;
  497.   UINT128 res;
  498.  
  499.   UINT128 sign;
  500.   UINT64 exp, comb;
  501.   UINT128 trailing;
  502.   UINT64 d0, d1, d2, d3, d4, d5, d6, d7, d8, d9, d10, d11;
  503.   UINT128 bcoeff;
  504.   UINT64 tl, th;
  505.   UINT64 nanb = 0;
  506.  
  507.   sign.w[1] = (in.w[HIGH_128W] & 0x8000000000000000ull);
  508.   sign.w[0] = 0;
  509.   comb = (in.w[HIGH_128W] & 0x7fffc00000000000ull) >> 46;
  510.   trailing.w[1] = (in.w[HIGH_128W] & 0x00003fffffffffffull);
  511.   trailing.w[0] = in.w[LOW_128W];
  512.   exp = 0;
  513.  
  514.   if ((comb & 0x1f000) == 0x1e000) {    // G0..G4 = 1111 -> Inf
  515.     res.w[HIGH_128W] = in.w[HIGH_128W] & 0xf800000000000000ull;
  516.     res.w[LOW_128W] = 0ull;
  517.     BID_RETURN (res);
  518.   } else if ((comb & 0x1f000) == 0x1f000) {     // G0..G4 = 11111 -> NaN
  519.     nanb = in.w[HIGH_128W] & 0xfe00000000000000ull;
  520.     exp = 0;
  521.     d0 = 0;
  522.   } else {      // Normal number
  523.     if ((comb & 0x18000) == 0x18000) {  // G0..G1 = 11 -> d0 = 8 + G4
  524.       d0 = 8 + (comb & 0x01000 ? 1 : 0);
  525.       exp =
  526.         (comb & 0x04000 ? 1 : 0) * 0x2000 +
  527.         (comb & 0x02000 ? 1 : 0) * 0x1000;
  528.       // exp leading bits are G2..G3
  529.     } else {
  530.       d0 =
  531.         4 * (comb & 0x04000 ? 1 : 0) + 2 * (comb & 0x2000 ? 1 : 0) +
  532.         (comb & 0x1000 ? 1 : 0);
  533.       exp =
  534.         (comb & 0x10000 ? 1 : 0) * 0x2000 +
  535.         (comb & 0x08000 ? 1 : 0) * 0x1000;
  536.       // exp loading bits are G0..G1
  537.     }
  538.   }
  539.  
  540.   d11 = d2b[(trailing.w[0]) & 0x3ff];
  541.   d10 = d2b[(trailing.w[0] >> 10) & 0x3ff];
  542.   d9 = d2b[(trailing.w[0] >> 20) & 0x3ff];
  543.   d8 = d2b[(trailing.w[0] >> 30) & 0x3ff];
  544.   d7 = d2b[(trailing.w[0] >> 40) & 0x3ff];
  545.   d6 = d2b[(trailing.w[0] >> 50) & 0x3ff];
  546.   d5 = d2b[(trailing.w[0] >> 60) | ((trailing.w[1] & 0x3f) << 4)];
  547.   d4 = d2b[(trailing.w[1] >> 6) & 0x3ff];
  548.   d3 = d2b[(trailing.w[1] >> 16) & 0x3ff];
  549.   d2 = d2b[(trailing.w[1] >> 26) & 0x3ff];
  550.   d1 = d2b[(trailing.w[1] >> 36) & 0x3ff];
  551.  
  552.   tl =
  553.     d11 + (d10 * 1000ull) + (d9 * 1000000ull) + (d8 * 1000000000ull) +
  554.     (d7 * 1000000000000ull) + (d6 * 1000000000000000ull);
  555.   th =
  556.     d5 + (d4 * 1000ull) + (d3 * 1000000ull) + (d2 * 1000000000ull) +
  557.     (d1 * 1000000000000ull) + (d0 * 1000000000000000ull);
  558.   __mul_64x64_to_128 (bcoeff, th, 1000000000000000000ull);
  559.   __add_128_64 (bcoeff, bcoeff, tl);
  560.  
  561.   if (!nanb)
  562.     exp += (comb & 0xfff);
  563.  
  564.   res.w[0] = bcoeff.w[0];
  565.   res.w[1] = (exp << 49) | sign.w[1] | bcoeff.w[1];
  566.  
  567.   res.w[1] |= nanb;
  568.  
  569.   BID_SWAP128 (res);
  570.   BID_RETURN (res);
  571. }
  572.  
  573. UINT128
  574. bid_to_bid128 (UINT128 bq) {
  575.   UINT128 res;
  576.   UINT64 sign, comb, exp;
  577.   UINT64 trailing;
  578.   UINT64 bcoeff;
  579.  
  580.   UINT128 rq;
  581.   UINT128 qcoeff;
  582.   UINT64 ba, bb;
  583.  
  584.   ba = *((UINT64 *) & bq + HIGH_128W);
  585.   bb = *((UINT64 *) & bq + LOW_128W);
  586.  
  587.   sign = (ba & 0x8000000000000000ull);
  588.   comb = (ba & 0x7fffc00000000000ull) >> 46;
  589.   trailing = (ba & 0x00003fffffffffffull);
  590.  
  591.   if ((comb & 0x18000) == 0x18000) {    // G0..G1 = 11 -> exp is G2..G11
  592.     exp = (comb >> 1) & 0x3fff;
  593.     bcoeff = ((8 + (comb & 1)) << 46) | trailing;
  594.   } else {
  595.     exp = (comb >> 3) & 0x3fff;
  596.     bcoeff = ((comb & 7) << 46) | trailing;
  597.   }
  598.  
  599.   if ((comb & 0x1f000) == 0x1f000) {    //NaN
  600.     ba &= 0xfe003fffffffffffULL;        // make exponent 0
  601.     bcoeff &= 0x00003fffffffffffull;    // NaN payloat is only T.  
  602.     if ((bcoeff > 0x0000314dc6448d93ULL) ||     // significand is non-canonical
  603.         ((bcoeff == 0x0000314dc6448d93ULL)
  604.          && (bb >= 0x38c15b0a00000000ULL))
  605.         // significand is non-canonical
  606.       ) {
  607.       bcoeff = 0ull;
  608.       ba &= ~0x00003fffffffffffull;
  609.       bb = 0ull;
  610.     }
  611.     *((UINT64 *) & rq + HIGH_128W) = ba;
  612.     *((UINT64 *) & rq + LOW_128W) = bb;
  613.     return rq;
  614.   } else if ((comb & 0x1e000) == 0x1e000) {     //Inf
  615.     ba &= 0xf800000000000000ULL;        // make exponent and significand 0
  616.     bb = 0;
  617.     *((UINT64 *) & rq + HIGH_128W) = ba;
  618.     *((UINT64 *) & rq + LOW_128W) = bb;
  619.     return rq;
  620.   }
  621.  
  622.   if ((bcoeff > 0x0001ed09bead87c0ull)
  623.       || ((bcoeff == 0x0001ed09bead87c0ull)
  624.           && (bb > 0x378d8e63ffffffffull))) {
  625.     // significand is non-canonical
  626.     bcoeff = 0ull;
  627.     bb = 0ull;
  628.   }
  629.  
  630.   *((UINT64 *) & qcoeff + 1) = bcoeff;
  631.   *((UINT64 *) & qcoeff + 0) = bb;
  632.  
  633.   get_BID128_fast (&res, sign, exp, qcoeff);
  634.  
  635.   BID_SWAP128 (res);
  636.   return res;
  637. }
  638.  
  639. UINT32
  640. bid32_canonize (UINT32 ba) {
  641.   FPSC bidrnd;
  642.   unsigned int rnd = 0;
  643.  
  644.   UINT32 res;
  645.   UINT32 sign, comb, exp;
  646.   UINT32 trailing;
  647.   UINT32 bcoeff;
  648.  
  649.   sign = (ba & 0x80000000);
  650.   comb = (ba & 0x7ff00000) >> 20;
  651.   trailing = (ba & 0x000fffff);
  652.  
  653.   if ((comb & 0x600) == 0x600) {        // G0..G1 = 11 -> exp is G2..G11
  654.     exp = (comb >> 1) & 0xff;
  655.     bcoeff = ((8 + (comb & 1)) << 20) | trailing;
  656.   } else {
  657.     exp = (comb >> 3) & 0xff;
  658.     bcoeff = ((comb & 7) << 20) | trailing;
  659.   }
  660.  
  661.   if ((comb & 0x7c0) == 0x7c0) {        //NaN
  662.     ba &= 0xfe0fffff;   // make exponent 0
  663.     bcoeff &= 0x000fffff;       // NaN payloat is only T.    
  664.     if (bcoeff >= 1000000)
  665.       ba &= 0xfff00000; //treat non-canonical significand
  666.     return ba;
  667.   } else if ((comb & 0x780) == 0x780) { //Inf
  668.     ba &= 0xf8000000;   // make exponent and significand 0
  669.     return ba;
  670.   }
  671.  
  672.   if (bcoeff >= 10000000)
  673.     bcoeff = 0;
  674.   rnd = bidrnd = ROUNDING_TO_NEAREST;
  675.   res = get_BID32 (sign, exp, bcoeff, rnd, &bidrnd);
  676.   return res;
  677. }
  678.  
  679. UINT64
  680. bid64_canonize (UINT64 ba) {
  681.   UINT64 res;
  682.   UINT64 sign, comb, exp;
  683.   UINT64 trailing;
  684.   UINT64 bcoeff;
  685.  
  686.   sign = (ba & 0x8000000000000000ull);
  687.   comb = (ba & 0x7ffc000000000000ull) >> 50;
  688.   trailing = (ba & 0x0003ffffffffffffull);
  689.  
  690.  
  691.   if ((comb & 0x1800) == 0x1800) {      // G0..G1 = 11 -> exp is G2..G11
  692.     exp = (comb >> 1) & 0x3ff;
  693.     bcoeff = ((8 + (comb & 1)) << 50) | trailing;
  694.   } else {
  695.     exp = (comb >> 3) & 0x3ff;
  696.     bcoeff = ((comb & 7) << 50) | trailing;
  697.   }
  698.  
  699.   if ((comb & 0x1f00) == 0x1f00) {      //NaN
  700.     ba &= 0xfe03ffffffffffffULL;        // make exponent 0
  701.     bcoeff &= 0x0003ffffffffffffull;    // NaN payloat is only T.  
  702.     if (bcoeff >= 1000000000000000ull)
  703.       ba &= 0xfe00000000000000ull;      // treat non canonical significand and zero G6-G12
  704.     return ba;
  705.   } else if ((comb & 0x1e00) == 0x1e00) {       //Inf
  706.     ba &= 0xf800000000000000ULL;        // make exponent and significand 0
  707.     return ba;
  708.   }
  709.  
  710.   if (bcoeff >= 10000000000000000ull) {
  711.     bcoeff = 0ull;
  712.   }
  713.   res = very_fast_get_BID64 (sign, exp, bcoeff);
  714.   return res;
  715. }
  716.  
  717. UINT128
  718. bid128_canonize (UINT128 bq) {
  719.   UINT128 res;
  720.   UINT64 sign, comb, exp;
  721.   UINT64 trailing;
  722.   UINT64 bcoeff;
  723.  
  724.   UINT128 rq;
  725.   UINT128 qcoeff;
  726.   UINT64 ba, bb;
  727.  
  728.   ba = *((UINT64 *) & bq + HIGH_128W);
  729.   bb = *((UINT64 *) & bq + LOW_128W);
  730.  
  731.   sign = (ba & 0x8000000000000000ull);
  732.   comb = (ba & 0x7fffc00000000000ull) >> 46;
  733.   trailing = (ba & 0x00003fffffffffffull);
  734.  
  735.   if ((comb & 0x18000) == 0x18000) {    // G0..G1 = 11 -> exp is G2..G11
  736.     exp = (comb >> 1) & 0x3fff;
  737.     bcoeff = ((8 + (comb & 1)) << 46) | trailing;
  738.   } else {
  739.     exp = (comb >> 3) & 0x3fff;
  740.     bcoeff = ((comb & 7) << 46) | trailing;
  741.   }
  742.  
  743.   if ((comb & 0x1f000) == 0x1f000) {    //NaN
  744.     ba &= 0xfe003fffffffffffULL;        // make exponent 0
  745.     bcoeff &= 0x00003fffffffffffull;    // NaN payload is only T.  
  746.  
  747.     if ((bcoeff > 0x0000314dc6448d93ULL) ||     // significand is non-canonical
  748.         ((bcoeff == 0x0000314dc6448d93ULL)
  749.          && (bb >= 0x38c15b0a00000000ULL))
  750.         // significand is non-canonical
  751.       ) {
  752.       bcoeff = 0ull;
  753.       ba &= ~0x00003fffffffffffull;
  754.       bb = 0ull;
  755.     }
  756.     *((UINT64 *) & rq + HIGH_128W) = ba;
  757.     *((UINT64 *) & rq + LOW_128W) = bb;
  758.     return rq;
  759.   } else if ((comb & 0x1e000) == 0x1e000) {     //Inf
  760.     ba &= 0xf800000000000000ULL;        // make exponent and significand 0
  761.     bb = 0;
  762.     *((UINT64 *) & rq + HIGH_128W) = ba;
  763.     *((UINT64 *) & rq + LOW_128W) = bb;
  764.     return rq;
  765.   }
  766.  
  767.   if ((bcoeff > 0x0001ed09bead87c0ull) ||       // significand is non-canonical
  768.       ((bcoeff == 0x0001ed09bead87c0ull)
  769.        && (bb > 0x378d8e63ffffffffull))
  770.       // significand is non-canonical
  771.     ) {
  772.     bcoeff = 0ull;
  773.     bb = 0ull;
  774.   }
  775.  
  776.   *((UINT64 *) & qcoeff + 1) = bcoeff;
  777.   *((UINT64 *) & qcoeff + 0) = bb;
  778.  
  779.   get_BID128_fast (&res, sign, exp, qcoeff);
  780.   BID_SWAP128 (res);
  781.   return res;
  782. }
  783.