1 /* Software floating-point emulation. 2 Convert a _BitInt to _Decimal128. 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 _Decimal128 __bid_floatbitinttd (const UBILtype *, SItype); 32 33 _Decimal128 34 __bid_floatbitinttd (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 UDItype mantissahi, mantissalo; 42 SItype exponent = 0; 43 UBILtype inexact = 0; 44 union { _Decimal128 d; UDItype u[2]; } 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 65 (_BitInt(32768)) 9999999999999999999999999999999999e+6111DL. */ 66 if (aiprec > 20414 + (iprec < 0)) 67 { 68 ovf: 69 if (iprec < 0) 70 u.d = -9000000000000000000000000000000000e+6111DL; 71 else 72 u.d = 9000000000000000000000000000000000e+6111DL; 73 __asm ("" : "+g" (u.d)); 74 u.d += u.d; 75 __asm ("" : "+g" (u.d)); 76 goto done; 77 } 78 /* Bit precision of 9999999999999999999999999999999999uwb. */ 79 if (aiprec >= 113) 80 { 81 USItype pow10_limbs, q_limbs, q2_limbs, j, k; 82 USItype exp_bits = 0, e; 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 >= 113 + 4 + 10) 88 { 89 exp_bits = ((aiprec - 113 - 4) * (UDItype) 30) / 299; 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 /* 127 is the highest number of quotient bits needed on 96 aiprec range of [127, 20414]. E.g. if aiprec is 20409, 97 exponent will be 6105 and exp_bits 20283. 20409 - 20283 + 1 98 is 127. */ 99 q_limbs = (127 + BIL_TYPE_SIZE - 1) / BIL_TYPE_SIZE; 100 q2_limbs = 128 / 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 + j]; 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 for (j = 3; j; ) 132 { 133 USItype eprev = e; 134 __bid_pow10bitint (buf + q_limbs + pow10_limbs * 2 + 1, 135 128, 33 + e + j); 136 for (k = BITINT_END (0, q_limbs - 1); 137 k != BITINT_END (q_limbs, (USItype) -1); k -= BITINT_INC) 138 if (buf[k] > buf[q_limbs + pow10_limbs * 2 + 1 + k]) 139 { 140 e += j; 141 break; 142 } 143 else if (buf[k] < buf[q_limbs + pow10_limbs * 2 + 1 + k]) 144 break; 145 if (k == BITINT_END (q_limbs, (USItype) -1)) 146 e += j; 147 if (j == 2 && e != eprev) 148 break; 149 else 150 --j; 151 } 152 exponent += e; 153 if (exponent > 6111) 154 goto ovf; 155 if (e) 156 { 157 UBILtype rem, half; 158 __bid_pow10bitint (buf + q_limbs + pow10_limbs * 2, 159 BIL_TYPE_SIZE, e); 160 __divmodbitint4 (buf + q_limbs + pow10_limbs * 2 + 1, 161 q2_limbs * BIL_TYPE_SIZE, 162 buf + q_limbs + pow10_limbs * 2 + 1 + q2_limbs, 163 BIL_TYPE_SIZE, 164 buf, q_limbs * BIL_TYPE_SIZE, 165 buf + q_limbs + pow10_limbs * 2, BIL_TYPE_SIZE); 166 half = buf[q_limbs + pow10_limbs * 2] / 2; 167 rem = buf[q_limbs + pow10_limbs * 2 + 1 + q2_limbs]; 168 if (inexact) 169 { 170 /* If first division discovered some non-0 digits 171 and this second division is by 10, e.g. 172 for XXXXXX5499999999999 or XXXXXX5000000000001 173 if first division is by 10^12 and second by 10^1, 174 doing rem |= 1 wouldn't change the 5. Similarly 175 for rem 4 doing rem |= 1 would change it to 5, 176 but we don't want to change it in that case. */ 177 if (e == 1) 178 { 179 if (rem == 5) 180 rem = 6; 181 else if (rem != 4) 182 rem |= 1; 183 } 184 else 185 rem |= 1; 186 } 187 /* Set inexact to 0, 1, 2, 3 depending on if remainder 188 of the divisions is exact 0, smaller than 10^exponent / 2, 189 exactly 10^exponent / 2 or greater than that. */ 190 if (rem >= half) 191 inexact = 2 + (rem > half); 192 else 193 inexact = (rem != 0); 194 #if BIL_TYPE_SIZE == 64 195 mantissahi = buf[q_limbs + pow10_limbs * 2 + 1 + BITINT_END (0, 1)]; 196 mantissalo = buf[q_limbs + pow10_limbs * 2 + 1 + BITINT_END (1, 0)]; 197 #else 198 mantissahi 199 = ((UDItype) 200 buf[q_limbs + pow10_limbs * 2 + 1 + BITINT_END (0, 3)] << 32 201 | buf[q_limbs + pow10_limbs * 2 + 1 + BITINT_END (1, 2)]); 202 mantissalo 203 = ((UDItype) 204 buf[q_limbs + pow10_limbs * 2 + 1 + BITINT_END (2, 1)] << 32 205 | buf[q_limbs + pow10_limbs * 2 + 1 + BITINT_END (3, 0)]); 206 #endif 207 } 208 else 209 { 210 #if BIL_TYPE_SIZE == 64 211 mantissahi = buf[BITINT_END (0, 1)]; 212 mantissalo = buf[BITINT_END (1, 0)]; 213 #else 214 mantissahi = ((UDItype) buf[BITINT_END (0, 3)] << 32 215 | buf[BITINT_END (1, 2)]); 216 mantissalo = ((UDItype) buf[BITINT_END (2, 1)] << 32 217 | buf[BITINT_END (3, 0)]); 218 #endif 219 } 220 } 221 else 222 { 223 mantissahi = iprec < 0 ? -1 : 0; 224 #if BIL_TYPE_SIZE == 64 225 if (in == 1) 226 mantissalo = msb; 227 else 228 { 229 mantissahi = msb; 230 mantissalo = i[BITINT_END (1, 0)]; 231 } 232 #else 233 if (in <= 2) 234 { 235 if (in == 1) 236 mantissalo = iprec < 0 ? (UDItype) (BILtype) msb : (UDItype) msb; 237 else 238 mantissalo = (UDItype) msb << 32 | i[BITINT_END (1, 0)]; 239 } 240 else 241 { 242 if (in == 3) 243 mantissahi = iprec < 0 ? (UDItype) (BILtype) msb : (UDItype) msb; 244 else 245 mantissahi = (UDItype) msb << 32 | i[BITINT_END (1, 2)]; 246 mantissalo = ((UDItype) i[BITINT_END (in - 2, 1)] << 32 247 | i[BITINT_END (in - 1, 0)]); 248 } 249 #endif 250 if (iprec < 0) 251 mantissahi 252 = ~mantissahi + __builtin_add_overflow (~mantissalo, 1, &mantissalo); 253 } 254 255 exponent += 6176; 256 u.u[__FLOAT_WORD_ORDER__ != __ORDER_BIG_ENDIAN__] 257 = ((((UDItype) (iprec < 0)) << 63) 258 | (((UDItype) exponent) << 49) 259 | mantissahi); 260 u.u[__FLOAT_WORD_ORDER__ == __ORDER_BIG_ENDIAN__] = mantissalo; 261 if (inexact) 262 { 263 ui.u[__FLOAT_WORD_ORDER__ != __ORDER_BIG_ENDIAN__] 264 = (((UDItype) (iprec < 0)) << 63) | (((UDItype) exponent - 1) << 49); 265 ui.u[__FLOAT_WORD_ORDER__ == __ORDER_BIG_ENDIAN__] = inexact + 3; 266 __asm ("" : "+g" (u.d)); 267 __asm ("" : "+g" (ui.d)); 268 u.d += ui.d; 269 __asm ("" : "+g" (u.d)); 270 } 271 272 done: 273 return u.d; 274 } 275 #endif 276