/* Software floating-point emulation. Convert a _BitInt to _Decimal128. Copyright (C) 2023 Free Software Foundation, Inc. This file is part of GCC. GCC is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation; either version 3, or (at your option) any later version. GCC is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details. Under Section 7 of GPL version 3, you are granted additional permissions described in the GCC Runtime Library Exception, version 3.1, as published by the Free Software Foundation. You should have received a copy of the GNU General Public License and a copy of the GCC Runtime Library Exception along with this program; see the files COPYING3 and COPYING.RUNTIME respectively. If not, see . */ #include "soft-fp.h" #include "bitint.h" #ifdef __BITINT_MAXWIDTH__ extern _Decimal128 __bid_floatbitinttd (const UBILtype *, SItype); _Decimal128 __bid_floatbitinttd (const UBILtype *i, SItype iprec) { iprec = bitint_reduce_prec (&i, iprec); USItype aiprec = iprec < 0 ? -iprec : iprec; USItype in = (aiprec + BIL_TYPE_SIZE - 1) / BIL_TYPE_SIZE; USItype idx = BITINT_END (0, in - 1); UBILtype msb = i[idx]; UDItype mantissahi, mantissalo; SItype exponent = 0; UBILtype inexact = 0; union { _Decimal128 d; UDItype u[2]; } u, ui; if (aiprec % BIL_TYPE_SIZE) { if (iprec > 0) msb &= ((UBILtype) 1 << (aiprec % BIL_TYPE_SIZE)) - 1; else msb |= (UBILtype) -1 << (aiprec % BIL_TYPE_SIZE); } if (iprec < 0) { SItype n = sizeof (0ULL) * __CHAR_BIT__ + 1 - __builtin_clzll (~msb); aiprec = (in - 1) * BIL_TYPE_SIZE + n; } else if (msb == 0) aiprec = 1; else { SItype n = sizeof (0ULL) * __CHAR_BIT__ - __builtin_clzll (msb); aiprec = (in - 1) * BIL_TYPE_SIZE + n; } /* Number of bits in (_BitInt(32768)) 9999999999999999999999999999999999e+6111DL. */ if (aiprec > 20414 + (iprec < 0)) { ovf: if (iprec < 0) u.d = -9000000000000000000000000000000000e+6111DL; else u.d = 9000000000000000000000000000000000e+6111DL; __asm ("" : "+g" (u.d)); u.d += u.d; __asm ("" : "+g" (u.d)); goto done; } /* Bit precision of 9999999999999999999999999999999999uwb. */ if (aiprec >= 113) { USItype pow10_limbs, q_limbs, q2_limbs, j, k; USItype exp_bits = 0, e; UBILtype *buf; /* First do a possibly large divide smaller enough such that we only need to check remainder for 0 or non-0 and then we'll do further division. */ if (aiprec >= 113 + 4 + 10) { exp_bits = ((aiprec - 113 - 4) * (UDItype) 30) / 299; exponent = exp_bits * 3; /* Upper estimate for pow10 (exponent) bits. */ exp_bits = exp_bits * 10 - exp_bits / 30; } pow10_limbs = (exp_bits + BIL_TYPE_SIZE - 1) / BIL_TYPE_SIZE; /* 127 is the highest number of quotient bits needed on aiprec range of [127, 20414]. E.g. if aiprec is 20409, exponent will be 6105 and exp_bits 20283. 20409 - 20283 + 1 is 127. */ q_limbs = (127 + BIL_TYPE_SIZE - 1) / BIL_TYPE_SIZE; q2_limbs = 128 / BIL_TYPE_SIZE; buf = __builtin_alloca ((q_limbs + pow10_limbs * 2 + q2_limbs + 2) * sizeof (UBILtype)); if (exponent) { __bid_pow10bitint (buf + q_limbs, exp_bits, exponent); __divmodbitint4 (buf, q_limbs * BIL_TYPE_SIZE, buf + q_limbs + pow10_limbs, pow10_limbs * BIL_TYPE_SIZE, i, iprec < 0 ? -aiprec : aiprec, buf + q_limbs, exp_bits); if (iprec < 0) bitint_negate (buf + BITINT_END (q_limbs - 1, 0), buf + BITINT_END (q_limbs - 1, 0), q_limbs); inexact = buf[q_limbs + pow10_limbs]; for (j = 1; j < pow10_limbs; ++j) inexact |= buf[q_limbs + pow10_limbs + j]; } else { __builtin_memcpy (buf + BITINT_END (q_limbs - in + 1, 0), i, (in - 1) * sizeof (UBILtype)); buf[BITINT_END (q_limbs - in, in - 1)] = msb; if (iprec < 0) bitint_negate (buf + BITINT_END (q_limbs - 1, 0), buf + BITINT_END (q_limbs - 1, 0), in); if (q_limbs > in) __builtin_memset (buf + BITINT_END (0, in), '\0', (q_limbs - in) * sizeof (UBILtype)); } e = 0; for (j = 3; j; ) { USItype eprev = e; __bid_pow10bitint (buf + q_limbs + pow10_limbs * 2 + 1, 128, 33 + e + j); for (k = BITINT_END (0, q_limbs - 1); k != BITINT_END (q_limbs, (USItype) -1); k -= BITINT_INC) if (buf[k] > buf[q_limbs + pow10_limbs * 2 + 1 + k]) { e += j; break; } else if (buf[k] < buf[q_limbs + pow10_limbs * 2 + 1 + k]) break; if (k == BITINT_END (q_limbs, (USItype) -1)) e += j; if (j == 2 && e != eprev) break; else --j; } exponent += e; if (exponent > 6111) goto ovf; if (e) { UBILtype rem, half; __bid_pow10bitint (buf + q_limbs + pow10_limbs * 2, BIL_TYPE_SIZE, e); __divmodbitint4 (buf + q_limbs + pow10_limbs * 2 + 1, q2_limbs * BIL_TYPE_SIZE, buf + q_limbs + pow10_limbs * 2 + 1 + q2_limbs, BIL_TYPE_SIZE, buf, q_limbs * BIL_TYPE_SIZE, buf + q_limbs + pow10_limbs * 2, BIL_TYPE_SIZE); half = buf[q_limbs + pow10_limbs * 2] / 2; rem = buf[q_limbs + pow10_limbs * 2 + 1 + q2_limbs]; if (inexact) { /* If first division discovered some non-0 digits and this second division is by 10, e.g. for XXXXXX5499999999999 or XXXXXX5000000000001 if first division is by 10^12 and second by 10^1, doing rem |= 1 wouldn't change the 5. Similarly for rem 4 doing rem |= 1 would change it to 5, but we don't want to change it in that case. */ if (e == 1) { if (rem == 5) rem = 6; else if (rem != 4) rem |= 1; } else rem |= 1; } /* Set inexact to 0, 1, 2, 3 depending on if remainder of the divisions is exact 0, smaller than 10^exponent / 2, exactly 10^exponent / 2 or greater than that. */ if (rem >= half) inexact = 2 + (rem > half); else inexact = (rem != 0); #if BIL_TYPE_SIZE == 64 mantissahi = buf[q_limbs + pow10_limbs * 2 + 1 + BITINT_END (0, 1)]; mantissalo = buf[q_limbs + pow10_limbs * 2 + 1 + BITINT_END (1, 0)]; #else mantissahi = ((UDItype) buf[q_limbs + pow10_limbs * 2 + 1 + BITINT_END (0, 3)] << 32 | buf[q_limbs + pow10_limbs * 2 + 1 + BITINT_END (1, 2)]); mantissalo = ((UDItype) buf[q_limbs + pow10_limbs * 2 + 1 + BITINT_END (2, 1)] << 32 | buf[q_limbs + pow10_limbs * 2 + 1 + BITINT_END (3, 0)]); #endif } else { #if BIL_TYPE_SIZE == 64 mantissahi = buf[BITINT_END (0, 1)]; mantissalo = buf[BITINT_END (1, 0)]; #else mantissahi = ((UDItype) buf[BITINT_END (0, 3)] << 32 | buf[BITINT_END (1, 2)]); mantissalo = ((UDItype) buf[BITINT_END (2, 1)] << 32 | buf[BITINT_END (3, 0)]); #endif } } else { mantissahi = iprec < 0 ? -1 : 0; #if BIL_TYPE_SIZE == 64 if (in == 1) mantissalo = msb; else { mantissahi = msb; mantissalo = i[BITINT_END (1, 0)]; } #else if (in <= 2) { if (in == 1) mantissalo = iprec < 0 ? (UDItype) (BILtype) msb : (UDItype) msb; else mantissalo = (UDItype) msb << 32 | i[BITINT_END (1, 0)]; } else { if (in == 3) mantissahi = iprec < 0 ? (UDItype) (BILtype) msb : (UDItype) msb; else mantissahi = (UDItype) msb << 32 | i[BITINT_END (1, 2)]; mantissalo = ((UDItype) i[BITINT_END (in - 2, 1)] << 32 | i[BITINT_END (in - 1, 0)]); } #endif if (iprec < 0) mantissahi = ~mantissahi + __builtin_add_overflow (~mantissalo, 1, &mantissalo); } exponent += 6176; u.u[__FLOAT_WORD_ORDER__ != __ORDER_BIG_ENDIAN__] = ((((UDItype) (iprec < 0)) << 63) | (((UDItype) exponent) << 49) | mantissahi); u.u[__FLOAT_WORD_ORDER__ == __ORDER_BIG_ENDIAN__] = mantissalo; if (inexact) { ui.u[__FLOAT_WORD_ORDER__ != __ORDER_BIG_ENDIAN__] = (((UDItype) (iprec < 0)) << 63) | (((UDItype) exponent - 1) << 49); ui.u[__FLOAT_WORD_ORDER__ == __ORDER_BIG_ENDIAN__] = inexact + 3; __asm ("" : "+g" (u.d)); __asm ("" : "+g" (ui.d)); u.d += ui.d; __asm ("" : "+g" (u.d)); } done: return u.d; } #endif