libgcc/config/libbid/bid64_to_bid128.c

   1 /* Copyright (C) 2007  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 2, or (at your option) any later
   8 version.
   9
  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.)
  18
  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.
  23
  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.  */
  28
  29 #define BID_128RES
  30 #include "bid_internal.h"
  31
  32 /*
  33  * Takes a BID64 as input and converts it to a BID128 and returns it.
  34  */
  35 TYPE0_FUNCTION_ARGTYPE1_NORND (UINT128, bid64_to_bid128, UINT64, x)
  36
  37      UINT128 new_coeff, res;
  38      UINT64 sign_x;
  39      int exponent_x;
  40      UINT64 coefficient_x;
  41
  42 if (!unpack_BID64 (&sign_x, &exponent_x, &coefficient_x, x)) {
  43 if (((x) << 1) >= 0xf000000000000000ull) {
  44 #ifdef SET_STATUS_FLAGS
  45   if (((x) & SNAN_MASK64) == SNAN_MASK64)       // sNaN
  46     __set_status_flags (pfpsf, INVALID_EXCEPTION);
  47 #endif
  48   res.w[0] = (coefficient_x & 0x0003ffffffffffffull);
  49   __mul_64x64_to_128 (res, res.w[0], power10_table_128[18].w[0]);
  50   res.w[1] |= ((coefficient_x) & 0xfc00000000000000ull);
  51   BID_RETURN (res);
  52 }
  53 }
  54
  55 new_coeff.w[0] = coefficient_x;
  56 new_coeff.w[1] = 0;
  57 get_BID128_very_fast (&res, sign_x,
  58                       exponent_x + DECIMAL_EXPONENT_BIAS_128 -
  59                       DECIMAL_EXPONENT_BIAS, new_coeff);
  60 BID_RETURN (res);
  61 }       // convert_bid64_to_bid128
  62
  63
  64
  65 /*
  66  * Takes a BID128 as input and converts it to a BID64 and returns it.
  67  */
  68 #if DECIMAL_CALL_BY_REFERENCE
  69
  70 void
  71 bid128_to_bid64 (UINT64 * pres,
  72                  UINT128 *
  73                  px _RND_MODE_PARAM _EXC_FLAGS_PARAM _EXC_MASKS_PARAM
  74                  _EXC_INFO_PARAM) {
  75   UINT128 x = *px;
  76 #else
  77
  78 UINT64
  79 bid128_to_bid64 (UINT128 x _RND_MODE_PARAM _EXC_FLAGS_PARAM
  80                  _EXC_MASKS_PARAM _EXC_INFO_PARAM) {
  81 #endif
  82   UINT128 CX, T128, TP128, Qh, Ql, Qh1, Stemp, Tmp, Tmp1, CX1;
  83   UINT64 sign_x, carry, cy, res;
  84   SINT64 D;
  85   int_float f64, fx;
  86   int exponent_x, extra_digits, amount, bin_expon_cx;
  87   unsigned rmode, status, uf_check = 0;
  88
  89 #if DECIMAL_CALL_BY_REFERENCE
  90 #if !DECIMAL_GLOBAL_ROUNDING
  91   _IDEC_round rnd_mode = *prnd_mode;
  92 #endif
  93 #endif
  94
  95   BID_SWAP128 (x);
  96   // unpack arguments, check for NaN or Infinity or 0
  97   if (!unpack_BID128_value (&sign_x, &exponent_x, &CX, x)) {
  98     if ((x.w[1] << 1) >= 0xf000000000000000ull) {
  99       Tmp.w[1] = (CX.w[1] & 0x00003fffffffffffull);
 100       Tmp.w[0] = CX.w[0];
 101       TP128 = reciprocals10_128[18];
 102       __mul_128x128_full (Qh, Ql, Tmp, TP128);
 103       amount = recip_scale[18];
 104       __shr_128 (Tmp, Qh, amount);
 105       res = (CX.w[1] & 0xfc00000000000000ull) | Tmp.w[0];
 106 #ifdef SET_STATUS_FLAGS
 107       if ((x.w[1] & SNAN_MASK64) == SNAN_MASK64)        // sNaN
 108         __set_status_flags (pfpsf, INVALID_EXCEPTION);
 109 #endif
 110       BID_RETURN_VAL (res);
 111     }
 112     exponent_x =
 113       exponent_x - DECIMAL_EXPONENT_BIAS_128 + DECIMAL_EXPONENT_BIAS;
 114     if (exponent_x < 0) {
 115       res = sign_x;
 116       BID_RETURN_VAL (res);
 117     }
 118     if (exponent_x > DECIMAL_MAX_EXPON_64)
 119       exponent_x = DECIMAL_MAX_EXPON_64;
 120     res = sign_x | (((UINT64) exponent_x) << 53);
 121     BID_RETURN_VAL (res);
 122   }
 123
 124   if (CX.w[1] || (CX.w[0] >= 10000000000000000ull)) {
 125     // find number of digits in coefficient
 126     // 2^64
 127     f64.i = 0x5f800000;
 128     // fx ~ CX
 129     fx.d = (float) CX.w[1] * f64.d + (float) CX.w[0];
 130     bin_expon_cx = ((fx.i >> 23) & 0xff) - 0x7f;
 131     extra_digits = estimate_decimal_digits[bin_expon_cx] - 16;
 132     // scale = 38-estimate_decimal_digits[bin_expon_cx];
 133     D = CX.w[1] - power10_index_binexp_128[bin_expon_cx].w[1];
 134     if (D > 0
 135         || (!D
 136             && CX.w[0] >= power10_index_binexp_128[bin_expon_cx].w[0]))
 137       extra_digits++;
 138
 139     exponent_x += extra_digits;
 140
 141 #ifndef IEEE_ROUND_NEAREST_TIES_AWAY
 142 #ifndef IEEE_ROUND_NEAREST
 143     rmode = rnd_mode;
 144     if (sign_x && (unsigned) (rmode - 1) < 2)
 145       rmode = 3 - rmode;
 146 #else
 147     rmode = 0;
 148 #endif
 149 #else
 150     rmode = 0;
 151 #endif
 152     if (exponent_x < DECIMAL_EXPONENT_BIAS_128 - DECIMAL_EXPONENT_BIAS) {
 153       uf_check = 1;
 154       if (-extra_digits + exponent_x - DECIMAL_EXPONENT_BIAS_128 +
 155           DECIMAL_EXPONENT_BIAS + 35 >= 0) {
 156         if (exponent_x ==
 157             DECIMAL_EXPONENT_BIAS_128 - DECIMAL_EXPONENT_BIAS - 1) {
 158           T128 = round_const_table_128[rmode][extra_digits];
 159           __add_carry_out (CX1.w[0], carry, T128.w[0], CX.w[0]);
 160           CX1.w[1] = CX.w[1] + T128.w[1] + carry;
 161           if (__unsigned_compare_ge_128
 162               (CX1, power10_table_128[extra_digits + 16]))
 163             uf_check = 0;
 164         }
 165         extra_digits =
 166           extra_digits + DECIMAL_EXPONENT_BIAS_128 -
 167           DECIMAL_EXPONENT_BIAS - exponent_x;
 168         exponent_x = DECIMAL_EXPONENT_BIAS_128 - DECIMAL_EXPONENT_BIAS;
 169         //uf_check = 2;
 170       } else
 171         rmode = ROUNDING_TO_ZERO;
 172     }
 173
 174     T128 = round_const_table_128[rmode][extra_digits];
 175     __add_carry_out (CX.w[0], carry, T128.w[0], CX.w[0]);
 176     CX.w[1] = CX.w[1] + T128.w[1] + carry;
 177
 178     TP128 = reciprocals10_128[extra_digits];
 179     __mul_128x128_full (Qh, Ql, CX, TP128);
 180     amount = recip_scale[extra_digits];
 181
 182     if (amount >= 64) {
 183       CX.w[0] = Qh.w[1] >> (amount - 64);
 184       CX.w[1] = 0;
 185     } else {
 186       __shr_128 (CX, Qh, amount);
 187     }
 188
 189 #ifndef IEEE_ROUND_NEAREST_TIES_AWAY
 190 #ifndef IEEE_ROUND_NEAREST
 191     if (!(rmode))
 192 #endif
 193       if (CX.w[0] & 1) {
 194         // check whether fractional part of initial_P/10^ed1 is exactly .5
 195
 196         // get remainder
 197         __shl_128_long (Qh1, Qh, (128 - amount));
 198
 199         if (!Qh1.w[1] && !Qh1.w[0]
 200             && (Ql.w[1] < reciprocals10_128[extra_digits].w[1]
 201                 || (Ql.w[1] == reciprocals10_128[extra_digits].w[1]
 202                     && Ql.w[0] < reciprocals10_128[extra_digits].w[0]))) {
 203           CX.w[0]--;
 204         }
 205       }
 206 #endif
 207
 208     {
 209       status = INEXACT_EXCEPTION;
 210       // get remainder
 211       __shl_128_long (Qh1, Qh, (128 - amount));
 212
 213       switch (rmode) {
 214       case ROUNDING_TO_NEAREST:
 215       case ROUNDING_TIES_AWAY:
 216         // test whether fractional part is 0
 217         if (Qh1.w[1] == 0x8000000000000000ull && (!Qh1.w[0])
 218             && (Ql.w[1] < reciprocals10_128[extra_digits].w[1]
 219                 || (Ql.w[1] == reciprocals10_128[extra_digits].w[1]
 220                     && Ql.w[0] < reciprocals10_128[extra_digits].w[0])))
 221           status = EXACT_STATUS;
 222         break;
 223       case ROUNDING_DOWN:
 224       case ROUNDING_TO_ZERO:
 225         if ((!Qh1.w[1]) && (!Qh1.w[0])
 226             && (Ql.w[1] < reciprocals10_128[extra_digits].w[1]
 227                 || (Ql.w[1] == reciprocals10_128[extra_digits].w[1]
 228                     && Ql.w[0] < reciprocals10_128[extra_digits].w[0])))
 229           status = EXACT_STATUS;
 230         break;
 231       default:
 232         // round up
 233         __add_carry_out (Stemp.w[0], cy, Ql.w[0],
 234                          reciprocals10_128[extra_digits].w[0]);
 235         __add_carry_in_out (Stemp.w[1], carry, Ql.w[1],
 236                             reciprocals10_128[extra_digits].w[1], cy);
 237         __shr_128_long (Qh, Qh1, (128 - amount));
 238         Tmp.w[0] = 1;
 239         Tmp.w[1] = 0;
 240         __shl_128_long (Tmp1, Tmp, amount);
 241         Qh.w[0] += carry;
 242         if (Qh.w[0] < carry)
 243           Qh.w[1]++;
 244         if (__unsigned_compare_ge_128 (Qh, Tmp1))
 245           status = EXACT_STATUS;
 246       }
 247
 248       if (status != EXACT_STATUS) {
 249         if (uf_check)
 250           status |= UNDERFLOW_EXCEPTION;
 251 #ifdef SET_STATUS_FLAGS
 252         __set_status_flags (pfpsf, status);
 253 #endif
 254       }
 255
 256
 257     }
 258
 259   }
 260
 261   res =
 262     get_BID64 (sign_x,
 263                exponent_x - DECIMAL_EXPONENT_BIAS_128 +
 264                DECIMAL_EXPONENT_BIAS, CX.w[0], rnd_mode, pfpsf);
 265   BID_RETURN_VAL (res);
 266
 267 }