libgcc/soft-fp/floatbitintsd.c

   1 /* Software floating-point emulation.
   2    Convert a _BitInt to _Decimal32.
   3
   4    Copyright (C) 2023 Free Software Foundation, Inc.
   5
   6 This file is part of GCC.
   7
   8 GCC is free software; you can redistribute it and/or modify it under
   9 the terms of the GNU General Public License as published by the Free
  10 Software Foundation; either version 3, or (at your option) any later
  11 version.
  12
  13 GCC is distributed in the hope that it will be useful, but WITHOUT ANY
  14 WARRANTY; without even the implied warranty of MERCHANTABILITY or
  15 FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
  16 for more details.
  17
  18 Under Section 7 of GPL version 3, you are granted additional
  19 permissions described in the GCC Runtime Library Exception, version
  20 3.1, as published by the Free Software Foundation.
  21
  22 You should have received a copy of the GNU General Public License and
  23 a copy of the GCC Runtime Library Exception along with this program;
  24 see the files COPYING3 and COPYING.RUNTIME respectively.  If not, see
  25 <http://www.gnu.org/licenses/>.  */
  26
  27 #include "soft-fp.h"
  28 #include "bitint.h"
  29
  30 #ifdef __BITINT_MAXWIDTH__
  31 extern _Decimal32 __bid_floatbitintsd (const UBILtype *, SItype);
  32
  33 _Decimal32
  34 __bid_floatbitintsd (const UBILtype *i, SItype iprec)
  35 {
  36   iprec = bitint_reduce_prec (&i, iprec);
  37   USItype aiprec = iprec < 0 ? -iprec : iprec;
  38   USItype in = (aiprec + BIL_TYPE_SIZE - 1) / BIL_TYPE_SIZE;
  39   USItype idx = BITINT_END (0, in - 1);
  40   UBILtype msb = i[idx];
  41   USItype mantissa;
  42   SItype exponent = 0;
  43   UBILtype inexact = 0;
  44   union { _Decimal32 d; USItype u; } u, ui;
  45   if (aiprec % BIL_TYPE_SIZE)
  46     {
  47       if (iprec > 0)
  48         msb &= ((UBILtype) 1 << (aiprec % BIL_TYPE_SIZE)) - 1;
  49       else
  50         msb |= (UBILtype) -1 << (aiprec % BIL_TYPE_SIZE);
  51     }
  52   if (iprec < 0)
  53     {
  54       SItype n = sizeof (0ULL) * __CHAR_BIT__ + 1 - __builtin_clzll (~msb);
  55       aiprec = (in - 1) * BIL_TYPE_SIZE + n;
  56     }
  57   else if (msb == 0)
  58     aiprec = 1;
  59   else
  60     {
  61       SItype n = sizeof (0ULL) * __CHAR_BIT__ - __builtin_clzll (msb);
  62       aiprec = (in - 1) * BIL_TYPE_SIZE + n;
  63     }
  64   /* Number of bits in (_BitInt(2048)) 9999999e+90DF.  */
  65   if (aiprec > 323 + (iprec < 0))
  66     {
  67     ovf:
  68       if (iprec < 0)
  69         u.d = -9000000e+90DF;
  70       else
  71         u.d = 9000000e+90DF;
  72       __asm ("" : "+g" (u.d));
  73       u.d += u.d;
  74       __asm ("" : "+g" (u.d));
  75       goto done;
  76     }
  77   /* Bit precision of 9999999uwb.  */
  78   if (aiprec >= 24)
  79     {
  80       USItype pow10_limbs, q_limbs, q2_limbs, j;
  81       USItype exp_bits = 0, e;
  82       UDItype m;
  83       UBILtype *buf;
  84       /* First do a possibly large divide smaller enough such that
  85          we only need to check remainder for 0 or non-0 and then
  86          we'll do further division.  */
  87       if (aiprec >= 24 + 4 + 10)
  88         {
  89           exp_bits = (aiprec - 24 - 4) / 10;
  90           exponent = exp_bits * 3;
  91           /* Upper estimate for pow10 (exponent) bits.  */
  92           exp_bits = exp_bits * 10 - exp_bits / 30;
  93         }
  94       pow10_limbs = (exp_bits + BIL_TYPE_SIZE - 1) / BIL_TYPE_SIZE;
  95       /* 38 is the highest number of quotient bits needed on
  96          aiprec range of [38, 323].  E.g. if aiprec is 317,
  97          exponent will be 84 and exp_bits 280.  317 - 280 + 1
  98          is 38.  */
  99       q_limbs = (38 + BIL_TYPE_SIZE - 1) / BIL_TYPE_SIZE;
 100       q2_limbs = (32 + BIL_TYPE_SIZE - 1) / BIL_TYPE_SIZE;
 101       buf = __builtin_alloca ((q_limbs + pow10_limbs * 2 + q2_limbs + 2)
 102                               * sizeof (UBILtype));
 103       if (exponent)
 104         {
 105           __bid_pow10bitint (buf + q_limbs, exp_bits, exponent);
 106           __divmodbitint4 (buf, q_limbs * BIL_TYPE_SIZE,
 107                            buf + q_limbs + pow10_limbs,
 108                            pow10_limbs * BIL_TYPE_SIZE,
 109                            i, iprec < 0 ? -aiprec : aiprec,
 110                            buf + q_limbs, exp_bits);
 111           if (iprec < 0)
 112             bitint_negate (buf + BITINT_END (q_limbs - 1, 0),
 113                            buf + BITINT_END (q_limbs - 1, 0), q_limbs);
 114           inexact = buf[q_limbs + pow10_limbs];
 115           for (j = 1; j < pow10_limbs; ++j)
 116             inexact |= buf[q_limbs + pow10_limbs + 1];
 117         }
 118       else
 119         {
 120           __builtin_memcpy (buf + BITINT_END (q_limbs - in + 1, 0), i,
 121                             (in - 1) * sizeof (UBILtype));
 122           buf[BITINT_END (q_limbs - in, in - 1)] = msb;
 123           if (iprec < 0)
 124             bitint_negate (buf + BITINT_END (q_limbs - 1, 0),
 125                            buf + BITINT_END (q_limbs - 1, 0), in);
 126           if (q_limbs > in)
 127             __builtin_memset (buf + BITINT_END (0, in), '\0',
 128                               (q_limbs - in) * sizeof (UBILtype));
 129         }
 130       e = 0;
 131 #if BIL_TYPE_SIZE == 64
 132       m = buf[0];
 133 #elif BIL_TYPE_SIZE == 32
 134       m = ((UDItype) buf[BITINT_END (0, 1)] << 32) | buf[BITINT_END (1, 0)];
 135 #else
 136 # error Unsupported BIL_TYPE_SIZE
 137 #endif
 138       if (m >= (UDItype) 10000000000)
 139         {
 140           if (m >= (UDItype) 100000000000)
 141             e = 5;
 142           else
 143             e = 4;
 144         }
 145       else if (m >= (UDItype) 100000000)
 146         {
 147           if (m >= (UDItype) 1000000000)
 148             e = 3;
 149           else
 150             e = 2;
 151         }
 152       else if (m >= (UDItype) 10000000)
 153         e = 1;
 154       exponent += e;
 155       if (exponent > 90)
 156         goto ovf;
 157       if (e)
 158         {
 159           UBILtype rem, half;
 160           __bid_pow10bitint (buf + q_limbs + pow10_limbs * 2,
 161                              BIL_TYPE_SIZE, e);
 162           __divmodbitint4 (buf + q_limbs + pow10_limbs * 2 + 1,
 163                            q2_limbs * BIL_TYPE_SIZE,
 164                            buf + q_limbs + pow10_limbs * 2 + 1 + q2_limbs,
 165                            BIL_TYPE_SIZE,
 166                            buf, q_limbs * BIL_TYPE_SIZE,
 167                            buf + q_limbs + pow10_limbs * 2, BIL_TYPE_SIZE);
 168           half = buf[q_limbs + pow10_limbs * 2] / 2;
 169           rem = buf[q_limbs + pow10_limbs * 2 + 1 + q2_limbs];
 170           if (inexact)
 171             {
 172               /* If first division discovered some non-0 digits
 173                  and this second division is by 10, e.g.
 174                  for XXXXXX5499999999999 or XXXXXX5000000000001
 175                  if first division is by 10^12 and second by 10^1,
 176                  doing rem |= 1 wouldn't change the 5.  Similarly
 177                  for rem 4 doing rem |= 1 would change it to 5,
 178                  but we don't want to change it in that case.  */
 179               if (e == 1)
 180                 {
 181                   if (rem == 5)
 182                     rem = 6;
 183                   else if (rem != 4)
 184                     rem |= 1;
 185                 }
 186               else
 187                 rem |= 1;
 188             }
 189           /* Set inexact to 0, 1, 2, 3 depending on if remainder
 190              of the divisions is exact 0, smaller than 10^exponent / 2,
 191              exactly 10^exponent / 2 or greater than that.  */
 192           if (rem >= half)
 193             inexact = 2 + (rem > half);
 194           else
 195             inexact = (rem != 0);
 196           mantissa = buf[q_limbs + pow10_limbs * 2 + 1];
 197         }
 198       else
 199 #if BIL_TYPE_SIZE == 64
 200         mantissa = buf[0];
 201 #else
 202         mantissa = buf[BITINT_END (1, 0)];
 203 #endif
 204     }
 205   else
 206     {
 207       mantissa = msb;
 208       if (iprec < 0)
 209         mantissa = -mantissa;
 210     }
 211
 212   exponent += 101;
 213   if (mantissa >= (USItype) 0x800000)
 214     u.u = (((((iprec < 0) << 2) | (USItype) 3) << 29)
 215            | (((USItype) exponent) << 21)
 216            | (mantissa ^ (USItype) 0x800000));
 217   else
 218     u.u = ((((USItype) (iprec < 0)) << 31)
 219            | (((USItype) exponent) << 23)
 220            | mantissa);
 221   if (inexact)
 222     {
 223       ui.u = ((((USItype) (iprec < 0)) << 31)
 224               | (((USItype) (exponent - 1)) << 23)
 225               | (inexact + 3));
 226       __asm ("" : "+g" (u.d));
 227       __asm ("" : "+g" (ui.d));
 228       u.d += ui.d;
 229       __asm ("" : "+g" (u.d));
 230     }
 231
 232 done:
 233   return u.d;
 234 }
 235 #endif