1 //===-- APFloat.cpp - Implement APFloat class -----------------------------===// 2 // 3 // The LLVM Compiler Infrastructure 4 // 5 // This file is distributed under the University of Illinois Open Source 6 // License. See LICENSE.TXT for details. 7 // 8 //===----------------------------------------------------------------------===// 9 // 10 // This file implements a class to represent arbitrary precision floating 11 // point values and provide a variety of arithmetic operations on them. 12 // 13 //===----------------------------------------------------------------------===// 14 15 #include "llvm/ADT/APFloat.h" 16 #include "llvm/ADT/APSInt.h" 17 #include "llvm/ADT/FoldingSet.h" 18 #include "llvm/ADT/Hashing.h" 19 #include "llvm/ADT/StringExtras.h" 20 #include "llvm/ADT/StringRef.h" 21 #include "llvm/Support/ErrorHandling.h" 22 #include "llvm/Support/MathExtras.h" 23 #include <cstring> 24 #include <limits.h> 25 26 using namespace llvm; 27 28 /// A macro used to combine two fcCategory enums into one key which can be used 29 /// in a switch statement to classify how the interaction of two APFloat's 30 /// categories affects an operation. 31 /// 32 /// TODO: If clang source code is ever allowed to use constexpr in its own 33 /// codebase, change this into a static inline function. 34 #define PackCategoriesIntoKey(_lhs, _rhs) ((_lhs) * 4 + (_rhs)) 35 36 /* Assumed in hexadecimal significand parsing, and conversion to 37 hexadecimal strings. */ 38 static_assert(integerPartWidth % 4 == 0, "Part width must be divisible by 4!"); 39 40 namespace llvm { 41 42 /* Represents floating point arithmetic semantics. */ 43 struct fltSemantics { 44 /* The largest E such that 2^E is representable; this matches the 45 definition of IEEE 754. */ 46 APFloat::ExponentType maxExponent; 47 48 /* The smallest E such that 2^E is a normalized number; this 49 matches the definition of IEEE 754. */ 50 APFloat::ExponentType minExponent; 51 52 /* Number of bits in the significand. This includes the integer 53 bit. */ 54 unsigned int precision; 55 56 /* Number of bits actually used in the semantics. */ 57 unsigned int sizeInBits; 58 }; 59 60 const fltSemantics APFloat::IEEEhalf = { 15, -14, 11, 16 }; 61 const fltSemantics APFloat::IEEEsingle = { 127, -126, 24, 32 }; 62 const fltSemantics APFloat::IEEEdouble = { 1023, -1022, 53, 64 }; 63 const fltSemantics APFloat::IEEEquad = { 16383, -16382, 113, 128 }; 64 const fltSemantics APFloat::x87DoubleExtended = { 16383, -16382, 64, 80 }; 65 const fltSemantics APFloat::Bogus = { 0, 0, 0, 0 }; 66 67 /* The PowerPC format consists of two doubles. It does not map cleanly 68 onto the usual format above. It is approximated using twice the 69 mantissa bits. Note that for exponents near the double minimum, 70 we no longer can represent the full 106 mantissa bits, so those 71 will be treated as denormal numbers. 72 73 FIXME: While this approximation is equivalent to what GCC uses for 74 compile-time arithmetic on PPC double-double numbers, it is not able 75 to represent all possible values held by a PPC double-double number, 76 for example: (long double) 1.0 + (long double) 0x1p-106 77 Should this be replaced by a full emulation of PPC double-double? */ 78 const fltSemantics APFloat::PPCDoubleDouble = { 1023, -1022 + 53, 53 + 53, 128 }; 79 80 /* A tight upper bound on number of parts required to hold the value 81 pow(5, power) is 82 83 power * 815 / (351 * integerPartWidth) + 1 84 85 However, whilst the result may require only this many parts, 86 because we are multiplying two values to get it, the 87 multiplication may require an extra part with the excess part 88 being zero (consider the trivial case of 1 * 1, tcFullMultiply 89 requires two parts to hold the single-part result). So we add an 90 extra one to guarantee enough space whilst multiplying. */ 91 const unsigned int maxExponent = 16383; 92 const unsigned int maxPrecision = 113; 93 const unsigned int maxPowerOfFiveExponent = maxExponent + maxPrecision - 1; 94 const unsigned int maxPowerOfFiveParts = 2 + ((maxPowerOfFiveExponent * 815) 95 / (351 * integerPartWidth)); 96 } 97 98 /* A bunch of private, handy routines. */ 99 100 static inline unsigned int 101 partCountForBits(unsigned int bits) 102 { 103 return ((bits) + integerPartWidth - 1) / integerPartWidth; 104 } 105 106 /* Returns 0U-9U. Return values >= 10U are not digits. */ 107 static inline unsigned int 108 decDigitValue(unsigned int c) 109 { 110 return c - '0'; 111 } 112 113 /* Return the value of a decimal exponent of the form 114 [+-]ddddddd. 115 116 If the exponent overflows, returns a large exponent with the 117 appropriate sign. */ 118 static int 119 readExponent(StringRef::iterator begin, StringRef::iterator end) 120 { 121 bool isNegative; 122 unsigned int absExponent; 123 const unsigned int overlargeExponent = 24000; /* FIXME. */ 124 StringRef::iterator p = begin; 125 126 assert(p != end && "Exponent has no digits"); 127 128 isNegative = (*p == '-'); 129 if (*p == '-' || *p == '+') { 130 p++; 131 assert(p != end && "Exponent has no digits"); 132 } 133 134 absExponent = decDigitValue(*p++); 135 assert(absExponent < 10U && "Invalid character in exponent"); 136 137 for (; p != end; ++p) { 138 unsigned int value; 139 140 value = decDigitValue(*p); 141 assert(value < 10U && "Invalid character in exponent"); 142 143 value += absExponent * 10; 144 if (absExponent >= overlargeExponent) { 145 absExponent = overlargeExponent; 146 p = end; /* outwit assert below */ 147 break; 148 } 149 absExponent = value; 150 } 151 152 assert(p == end && "Invalid exponent in exponent"); 153 154 if (isNegative) 155 return -(int) absExponent; 156 else 157 return (int) absExponent; 158 } 159 160 /* This is ugly and needs cleaning up, but I don't immediately see 161 how whilst remaining safe. */ 162 static int 163 totalExponent(StringRef::iterator p, StringRef::iterator end, 164 int exponentAdjustment) 165 { 166 int unsignedExponent; 167 bool negative, overflow; 168 int exponent = 0; 169 170 assert(p != end && "Exponent has no digits"); 171 172 negative = *p == '-'; 173 if (*p == '-' || *p == '+') { 174 p++; 175 assert(p != end && "Exponent has no digits"); 176 } 177 178 unsignedExponent = 0; 179 overflow = false; 180 for (; p != end; ++p) { 181 unsigned int value; 182 183 value = decDigitValue(*p); 184 assert(value < 10U && "Invalid character in exponent"); 185 186 unsignedExponent = unsignedExponent * 10 + value; 187 if (unsignedExponent > 32767) { 188 overflow = true; 189 break; 190 } 191 } 192 193 if (exponentAdjustment > 32767 || exponentAdjustment < -32768) 194 overflow = true; 195 196 if (!overflow) { 197 exponent = unsignedExponent; 198 if (negative) 199 exponent = -exponent; 200 exponent += exponentAdjustment; 201 if (exponent > 32767 || exponent < -32768) 202 overflow = true; 203 } 204 205 if (overflow) 206 exponent = negative ? -32768: 32767; 207 208 return exponent; 209 } 210 211 static StringRef::iterator 212 skipLeadingZeroesAndAnyDot(StringRef::iterator begin, StringRef::iterator end, 213 StringRef::iterator *dot) 214 { 215 StringRef::iterator p = begin; 216 *dot = end; 217 while (p != end && *p == '0') 218 p++; 219 220 if (p != end && *p == '.') { 221 *dot = p++; 222 223 assert(end - begin != 1 && "Significand has no digits"); 224 225 while (p != end && *p == '0') 226 p++; 227 } 228 229 return p; 230 } 231 232 /* Given a normal decimal floating point number of the form 233 234 dddd.dddd[eE][+-]ddd 235 236 where the decimal point and exponent are optional, fill out the 237 structure D. Exponent is appropriate if the significand is 238 treated as an integer, and normalizedExponent if the significand 239 is taken to have the decimal point after a single leading 240 non-zero digit. 241 242 If the value is zero, V->firstSigDigit points to a non-digit, and 243 the return exponent is zero. 244 */ 245 struct decimalInfo { 246 const char *firstSigDigit; 247 const char *lastSigDigit; 248 int exponent; 249 int normalizedExponent; 250 }; 251 252 static void 253 interpretDecimal(StringRef::iterator begin, StringRef::iterator end, 254 decimalInfo *D) 255 { 256 StringRef::iterator dot = end; 257 StringRef::iterator p = skipLeadingZeroesAndAnyDot (begin, end, &dot); 258 259 D->firstSigDigit = p; 260 D->exponent = 0; 261 D->normalizedExponent = 0; 262 263 for (; p != end; ++p) { 264 if (*p == '.') { 265 assert(dot == end && "String contains multiple dots"); 266 dot = p++; 267 if (p == end) 268 break; 269 } 270 if (decDigitValue(*p) >= 10U) 271 break; 272 } 273 274 if (p != end) { 275 assert((*p == 'e' || *p == 'E') && "Invalid character in significand"); 276 assert(p != begin && "Significand has no digits"); 277 assert((dot == end || p - begin != 1) && "Significand has no digits"); 278 279 /* p points to the first non-digit in the string */ 280 D->exponent = readExponent(p + 1, end); 281 282 /* Implied decimal point? */ 283 if (dot == end) 284 dot = p; 285 } 286 287 /* If number is all zeroes accept any exponent. */ 288 if (p != D->firstSigDigit) { 289 /* Drop insignificant trailing zeroes. */ 290 if (p != begin) { 291 do 292 do 293 p--; 294 while (p != begin && *p == '0'); 295 while (p != begin && *p == '.'); 296 } 297 298 /* Adjust the exponents for any decimal point. */ 299 D->exponent += static_cast<APFloat::ExponentType>((dot - p) - (dot > p)); 300 D->normalizedExponent = (D->exponent + 301 static_cast<APFloat::ExponentType>((p - D->firstSigDigit) 302 - (dot > D->firstSigDigit && dot < p))); 303 } 304 305 D->lastSigDigit = p; 306 } 307 308 /* Return the trailing fraction of a hexadecimal number. 309 DIGITVALUE is the first hex digit of the fraction, P points to 310 the next digit. */ 311 static lostFraction 312 trailingHexadecimalFraction(StringRef::iterator p, StringRef::iterator end, 313 unsigned int digitValue) 314 { 315 unsigned int hexDigit; 316 317 /* If the first trailing digit isn't 0 or 8 we can work out the 318 fraction immediately. */ 319 if (digitValue > 8) 320 return lfMoreThanHalf; 321 else if (digitValue < 8 && digitValue > 0) 322 return lfLessThanHalf; 323 324 // Otherwise we need to find the first non-zero digit. 325 while (p != end && (*p == '0' || *p == '.')) 326 p++; 327 328 assert(p != end && "Invalid trailing hexadecimal fraction!"); 329 330 hexDigit = hexDigitValue(*p); 331 332 /* If we ran off the end it is exactly zero or one-half, otherwise 333 a little more. */ 334 if (hexDigit == -1U) 335 return digitValue == 0 ? lfExactlyZero: lfExactlyHalf; 336 else 337 return digitValue == 0 ? lfLessThanHalf: lfMoreThanHalf; 338 } 339 340 /* Return the fraction lost were a bignum truncated losing the least 341 significant BITS bits. */ 342 static lostFraction 343 lostFractionThroughTruncation(const integerPart *parts, 344 unsigned int partCount, 345 unsigned int bits) 346 { 347 unsigned int lsb; 348 349 lsb = APInt::tcLSB(parts, partCount); 350 351 /* Note this is guaranteed true if bits == 0, or LSB == -1U. */ 352 if (bits <= lsb) 353 return lfExactlyZero; 354 if (bits == lsb + 1) 355 return lfExactlyHalf; 356 if (bits <= partCount * integerPartWidth && 357 APInt::tcExtractBit(parts, bits - 1)) 358 return lfMoreThanHalf; 359 360 return lfLessThanHalf; 361 } 362 363 /* Shift DST right BITS bits noting lost fraction. */ 364 static lostFraction 365 shiftRight(integerPart *dst, unsigned int parts, unsigned int bits) 366 { 367 lostFraction lost_fraction; 368 369 lost_fraction = lostFractionThroughTruncation(dst, parts, bits); 370 371 APInt::tcShiftRight(dst, parts, bits); 372 373 return lost_fraction; 374 } 375 376 /* Combine the effect of two lost fractions. */ 377 static lostFraction 378 combineLostFractions(lostFraction moreSignificant, 379 lostFraction lessSignificant) 380 { 381 if (lessSignificant != lfExactlyZero) { 382 if (moreSignificant == lfExactlyZero) 383 moreSignificant = lfLessThanHalf; 384 else if (moreSignificant == lfExactlyHalf) 385 moreSignificant = lfMoreThanHalf; 386 } 387 388 return moreSignificant; 389 } 390 391 /* The error from the true value, in half-ulps, on multiplying two 392 floating point numbers, which differ from the value they 393 approximate by at most HUE1 and HUE2 half-ulps, is strictly less 394 than the returned value. 395 396 See "How to Read Floating Point Numbers Accurately" by William D 397 Clinger. */ 398 static unsigned int 399 HUerrBound(bool inexactMultiply, unsigned int HUerr1, unsigned int HUerr2) 400 { 401 assert(HUerr1 < 2 || HUerr2 < 2 || (HUerr1 + HUerr2 < 8)); 402 403 if (HUerr1 + HUerr2 == 0) 404 return inexactMultiply * 2; /* <= inexactMultiply half-ulps. */ 405 else 406 return inexactMultiply + 2 * (HUerr1 + HUerr2); 407 } 408 409 /* The number of ulps from the boundary (zero, or half if ISNEAREST) 410 when the least significant BITS are truncated. BITS cannot be 411 zero. */ 412 static integerPart 413 ulpsFromBoundary(const integerPart *parts, unsigned int bits, bool isNearest) 414 { 415 unsigned int count, partBits; 416 integerPart part, boundary; 417 418 assert(bits != 0); 419 420 bits--; 421 count = bits / integerPartWidth; 422 partBits = bits % integerPartWidth + 1; 423 424 part = parts[count] & (~(integerPart) 0 >> (integerPartWidth - partBits)); 425 426 if (isNearest) 427 boundary = (integerPart) 1 << (partBits - 1); 428 else 429 boundary = 0; 430 431 if (count == 0) { 432 if (part - boundary <= boundary - part) 433 return part - boundary; 434 else 435 return boundary - part; 436 } 437 438 if (part == boundary) { 439 while (--count) 440 if (parts[count]) 441 return ~(integerPart) 0; /* A lot. */ 442 443 return parts[0]; 444 } else if (part == boundary - 1) { 445 while (--count) 446 if (~parts[count]) 447 return ~(integerPart) 0; /* A lot. */ 448 449 return -parts[0]; 450 } 451 452 return ~(integerPart) 0; /* A lot. */ 453 } 454 455 /* Place pow(5, power) in DST, and return the number of parts used. 456 DST must be at least one part larger than size of the answer. */ 457 static unsigned int 458 powerOf5(integerPart *dst, unsigned int power) 459 { 460 static const integerPart firstEightPowers[] = { 1, 5, 25, 125, 625, 3125, 461 15625, 78125 }; 462 integerPart pow5s[maxPowerOfFiveParts * 2 + 5]; 463 pow5s[0] = 78125 * 5; 464 465 unsigned int partsCount[16] = { 1 }; 466 integerPart scratch[maxPowerOfFiveParts], *p1, *p2, *pow5; 467 unsigned int result; 468 assert(power <= maxExponent); 469 470 p1 = dst; 471 p2 = scratch; 472 473 *p1 = firstEightPowers[power & 7]; 474 power >>= 3; 475 476 result = 1; 477 pow5 = pow5s; 478 479 for (unsigned int n = 0; power; power >>= 1, n++) { 480 unsigned int pc; 481 482 pc = partsCount[n]; 483 484 /* Calculate pow(5,pow(2,n+3)) if we haven't yet. */ 485 if (pc == 0) { 486 pc = partsCount[n - 1]; 487 APInt::tcFullMultiply(pow5, pow5 - pc, pow5 - pc, pc, pc); 488 pc *= 2; 489 if (pow5[pc - 1] == 0) 490 pc--; 491 partsCount[n] = pc; 492 } 493 494 if (power & 1) { 495 integerPart *tmp; 496 497 APInt::tcFullMultiply(p2, p1, pow5, result, pc); 498 result += pc; 499 if (p2[result - 1] == 0) 500 result--; 501 502 /* Now result is in p1 with partsCount parts and p2 is scratch 503 space. */ 504 tmp = p1; 505 p1 = p2; 506 p2 = tmp; 507 } 508 509 pow5 += pc; 510 } 511 512 if (p1 != dst) 513 APInt::tcAssign(dst, p1, result); 514 515 return result; 516 } 517 518 /* Zero at the end to avoid modular arithmetic when adding one; used 519 when rounding up during hexadecimal output. */ 520 static const char hexDigitsLower[] = "0123456789abcdef0"; 521 static const char hexDigitsUpper[] = "0123456789ABCDEF0"; 522 static const char infinityL[] = "infinity"; 523 static const char infinityU[] = "INFINITY"; 524 static const char NaNL[] = "nan"; 525 static const char NaNU[] = "NAN"; 526 527 /* Write out an integerPart in hexadecimal, starting with the most 528 significant nibble. Write out exactly COUNT hexdigits, return 529 COUNT. */ 530 static unsigned int 531 partAsHex (char *dst, integerPart part, unsigned int count, 532 const char *hexDigitChars) 533 { 534 unsigned int result = count; 535 536 assert(count != 0 && count <= integerPartWidth / 4); 537 538 part >>= (integerPartWidth - 4 * count); 539 while (count--) { 540 dst[count] = hexDigitChars[part & 0xf]; 541 part >>= 4; 542 } 543 544 return result; 545 } 546 547 /* Write out an unsigned decimal integer. */ 548 static char * 549 writeUnsignedDecimal (char *dst, unsigned int n) 550 { 551 char buff[40], *p; 552 553 p = buff; 554 do 555 *p++ = '0' + n % 10; 556 while (n /= 10); 557 558 do 559 *dst++ = *--p; 560 while (p != buff); 561 562 return dst; 563 } 564 565 /* Write out a signed decimal integer. */ 566 static char * 567 writeSignedDecimal (char *dst, int value) 568 { 569 if (value < 0) { 570 *dst++ = '-'; 571 dst = writeUnsignedDecimal(dst, -(unsigned) value); 572 } else 573 dst = writeUnsignedDecimal(dst, value); 574 575 return dst; 576 } 577 578 /* Constructors. */ 579 void 580 APFloat::initialize(const fltSemantics *ourSemantics) 581 { 582 unsigned int count; 583 584 semantics = ourSemantics; 585 count = partCount(); 586 if (count > 1) 587 significand.parts = new integerPart[count]; 588 } 589 590 void 591 APFloat::freeSignificand() 592 { 593 if (needsCleanup()) 594 delete [] significand.parts; 595 } 596 597 void 598 APFloat::assign(const APFloat &rhs) 599 { 600 assert(semantics == rhs.semantics); 601 602 sign = rhs.sign; 603 category = rhs.category; 604 exponent = rhs.exponent; 605 if (isFiniteNonZero() || category == fcNaN) 606 copySignificand(rhs); 607 } 608 609 void 610 APFloat::copySignificand(const APFloat &rhs) 611 { 612 assert(isFiniteNonZero() || category == fcNaN); 613 assert(rhs.partCount() >= partCount()); 614 615 APInt::tcAssign(significandParts(), rhs.significandParts(), 616 partCount()); 617 } 618 619 /* Make this number a NaN, with an arbitrary but deterministic value 620 for the significand. If double or longer, this is a signalling NaN, 621 which may not be ideal. If float, this is QNaN(0). */ 622 void APFloat::makeNaN(bool SNaN, bool Negative, const APInt *fill) 623 { 624 category = fcNaN; 625 sign = Negative; 626 627 integerPart *significand = significandParts(); 628 unsigned numParts = partCount(); 629 630 // Set the significand bits to the fill. 631 if (!fill || fill->getNumWords() < numParts) 632 APInt::tcSet(significand, 0, numParts); 633 if (fill) { 634 APInt::tcAssign(significand, fill->getRawData(), 635 std::min(fill->getNumWords(), numParts)); 636 637 // Zero out the excess bits of the significand. 638 unsigned bitsToPreserve = semantics->precision - 1; 639 unsigned part = bitsToPreserve / 64; 640 bitsToPreserve %= 64; 641 significand[part] &= ((1ULL << bitsToPreserve) - 1); 642 for (part++; part != numParts; ++part) 643 significand[part] = 0; 644 } 645 646 unsigned QNaNBit = semantics->precision - 2; 647 648 if (SNaN) { 649 // We always have to clear the QNaN bit to make it an SNaN. 650 APInt::tcClearBit(significand, QNaNBit); 651 652 // If there are no bits set in the payload, we have to set 653 // *something* to make it a NaN instead of an infinity; 654 // conventionally, this is the next bit down from the QNaN bit. 655 if (APInt::tcIsZero(significand, numParts)) 656 APInt::tcSetBit(significand, QNaNBit - 1); 657 } else { 658 // We always have to set the QNaN bit to make it a QNaN. 659 APInt::tcSetBit(significand, QNaNBit); 660 } 661 662 // For x87 extended precision, we want to make a NaN, not a 663 // pseudo-NaN. Maybe we should expose the ability to make 664 // pseudo-NaNs? 665 if (semantics == &APFloat::x87DoubleExtended) 666 APInt::tcSetBit(significand, QNaNBit + 1); 667 } 668 669 APFloat APFloat::makeNaN(const fltSemantics &Sem, bool SNaN, bool Negative, 670 const APInt *fill) { 671 APFloat value(Sem, uninitialized); 672 value.makeNaN(SNaN, Negative, fill); 673 return value; 674 } 675 676 APFloat & 677 APFloat::operator=(const APFloat &rhs) 678 { 679 if (this != &rhs) { 680 if (semantics != rhs.semantics) { 681 freeSignificand(); 682 initialize(rhs.semantics); 683 } 684 assign(rhs); 685 } 686 687 return *this; 688 } 689 690 APFloat & 691 APFloat::operator=(APFloat &&rhs) { 692 freeSignificand(); 693 694 semantics = rhs.semantics; 695 significand = rhs.significand; 696 exponent = rhs.exponent; 697 category = rhs.category; 698 sign = rhs.sign; 699 700 rhs.semantics = &Bogus; 701 return *this; 702 } 703 704 bool 705 APFloat::isDenormal() const { 706 return isFiniteNonZero() && (exponent == semantics->minExponent) && 707 (APInt::tcExtractBit(significandParts(), 708 semantics->precision - 1) == 0); 709 } 710 711 bool 712 APFloat::isSmallest() const { 713 // The smallest number by magnitude in our format will be the smallest 714 // denormal, i.e. the floating point number with exponent being minimum 715 // exponent and significand bitwise equal to 1 (i.e. with MSB equal to 0). 716 return isFiniteNonZero() && exponent == semantics->minExponent && 717 significandMSB() == 0; 718 } 719 720 bool APFloat::isSignificandAllOnes() const { 721 // Test if the significand excluding the integral bit is all ones. This allows 722 // us to test for binade boundaries. 723 const integerPart *Parts = significandParts(); 724 const unsigned PartCount = partCount(); 725 for (unsigned i = 0; i < PartCount - 1; i++) 726 if (~Parts[i]) 727 return false; 728 729 // Set the unused high bits to all ones when we compare. 730 const unsigned NumHighBits = 731 PartCount*integerPartWidth - semantics->precision + 1; 732 assert(NumHighBits <= integerPartWidth && "Can not have more high bits to " 733 "fill than integerPartWidth"); 734 const integerPart HighBitFill = 735 ~integerPart(0) << (integerPartWidth - NumHighBits); 736 if (~(Parts[PartCount - 1] | HighBitFill)) 737 return false; 738 739 return true; 740 } 741 742 bool APFloat::isSignificandAllZeros() const { 743 // Test if the significand excluding the integral bit is all zeros. This 744 // allows us to test for binade boundaries. 745 const integerPart *Parts = significandParts(); 746 const unsigned PartCount = partCount(); 747 748 for (unsigned i = 0; i < PartCount - 1; i++) 749 if (Parts[i]) 750 return false; 751 752 const unsigned NumHighBits = 753 PartCount*integerPartWidth - semantics->precision + 1; 754 assert(NumHighBits <= integerPartWidth && "Can not have more high bits to " 755 "clear than integerPartWidth"); 756 const integerPart HighBitMask = ~integerPart(0) >> NumHighBits; 757 758 if (Parts[PartCount - 1] & HighBitMask) 759 return false; 760 761 return true; 762 } 763 764 bool 765 APFloat::isLargest() const { 766 // The largest number by magnitude in our format will be the floating point 767 // number with maximum exponent and with significand that is all ones. 768 return isFiniteNonZero() && exponent == semantics->maxExponent 769 && isSignificandAllOnes(); 770 } 771 772 bool 773 APFloat::isInteger() const { 774 // This could be made more efficient; I'm going for obviously correct. 775 if (!isFinite()) return false; 776 APFloat truncated = *this; 777 truncated.roundToIntegral(rmTowardZero); 778 return compare(truncated) == cmpEqual; 779 } 780 781 bool 782 APFloat::bitwiseIsEqual(const APFloat &rhs) const { 783 if (this == &rhs) 784 return true; 785 if (semantics != rhs.semantics || 786 category != rhs.category || 787 sign != rhs.sign) 788 return false; 789 if (category==fcZero || category==fcInfinity) 790 return true; 791 792 if (isFiniteNonZero() && exponent != rhs.exponent) 793 return false; 794 795 return std::equal(significandParts(), significandParts() + partCount(), 796 rhs.significandParts()); 797 } 798 799 APFloat::APFloat(const fltSemantics &ourSemantics, integerPart value) { 800 initialize(&ourSemantics); 801 sign = 0; 802 category = fcNormal; 803 zeroSignificand(); 804 exponent = ourSemantics.precision - 1; 805 significandParts()[0] = value; 806 normalize(rmNearestTiesToEven, lfExactlyZero); 807 } 808 809 APFloat::APFloat(const fltSemantics &ourSemantics) { 810 initialize(&ourSemantics); 811 category = fcZero; 812 sign = false; 813 } 814 815 APFloat::APFloat(const fltSemantics &ourSemantics, uninitializedTag tag) { 816 // Allocates storage if necessary but does not initialize it. 817 initialize(&ourSemantics); 818 } 819 820 APFloat::APFloat(const fltSemantics &ourSemantics, StringRef text) { 821 initialize(&ourSemantics); 822 convertFromString(text, rmNearestTiesToEven); 823 } 824 825 APFloat::APFloat(const APFloat &rhs) { 826 initialize(rhs.semantics); 827 assign(rhs); 828 } 829 830 APFloat::APFloat(APFloat &&rhs) : semantics(&Bogus) { 831 *this = std::move(rhs); 832 } 833 834 APFloat::~APFloat() 835 { 836 freeSignificand(); 837 } 838 839 // Profile - This method 'profiles' an APFloat for use with FoldingSet. 840 void APFloat::Profile(FoldingSetNodeID& ID) const { 841 ID.Add(bitcastToAPInt()); 842 } 843 844 unsigned int 845 APFloat::partCount() const 846 { 847 return partCountForBits(semantics->precision + 1); 848 } 849 850 unsigned int 851 APFloat::semanticsPrecision(const fltSemantics &semantics) 852 { 853 return semantics.precision; 854 } 855 APFloat::ExponentType 856 APFloat::semanticsMaxExponent(const fltSemantics &semantics) 857 { 858 return semantics.maxExponent; 859 } 860 APFloat::ExponentType 861 APFloat::semanticsMinExponent(const fltSemantics &semantics) 862 { 863 return semantics.minExponent; 864 } 865 unsigned int 866 APFloat::semanticsSizeInBits(const fltSemantics &semantics) 867 { 868 return semantics.sizeInBits; 869 } 870 871 const integerPart * 872 APFloat::significandParts() const 873 { 874 return const_cast<APFloat *>(this)->significandParts(); 875 } 876 877 integerPart * 878 APFloat::significandParts() 879 { 880 if (partCount() > 1) 881 return significand.parts; 882 else 883 return &significand.part; 884 } 885 886 void 887 APFloat::zeroSignificand() 888 { 889 APInt::tcSet(significandParts(), 0, partCount()); 890 } 891 892 /* Increment an fcNormal floating point number's significand. */ 893 void 894 APFloat::incrementSignificand() 895 { 896 integerPart carry; 897 898 carry = APInt::tcIncrement(significandParts(), partCount()); 899 900 /* Our callers should never cause us to overflow. */ 901 assert(carry == 0); 902 (void)carry; 903 } 904 905 /* Add the significand of the RHS. Returns the carry flag. */ 906 integerPart 907 APFloat::addSignificand(const APFloat &rhs) 908 { 909 integerPart *parts; 910 911 parts = significandParts(); 912 913 assert(semantics == rhs.semantics); 914 assert(exponent == rhs.exponent); 915 916 return APInt::tcAdd(parts, rhs.significandParts(), 0, partCount()); 917 } 918 919 /* Subtract the significand of the RHS with a borrow flag. Returns 920 the borrow flag. */ 921 integerPart 922 APFloat::subtractSignificand(const APFloat &rhs, integerPart borrow) 923 { 924 integerPart *parts; 925 926 parts = significandParts(); 927 928 assert(semantics == rhs.semantics); 929 assert(exponent == rhs.exponent); 930 931 return APInt::tcSubtract(parts, rhs.significandParts(), borrow, 932 partCount()); 933 } 934 935 /* Multiply the significand of the RHS. If ADDEND is non-NULL, add it 936 on to the full-precision result of the multiplication. Returns the 937 lost fraction. */ 938 lostFraction 939 APFloat::multiplySignificand(const APFloat &rhs, const APFloat *addend) 940 { 941 unsigned int omsb; // One, not zero, based MSB. 942 unsigned int partsCount, newPartsCount, precision; 943 integerPart *lhsSignificand; 944 integerPart scratch[4]; 945 integerPart *fullSignificand; 946 lostFraction lost_fraction; 947 bool ignored; 948 949 assert(semantics == rhs.semantics); 950 951 precision = semantics->precision; 952 953 // Allocate space for twice as many bits as the original significand, plus one 954 // extra bit for the addition to overflow into. 955 newPartsCount = partCountForBits(precision * 2 + 1); 956 957 if (newPartsCount > 4) 958 fullSignificand = new integerPart[newPartsCount]; 959 else 960 fullSignificand = scratch; 961 962 lhsSignificand = significandParts(); 963 partsCount = partCount(); 964 965 APInt::tcFullMultiply(fullSignificand, lhsSignificand, 966 rhs.significandParts(), partsCount, partsCount); 967 968 lost_fraction = lfExactlyZero; 969 omsb = APInt::tcMSB(fullSignificand, newPartsCount) + 1; 970 exponent += rhs.exponent; 971 972 // Assume the operands involved in the multiplication are single-precision 973 // FP, and the two multiplicants are: 974 // *this = a23 . a22 ... a0 * 2^e1 975 // rhs = b23 . b22 ... b0 * 2^e2 976 // the result of multiplication is: 977 // *this = c48 c47 c46 . c45 ... c0 * 2^(e1+e2) 978 // Note that there are three significant bits at the left-hand side of the 979 // radix point: two for the multiplication, and an overflow bit for the 980 // addition (that will always be zero at this point). Move the radix point 981 // toward left by two bits, and adjust exponent accordingly. 982 exponent += 2; 983 984 if (addend && addend->isNonZero()) { 985 // The intermediate result of the multiplication has "2 * precision" 986 // signicant bit; adjust the addend to be consistent with mul result. 987 // 988 Significand savedSignificand = significand; 989 const fltSemantics *savedSemantics = semantics; 990 fltSemantics extendedSemantics; 991 opStatus status; 992 unsigned int extendedPrecision; 993 994 // Normalize our MSB to one below the top bit to allow for overflow. 995 extendedPrecision = 2 * precision + 1; 996 if (omsb != extendedPrecision - 1) { 997 assert(extendedPrecision > omsb); 998 APInt::tcShiftLeft(fullSignificand, newPartsCount, 999 (extendedPrecision - 1) - omsb); 1000 exponent -= (extendedPrecision - 1) - omsb; 1001 } 1002 1003 /* Create new semantics. */ 1004 extendedSemantics = *semantics; 1005 extendedSemantics.precision = extendedPrecision; 1006 1007 if (newPartsCount == 1) 1008 significand.part = fullSignificand[0]; 1009 else 1010 significand.parts = fullSignificand; 1011 semantics = &extendedSemantics; 1012 1013 APFloat extendedAddend(*addend); 1014 status = extendedAddend.convert(extendedSemantics, rmTowardZero, &ignored); 1015 assert(status == opOK); 1016 (void)status; 1017 1018 // Shift the significand of the addend right by one bit. This guarantees 1019 // that the high bit of the significand is zero (same as fullSignificand), 1020 // so the addition will overflow (if it does overflow at all) into the top bit. 1021 lost_fraction = extendedAddend.shiftSignificandRight(1); 1022 assert(lost_fraction == lfExactlyZero && 1023 "Lost precision while shifting addend for fused-multiply-add."); 1024 1025 lost_fraction = addOrSubtractSignificand(extendedAddend, false); 1026 1027 /* Restore our state. */ 1028 if (newPartsCount == 1) 1029 fullSignificand[0] = significand.part; 1030 significand = savedSignificand; 1031 semantics = savedSemantics; 1032 1033 omsb = APInt::tcMSB(fullSignificand, newPartsCount) + 1; 1034 } 1035 1036 // Convert the result having "2 * precision" significant-bits back to the one 1037 // having "precision" significant-bits. First, move the radix point from 1038 // poision "2*precision - 1" to "precision - 1". The exponent need to be 1039 // adjusted by "2*precision - 1" - "precision - 1" = "precision". 1040 exponent -= precision + 1; 1041 1042 // In case MSB resides at the left-hand side of radix point, shift the 1043 // mantissa right by some amount to make sure the MSB reside right before 1044 // the radix point (i.e. "MSB . rest-significant-bits"). 1045 // 1046 // Note that the result is not normalized when "omsb < precision". So, the 1047 // caller needs to call APFloat::normalize() if normalized value is expected. 1048 if (omsb > precision) { 1049 unsigned int bits, significantParts; 1050 lostFraction lf; 1051 1052 bits = omsb - precision; 1053 significantParts = partCountForBits(omsb); 1054 lf = shiftRight(fullSignificand, significantParts, bits); 1055 lost_fraction = combineLostFractions(lf, lost_fraction); 1056 exponent += bits; 1057 } 1058 1059 APInt::tcAssign(lhsSignificand, fullSignificand, partsCount); 1060 1061 if (newPartsCount > 4) 1062 delete [] fullSignificand; 1063 1064 return lost_fraction; 1065 } 1066 1067 /* Multiply the significands of LHS and RHS to DST. */ 1068 lostFraction 1069 APFloat::divideSignificand(const APFloat &rhs) 1070 { 1071 unsigned int bit, i, partsCount; 1072 const integerPart *rhsSignificand; 1073 integerPart *lhsSignificand, *dividend, *divisor; 1074 integerPart scratch[4]; 1075 lostFraction lost_fraction; 1076 1077 assert(semantics == rhs.semantics); 1078 1079 lhsSignificand = significandParts(); 1080 rhsSignificand = rhs.significandParts(); 1081 partsCount = partCount(); 1082 1083 if (partsCount > 2) 1084 dividend = new integerPart[partsCount * 2]; 1085 else 1086 dividend = scratch; 1087 1088 divisor = dividend + partsCount; 1089 1090 /* Copy the dividend and divisor as they will be modified in-place. */ 1091 for (i = 0; i < partsCount; i++) { 1092 dividend[i] = lhsSignificand[i]; 1093 divisor[i] = rhsSignificand[i]; 1094 lhsSignificand[i] = 0; 1095 } 1096 1097 exponent -= rhs.exponent; 1098 1099 unsigned int precision = semantics->precision; 1100 1101 /* Normalize the divisor. */ 1102 bit = precision - APInt::tcMSB(divisor, partsCount) - 1; 1103 if (bit) { 1104 exponent += bit; 1105 APInt::tcShiftLeft(divisor, partsCount, bit); 1106 } 1107 1108 /* Normalize the dividend. */ 1109 bit = precision - APInt::tcMSB(dividend, partsCount) - 1; 1110 if (bit) { 1111 exponent -= bit; 1112 APInt::tcShiftLeft(dividend, partsCount, bit); 1113 } 1114 1115 /* Ensure the dividend >= divisor initially for the loop below. 1116 Incidentally, this means that the division loop below is 1117 guaranteed to set the integer bit to one. */ 1118 if (APInt::tcCompare(dividend, divisor, partsCount) < 0) { 1119 exponent--; 1120 APInt::tcShiftLeft(dividend, partsCount, 1); 1121 assert(APInt::tcCompare(dividend, divisor, partsCount) >= 0); 1122 } 1123 1124 /* Long division. */ 1125 for (bit = precision; bit; bit -= 1) { 1126 if (APInt::tcCompare(dividend, divisor, partsCount) >= 0) { 1127 APInt::tcSubtract(dividend, divisor, 0, partsCount); 1128 APInt::tcSetBit(lhsSignificand, bit - 1); 1129 } 1130 1131 APInt::tcShiftLeft(dividend, partsCount, 1); 1132 } 1133 1134 /* Figure out the lost fraction. */ 1135 int cmp = APInt::tcCompare(dividend, divisor, partsCount); 1136 1137 if (cmp > 0) 1138 lost_fraction = lfMoreThanHalf; 1139 else if (cmp == 0) 1140 lost_fraction = lfExactlyHalf; 1141 else if (APInt::tcIsZero(dividend, partsCount)) 1142 lost_fraction = lfExactlyZero; 1143 else 1144 lost_fraction = lfLessThanHalf; 1145 1146 if (partsCount > 2) 1147 delete [] dividend; 1148 1149 return lost_fraction; 1150 } 1151 1152 unsigned int 1153 APFloat::significandMSB() const 1154 { 1155 return APInt::tcMSB(significandParts(), partCount()); 1156 } 1157 1158 unsigned int 1159 APFloat::significandLSB() const 1160 { 1161 return APInt::tcLSB(significandParts(), partCount()); 1162 } 1163 1164 /* Note that a zero result is NOT normalized to fcZero. */ 1165 lostFraction 1166 APFloat::shiftSignificandRight(unsigned int bits) 1167 { 1168 /* Our exponent should not overflow. */ 1169 assert((ExponentType) (exponent + bits) >= exponent); 1170 1171 exponent += bits; 1172 1173 return shiftRight(significandParts(), partCount(), bits); 1174 } 1175 1176 /* Shift the significand left BITS bits, subtract BITS from its exponent. */ 1177 void 1178 APFloat::shiftSignificandLeft(unsigned int bits) 1179 { 1180 assert(bits < semantics->precision); 1181 1182 if (bits) { 1183 unsigned int partsCount = partCount(); 1184 1185 APInt::tcShiftLeft(significandParts(), partsCount, bits); 1186 exponent -= bits; 1187 1188 assert(!APInt::tcIsZero(significandParts(), partsCount)); 1189 } 1190 } 1191 1192 APFloat::cmpResult 1193 APFloat::compareAbsoluteValue(const APFloat &rhs) const 1194 { 1195 int compare; 1196 1197 assert(semantics == rhs.semantics); 1198 assert(isFiniteNonZero()); 1199 assert(rhs.isFiniteNonZero()); 1200 1201 compare = exponent - rhs.exponent; 1202 1203 /* If exponents are equal, do an unsigned bignum comparison of the 1204 significands. */ 1205 if (compare == 0) 1206 compare = APInt::tcCompare(significandParts(), rhs.significandParts(), 1207 partCount()); 1208 1209 if (compare > 0) 1210 return cmpGreaterThan; 1211 else if (compare < 0) 1212 return cmpLessThan; 1213 else 1214 return cmpEqual; 1215 } 1216 1217 /* Handle overflow. Sign is preserved. We either become infinity or 1218 the largest finite number. */ 1219 APFloat::opStatus 1220 APFloat::handleOverflow(roundingMode rounding_mode) 1221 { 1222 /* Infinity? */ 1223 if (rounding_mode == rmNearestTiesToEven || 1224 rounding_mode == rmNearestTiesToAway || 1225 (rounding_mode == rmTowardPositive && !sign) || 1226 (rounding_mode == rmTowardNegative && sign)) { 1227 category = fcInfinity; 1228 return (opStatus) (opOverflow | opInexact); 1229 } 1230 1231 /* Otherwise we become the largest finite number. */ 1232 category = fcNormal; 1233 exponent = semantics->maxExponent; 1234 APInt::tcSetLeastSignificantBits(significandParts(), partCount(), 1235 semantics->precision); 1236 1237 return opInexact; 1238 } 1239 1240 /* Returns TRUE if, when truncating the current number, with BIT the 1241 new LSB, with the given lost fraction and rounding mode, the result 1242 would need to be rounded away from zero (i.e., by increasing the 1243 signficand). This routine must work for fcZero of both signs, and 1244 fcNormal numbers. */ 1245 bool 1246 APFloat::roundAwayFromZero(roundingMode rounding_mode, 1247 lostFraction lost_fraction, 1248 unsigned int bit) const 1249 { 1250 /* NaNs and infinities should not have lost fractions. */ 1251 assert(isFiniteNonZero() || category == fcZero); 1252 1253 /* Current callers never pass this so we don't handle it. */ 1254 assert(lost_fraction != lfExactlyZero); 1255 1256 switch (rounding_mode) { 1257 case rmNearestTiesToAway: 1258 return lost_fraction == lfExactlyHalf || lost_fraction == lfMoreThanHalf; 1259 1260 case rmNearestTiesToEven: 1261 if (lost_fraction == lfMoreThanHalf) 1262 return true; 1263 1264 /* Our zeroes don't have a significand to test. */ 1265 if (lost_fraction == lfExactlyHalf && category != fcZero) 1266 return APInt::tcExtractBit(significandParts(), bit); 1267 1268 return false; 1269 1270 case rmTowardZero: 1271 return false; 1272 1273 case rmTowardPositive: 1274 return !sign; 1275 1276 case rmTowardNegative: 1277 return sign; 1278 } 1279 llvm_unreachable("Invalid rounding mode found"); 1280 } 1281 1282 APFloat::opStatus 1283 APFloat::normalize(roundingMode rounding_mode, 1284 lostFraction lost_fraction) 1285 { 1286 unsigned int omsb; /* One, not zero, based MSB. */ 1287 int exponentChange; 1288 1289 if (!isFiniteNonZero()) 1290 return opOK; 1291 1292 /* Before rounding normalize the exponent of fcNormal numbers. */ 1293 omsb = significandMSB() + 1; 1294 1295 if (omsb) { 1296 /* OMSB is numbered from 1. We want to place it in the integer 1297 bit numbered PRECISION if possible, with a compensating change in 1298 the exponent. */ 1299 exponentChange = omsb - semantics->precision; 1300 1301 /* If the resulting exponent is too high, overflow according to 1302 the rounding mode. */ 1303 if (exponent + exponentChange > semantics->maxExponent) 1304 return handleOverflow(rounding_mode); 1305 1306 /* Subnormal numbers have exponent minExponent, and their MSB 1307 is forced based on that. */ 1308 if (exponent + exponentChange < semantics->minExponent) 1309 exponentChange = semantics->minExponent - exponent; 1310 1311 /* Shifting left is easy as we don't lose precision. */ 1312 if (exponentChange < 0) { 1313 assert(lost_fraction == lfExactlyZero); 1314 1315 shiftSignificandLeft(-exponentChange); 1316 1317 return opOK; 1318 } 1319 1320 if (exponentChange > 0) { 1321 lostFraction lf; 1322 1323 /* Shift right and capture any new lost fraction. */ 1324 lf = shiftSignificandRight(exponentChange); 1325 1326 lost_fraction = combineLostFractions(lf, lost_fraction); 1327 1328 /* Keep OMSB up-to-date. */ 1329 if (omsb > (unsigned) exponentChange) 1330 omsb -= exponentChange; 1331 else 1332 omsb = 0; 1333 } 1334 } 1335 1336 /* Now round the number according to rounding_mode given the lost 1337 fraction. */ 1338 1339 /* As specified in IEEE 754, since we do not trap we do not report 1340 underflow for exact results. */ 1341 if (lost_fraction == lfExactlyZero) { 1342 /* Canonicalize zeroes. */ 1343 if (omsb == 0) 1344 category = fcZero; 1345 1346 return opOK; 1347 } 1348 1349 /* Increment the significand if we're rounding away from zero. */ 1350 if (roundAwayFromZero(rounding_mode, lost_fraction, 0)) { 1351 if (omsb == 0) 1352 exponent = semantics->minExponent; 1353 1354 incrementSignificand(); 1355 omsb = significandMSB() + 1; 1356 1357 /* Did the significand increment overflow? */ 1358 if (omsb == (unsigned) semantics->precision + 1) { 1359 /* Renormalize by incrementing the exponent and shifting our 1360 significand right one. However if we already have the 1361 maximum exponent we overflow to infinity. */ 1362 if (exponent == semantics->maxExponent) { 1363 category = fcInfinity; 1364 1365 return (opStatus) (opOverflow | opInexact); 1366 } 1367 1368 shiftSignificandRight(1); 1369 1370 return opInexact; 1371 } 1372 } 1373 1374 /* The normal case - we were and are not denormal, and any 1375 significand increment above didn't overflow. */ 1376 if (omsb == semantics->precision) 1377 return opInexact; 1378 1379 /* We have a non-zero denormal. */ 1380 assert(omsb < semantics->precision); 1381 1382 /* Canonicalize zeroes. */ 1383 if (omsb == 0) 1384 category = fcZero; 1385 1386 /* The fcZero case is a denormal that underflowed to zero. */ 1387 return (opStatus) (opUnderflow | opInexact); 1388 } 1389 1390 APFloat::opStatus 1391 APFloat::addOrSubtractSpecials(const APFloat &rhs, bool subtract) 1392 { 1393 switch (PackCategoriesIntoKey(category, rhs.category)) { 1394 default: 1395 llvm_unreachable(nullptr); 1396 1397 case PackCategoriesIntoKey(fcNaN, fcZero): 1398 case PackCategoriesIntoKey(fcNaN, fcNormal): 1399 case PackCategoriesIntoKey(fcNaN, fcInfinity): 1400 case PackCategoriesIntoKey(fcNaN, fcNaN): 1401 case PackCategoriesIntoKey(fcNormal, fcZero): 1402 case PackCategoriesIntoKey(fcInfinity, fcNormal): 1403 case PackCategoriesIntoKey(fcInfinity, fcZero): 1404 return opOK; 1405 1406 case PackCategoriesIntoKey(fcZero, fcNaN): 1407 case PackCategoriesIntoKey(fcNormal, fcNaN): 1408 case PackCategoriesIntoKey(fcInfinity, fcNaN): 1409 // We need to be sure to flip the sign here for subtraction because we 1410 // don't have a separate negate operation so -NaN becomes 0 - NaN here. 1411 sign = rhs.sign ^ subtract; 1412 category = fcNaN; 1413 copySignificand(rhs); 1414 return opOK; 1415 1416 case PackCategoriesIntoKey(fcNormal, fcInfinity): 1417 case PackCategoriesIntoKey(fcZero, fcInfinity): 1418 category = fcInfinity; 1419 sign = rhs.sign ^ subtract; 1420 return opOK; 1421 1422 case PackCategoriesIntoKey(fcZero, fcNormal): 1423 assign(rhs); 1424 sign = rhs.sign ^ subtract; 1425 return opOK; 1426 1427 case PackCategoriesIntoKey(fcZero, fcZero): 1428 /* Sign depends on rounding mode; handled by caller. */ 1429 return opOK; 1430 1431 case PackCategoriesIntoKey(fcInfinity, fcInfinity): 1432 /* Differently signed infinities can only be validly 1433 subtracted. */ 1434 if (((sign ^ rhs.sign)!=0) != subtract) { 1435 makeNaN(); 1436 return opInvalidOp; 1437 } 1438 1439 return opOK; 1440 1441 case PackCategoriesIntoKey(fcNormal, fcNormal): 1442 return opDivByZero; 1443 } 1444 } 1445 1446 /* Add or subtract two normal numbers. */ 1447 lostFraction 1448 APFloat::addOrSubtractSignificand(const APFloat &rhs, bool subtract) 1449 { 1450 integerPart carry; 1451 lostFraction lost_fraction; 1452 int bits; 1453 1454 /* Determine if the operation on the absolute values is effectively 1455 an addition or subtraction. */ 1456 subtract ^= static_cast<bool>(sign ^ rhs.sign); 1457 1458 /* Are we bigger exponent-wise than the RHS? */ 1459 bits = exponent - rhs.exponent; 1460 1461 /* Subtraction is more subtle than one might naively expect. */ 1462 if (subtract) { 1463 APFloat temp_rhs(rhs); 1464 bool reverse; 1465 1466 if (bits == 0) { 1467 reverse = compareAbsoluteValue(temp_rhs) == cmpLessThan; 1468 lost_fraction = lfExactlyZero; 1469 } else if (bits > 0) { 1470 lost_fraction = temp_rhs.shiftSignificandRight(bits - 1); 1471 shiftSignificandLeft(1); 1472 reverse = false; 1473 } else { 1474 lost_fraction = shiftSignificandRight(-bits - 1); 1475 temp_rhs.shiftSignificandLeft(1); 1476 reverse = true; 1477 } 1478 1479 if (reverse) { 1480 carry = temp_rhs.subtractSignificand 1481 (*this, lost_fraction != lfExactlyZero); 1482 copySignificand(temp_rhs); 1483 sign = !sign; 1484 } else { 1485 carry = subtractSignificand 1486 (temp_rhs, lost_fraction != lfExactlyZero); 1487 } 1488 1489 /* Invert the lost fraction - it was on the RHS and 1490 subtracted. */ 1491 if (lost_fraction == lfLessThanHalf) 1492 lost_fraction = lfMoreThanHalf; 1493 else if (lost_fraction == lfMoreThanHalf) 1494 lost_fraction = lfLessThanHalf; 1495 1496 /* The code above is intended to ensure that no borrow is 1497 necessary. */ 1498 assert(!carry); 1499 (void)carry; 1500 } else { 1501 if (bits > 0) { 1502 APFloat temp_rhs(rhs); 1503 1504 lost_fraction = temp_rhs.shiftSignificandRight(bits); 1505 carry = addSignificand(temp_rhs); 1506 } else { 1507 lost_fraction = shiftSignificandRight(-bits); 1508 carry = addSignificand(rhs); 1509 } 1510 1511 /* We have a guard bit; generating a carry cannot happen. */ 1512 assert(!carry); 1513 (void)carry; 1514 } 1515 1516 return lost_fraction; 1517 } 1518 1519 APFloat::opStatus 1520 APFloat::multiplySpecials(const APFloat &rhs) 1521 { 1522 switch (PackCategoriesIntoKey(category, rhs.category)) { 1523 default: 1524 llvm_unreachable(nullptr); 1525 1526 case PackCategoriesIntoKey(fcNaN, fcZero): 1527 case PackCategoriesIntoKey(fcNaN, fcNormal): 1528 case PackCategoriesIntoKey(fcNaN, fcInfinity): 1529 case PackCategoriesIntoKey(fcNaN, fcNaN): 1530 sign = false; 1531 return opOK; 1532 1533 case PackCategoriesIntoKey(fcZero, fcNaN): 1534 case PackCategoriesIntoKey(fcNormal, fcNaN): 1535 case PackCategoriesIntoKey(fcInfinity, fcNaN): 1536 sign = false; 1537 category = fcNaN; 1538 copySignificand(rhs); 1539 return opOK; 1540 1541 case PackCategoriesIntoKey(fcNormal, fcInfinity): 1542 case PackCategoriesIntoKey(fcInfinity, fcNormal): 1543 case PackCategoriesIntoKey(fcInfinity, fcInfinity): 1544 category = fcInfinity; 1545 return opOK; 1546 1547 case PackCategoriesIntoKey(fcZero, fcNormal): 1548 case PackCategoriesIntoKey(fcNormal, fcZero): 1549 case PackCategoriesIntoKey(fcZero, fcZero): 1550 category = fcZero; 1551 return opOK; 1552 1553 case PackCategoriesIntoKey(fcZero, fcInfinity): 1554 case PackCategoriesIntoKey(fcInfinity, fcZero): 1555 makeNaN(); 1556 return opInvalidOp; 1557 1558 case PackCategoriesIntoKey(fcNormal, fcNormal): 1559 return opOK; 1560 } 1561 } 1562 1563 APFloat::opStatus 1564 APFloat::divideSpecials(const APFloat &rhs) 1565 { 1566 switch (PackCategoriesIntoKey(category, rhs.category)) { 1567 default: 1568 llvm_unreachable(nullptr); 1569 1570 case PackCategoriesIntoKey(fcZero, fcNaN): 1571 case PackCategoriesIntoKey(fcNormal, fcNaN): 1572 case PackCategoriesIntoKey(fcInfinity, fcNaN): 1573 category = fcNaN; 1574 copySignificand(rhs); 1575 case PackCategoriesIntoKey(fcNaN, fcZero): 1576 case PackCategoriesIntoKey(fcNaN, fcNormal): 1577 case PackCategoriesIntoKey(fcNaN, fcInfinity): 1578 case PackCategoriesIntoKey(fcNaN, fcNaN): 1579 sign = false; 1580 case PackCategoriesIntoKey(fcInfinity, fcZero): 1581 case PackCategoriesIntoKey(fcInfinity, fcNormal): 1582 case PackCategoriesIntoKey(fcZero, fcInfinity): 1583 case PackCategoriesIntoKey(fcZero, fcNormal): 1584 return opOK; 1585 1586 case PackCategoriesIntoKey(fcNormal, fcInfinity): 1587 category = fcZero; 1588 return opOK; 1589 1590 case PackCategoriesIntoKey(fcNormal, fcZero): 1591 category = fcInfinity; 1592 return opDivByZero; 1593 1594 case PackCategoriesIntoKey(fcInfinity, fcInfinity): 1595 case PackCategoriesIntoKey(fcZero, fcZero): 1596 makeNaN(); 1597 return opInvalidOp; 1598 1599 case PackCategoriesIntoKey(fcNormal, fcNormal): 1600 return opOK; 1601 } 1602 } 1603 1604 APFloat::opStatus 1605 APFloat::modSpecials(const APFloat &rhs) 1606 { 1607 switch (PackCategoriesIntoKey(category, rhs.category)) { 1608 default: 1609 llvm_unreachable(nullptr); 1610 1611 case PackCategoriesIntoKey(fcNaN, fcZero): 1612 case PackCategoriesIntoKey(fcNaN, fcNormal): 1613 case PackCategoriesIntoKey(fcNaN, fcInfinity): 1614 case PackCategoriesIntoKey(fcNaN, fcNaN): 1615 case PackCategoriesIntoKey(fcZero, fcInfinity): 1616 case PackCategoriesIntoKey(fcZero, fcNormal): 1617 case PackCategoriesIntoKey(fcNormal, fcInfinity): 1618 return opOK; 1619 1620 case PackCategoriesIntoKey(fcZero, fcNaN): 1621 case PackCategoriesIntoKey(fcNormal, fcNaN): 1622 case PackCategoriesIntoKey(fcInfinity, fcNaN): 1623 sign = false; 1624 category = fcNaN; 1625 copySignificand(rhs); 1626 return opOK; 1627 1628 case PackCategoriesIntoKey(fcNormal, fcZero): 1629 case PackCategoriesIntoKey(fcInfinity, fcZero): 1630 case PackCategoriesIntoKey(fcInfinity, fcNormal): 1631 case PackCategoriesIntoKey(fcInfinity, fcInfinity): 1632 case PackCategoriesIntoKey(fcZero, fcZero): 1633 makeNaN(); 1634 return opInvalidOp; 1635 1636 case PackCategoriesIntoKey(fcNormal, fcNormal): 1637 return opOK; 1638 } 1639 } 1640 1641 /* Change sign. */ 1642 void 1643 APFloat::changeSign() 1644 { 1645 /* Look mummy, this one's easy. */ 1646 sign = !sign; 1647 } 1648 1649 void 1650 APFloat::clearSign() 1651 { 1652 /* So is this one. */ 1653 sign = 0; 1654 } 1655 1656 void 1657 APFloat::copySign(const APFloat &rhs) 1658 { 1659 /* And this one. */ 1660 sign = rhs.sign; 1661 } 1662 1663 /* Normalized addition or subtraction. */ 1664 APFloat::opStatus 1665 APFloat::addOrSubtract(const APFloat &rhs, roundingMode rounding_mode, 1666 bool subtract) 1667 { 1668 opStatus fs; 1669 1670 fs = addOrSubtractSpecials(rhs, subtract); 1671 1672 /* This return code means it was not a simple case. */ 1673 if (fs == opDivByZero) { 1674 lostFraction lost_fraction; 1675 1676 lost_fraction = addOrSubtractSignificand(rhs, subtract); 1677 fs = normalize(rounding_mode, lost_fraction); 1678 1679 /* Can only be zero if we lost no fraction. */ 1680 assert(category != fcZero || lost_fraction == lfExactlyZero); 1681 } 1682 1683 /* If two numbers add (exactly) to zero, IEEE 754 decrees it is a 1684 positive zero unless rounding to minus infinity, except that 1685 adding two like-signed zeroes gives that zero. */ 1686 if (category == fcZero) { 1687 if (rhs.category != fcZero || (sign == rhs.sign) == subtract) 1688 sign = (rounding_mode == rmTowardNegative); 1689 } 1690 1691 return fs; 1692 } 1693 1694 /* Normalized addition. */ 1695 APFloat::opStatus 1696 APFloat::add(const APFloat &rhs, roundingMode rounding_mode) 1697 { 1698 return addOrSubtract(rhs, rounding_mode, false); 1699 } 1700 1701 /* Normalized subtraction. */ 1702 APFloat::opStatus 1703 APFloat::subtract(const APFloat &rhs, roundingMode rounding_mode) 1704 { 1705 return addOrSubtract(rhs, rounding_mode, true); 1706 } 1707 1708 /* Normalized multiply. */ 1709 APFloat::opStatus 1710 APFloat::multiply(const APFloat &rhs, roundingMode rounding_mode) 1711 { 1712 opStatus fs; 1713 1714 sign ^= rhs.sign; 1715 fs = multiplySpecials(rhs); 1716 1717 if (isFiniteNonZero()) { 1718 lostFraction lost_fraction = multiplySignificand(rhs, nullptr); 1719 fs = normalize(rounding_mode, lost_fraction); 1720 if (lost_fraction != lfExactlyZero) 1721 fs = (opStatus) (fs | opInexact); 1722 } 1723 1724 return fs; 1725 } 1726 1727 /* Normalized divide. */ 1728 APFloat::opStatus 1729 APFloat::divide(const APFloat &rhs, roundingMode rounding_mode) 1730 { 1731 opStatus fs; 1732 1733 sign ^= rhs.sign; 1734 fs = divideSpecials(rhs); 1735 1736 if (isFiniteNonZero()) { 1737 lostFraction lost_fraction = divideSignificand(rhs); 1738 fs = normalize(rounding_mode, lost_fraction); 1739 if (lost_fraction != lfExactlyZero) 1740 fs = (opStatus) (fs | opInexact); 1741 } 1742 1743 return fs; 1744 } 1745 1746 /* Normalized remainder. This is not currently correct in all cases. */ 1747 APFloat::opStatus 1748 APFloat::remainder(const APFloat &rhs) 1749 { 1750 opStatus fs; 1751 APFloat V = *this; 1752 unsigned int origSign = sign; 1753 1754 fs = V.divide(rhs, rmNearestTiesToEven); 1755 if (fs == opDivByZero) 1756 return fs; 1757 1758 int parts = partCount(); 1759 integerPart *x = new integerPart[parts]; 1760 bool ignored; 1761 fs = V.convertToInteger(x, parts * integerPartWidth, true, 1762 rmNearestTiesToEven, &ignored); 1763 if (fs==opInvalidOp) 1764 return fs; 1765 1766 fs = V.convertFromZeroExtendedInteger(x, parts * integerPartWidth, true, 1767 rmNearestTiesToEven); 1768 assert(fs==opOK); // should always work 1769 1770 fs = V.multiply(rhs, rmNearestTiesToEven); 1771 assert(fs==opOK || fs==opInexact); // should not overflow or underflow 1772 1773 fs = subtract(V, rmNearestTiesToEven); 1774 assert(fs==opOK || fs==opInexact); // likewise 1775 1776 if (isZero()) 1777 sign = origSign; // IEEE754 requires this 1778 delete[] x; 1779 return fs; 1780 } 1781 1782 /* Normalized llvm frem (C fmod). 1783 This is not currently correct in all cases. */ 1784 APFloat::opStatus 1785 APFloat::mod(const APFloat &rhs) 1786 { 1787 opStatus fs; 1788 fs = modSpecials(rhs); 1789 1790 if (isFiniteNonZero() && rhs.isFiniteNonZero()) { 1791 APFloat V = *this; 1792 unsigned int origSign = sign; 1793 1794 fs = V.divide(rhs, rmNearestTiesToEven); 1795 if (fs == opDivByZero) 1796 return fs; 1797 1798 int parts = partCount(); 1799 integerPart *x = new integerPart[parts]; 1800 bool ignored; 1801 fs = V.convertToInteger(x, parts * integerPartWidth, true, 1802 rmTowardZero, &ignored); 1803 if (fs==opInvalidOp) 1804 return fs; 1805 1806 fs = V.convertFromZeroExtendedInteger(x, parts * integerPartWidth, true, 1807 rmNearestTiesToEven); 1808 assert(fs==opOK); // should always work 1809 1810 fs = V.multiply(rhs, rmNearestTiesToEven); 1811 assert(fs==opOK || fs==opInexact); // should not overflow or underflow 1812 1813 fs = subtract(V, rmNearestTiesToEven); 1814 assert(fs==opOK || fs==opInexact); // likewise 1815 1816 if (isZero()) 1817 sign = origSign; // IEEE754 requires this 1818 delete[] x; 1819 } 1820 return fs; 1821 } 1822 1823 /* Normalized fused-multiply-add. */ 1824 APFloat::opStatus 1825 APFloat::fusedMultiplyAdd(const APFloat &multiplicand, 1826 const APFloat &addend, 1827 roundingMode rounding_mode) 1828 { 1829 opStatus fs; 1830 1831 /* Post-multiplication sign, before addition. */ 1832 sign ^= multiplicand.sign; 1833 1834 /* If and only if all arguments are normal do we need to do an 1835 extended-precision calculation. */ 1836 if (isFiniteNonZero() && 1837 multiplicand.isFiniteNonZero() && 1838 addend.isFinite()) { 1839 lostFraction lost_fraction; 1840 1841 lost_fraction = multiplySignificand(multiplicand, &addend); 1842 fs = normalize(rounding_mode, lost_fraction); 1843 if (lost_fraction != lfExactlyZero) 1844 fs = (opStatus) (fs | opInexact); 1845 1846 /* If two numbers add (exactly) to zero, IEEE 754 decrees it is a 1847 positive zero unless rounding to minus infinity, except that 1848 adding two like-signed zeroes gives that zero. */ 1849 if (category == fcZero && !(fs & opUnderflow) && sign != addend.sign) 1850 sign = (rounding_mode == rmTowardNegative); 1851 } else { 1852 fs = multiplySpecials(multiplicand); 1853 1854 /* FS can only be opOK or opInvalidOp. There is no more work 1855 to do in the latter case. The IEEE-754R standard says it is 1856 implementation-defined in this case whether, if ADDEND is a 1857 quiet NaN, we raise invalid op; this implementation does so. 1858 1859 If we need to do the addition we can do so with normal 1860 precision. */ 1861 if (fs == opOK) 1862 fs = addOrSubtract(addend, rounding_mode, false); 1863 } 1864 1865 return fs; 1866 } 1867 1868 /* Rounding-mode corrrect round to integral value. */ 1869 APFloat::opStatus APFloat::roundToIntegral(roundingMode rounding_mode) { 1870 opStatus fs; 1871 1872 // If the exponent is large enough, we know that this value is already 1873 // integral, and the arithmetic below would potentially cause it to saturate 1874 // to +/-Inf. Bail out early instead. 1875 if (isFiniteNonZero() && exponent+1 >= (int)semanticsPrecision(*semantics)) 1876 return opOK; 1877 1878 // The algorithm here is quite simple: we add 2^(p-1), where p is the 1879 // precision of our format, and then subtract it back off again. The choice 1880 // of rounding modes for the addition/subtraction determines the rounding mode 1881 // for our integral rounding as well. 1882 // NOTE: When the input value is negative, we do subtraction followed by 1883 // addition instead. 1884 APInt IntegerConstant(NextPowerOf2(semanticsPrecision(*semantics)), 1); 1885 IntegerConstant <<= semanticsPrecision(*semantics)-1; 1886 APFloat MagicConstant(*semantics); 1887 fs = MagicConstant.convertFromAPInt(IntegerConstant, false, 1888 rmNearestTiesToEven); 1889 MagicConstant.copySign(*this); 1890 1891 if (fs != opOK) 1892 return fs; 1893 1894 // Preserve the input sign so that we can handle 0.0/-0.0 cases correctly. 1895 bool inputSign = isNegative(); 1896 1897 fs = add(MagicConstant, rounding_mode); 1898 if (fs != opOK && fs != opInexact) 1899 return fs; 1900 1901 fs = subtract(MagicConstant, rounding_mode); 1902 1903 // Restore the input sign. 1904 if (inputSign != isNegative()) 1905 changeSign(); 1906 1907 return fs; 1908 } 1909 1910 1911 /* Comparison requires normalized numbers. */ 1912 APFloat::cmpResult 1913 APFloat::compare(const APFloat &rhs) const 1914 { 1915 cmpResult result; 1916 1917 assert(semantics == rhs.semantics); 1918 1919 switch (PackCategoriesIntoKey(category, rhs.category)) { 1920 default: 1921 llvm_unreachable(nullptr); 1922 1923 case PackCategoriesIntoKey(fcNaN, fcZero): 1924 case PackCategoriesIntoKey(fcNaN, fcNormal): 1925 case PackCategoriesIntoKey(fcNaN, fcInfinity): 1926 case PackCategoriesIntoKey(fcNaN, fcNaN): 1927 case PackCategoriesIntoKey(fcZero, fcNaN): 1928 case PackCategoriesIntoKey(fcNormal, fcNaN): 1929 case PackCategoriesIntoKey(fcInfinity, fcNaN): 1930 return cmpUnordered; 1931 1932 case PackCategoriesIntoKey(fcInfinity, fcNormal): 1933 case PackCategoriesIntoKey(fcInfinity, fcZero): 1934 case PackCategoriesIntoKey(fcNormal, fcZero): 1935 if (sign) 1936 return cmpLessThan; 1937 else 1938 return cmpGreaterThan; 1939 1940 case PackCategoriesIntoKey(fcNormal, fcInfinity): 1941 case PackCategoriesIntoKey(fcZero, fcInfinity): 1942 case PackCategoriesIntoKey(fcZero, fcNormal): 1943 if (rhs.sign) 1944 return cmpGreaterThan; 1945 else 1946 return cmpLessThan; 1947 1948 case PackCategoriesIntoKey(fcInfinity, fcInfinity): 1949 if (sign == rhs.sign) 1950 return cmpEqual; 1951 else if (sign) 1952 return cmpLessThan; 1953 else 1954 return cmpGreaterThan; 1955 1956 case PackCategoriesIntoKey(fcZero, fcZero): 1957 return cmpEqual; 1958 1959 case PackCategoriesIntoKey(fcNormal, fcNormal): 1960 break; 1961 } 1962 1963 /* Two normal numbers. Do they have the same sign? */ 1964 if (sign != rhs.sign) { 1965 if (sign) 1966 result = cmpLessThan; 1967 else 1968 result = cmpGreaterThan; 1969 } else { 1970 /* Compare absolute values; invert result if negative. */ 1971 result = compareAbsoluteValue(rhs); 1972 1973 if (sign) { 1974 if (result == cmpLessThan) 1975 result = cmpGreaterThan; 1976 else if (result == cmpGreaterThan) 1977 result = cmpLessThan; 1978 } 1979 } 1980 1981 return result; 1982 } 1983 1984 /// APFloat::convert - convert a value of one floating point type to another. 1985 /// The return value corresponds to the IEEE754 exceptions. *losesInfo 1986 /// records whether the transformation lost information, i.e. whether 1987 /// converting the result back to the original type will produce the 1988 /// original value (this is almost the same as return value==fsOK, but there 1989 /// are edge cases where this is not so). 1990 1991 APFloat::opStatus 1992 APFloat::convert(const fltSemantics &toSemantics, 1993 roundingMode rounding_mode, bool *losesInfo) 1994 { 1995 lostFraction lostFraction; 1996 unsigned int newPartCount, oldPartCount; 1997 opStatus fs; 1998 int shift; 1999 const fltSemantics &fromSemantics = *semantics; 2000 2001 lostFraction = lfExactlyZero; 2002 newPartCount = partCountForBits(toSemantics.precision + 1); 2003 oldPartCount = partCount(); 2004 shift = toSemantics.precision - fromSemantics.precision; 2005 2006 bool X86SpecialNan = false; 2007 if (&fromSemantics == &APFloat::x87DoubleExtended && 2008 &toSemantics != &APFloat::x87DoubleExtended && category == fcNaN && 2009 (!(*significandParts() & 0x8000000000000000ULL) || 2010 !(*significandParts() & 0x4000000000000000ULL))) { 2011 // x86 has some unusual NaNs which cannot be represented in any other 2012 // format; note them here. 2013 X86SpecialNan = true; 2014 } 2015 2016 // If this is a truncation of a denormal number, and the target semantics 2017 // has larger exponent range than the source semantics (this can happen 2018 // when truncating from PowerPC double-double to double format), the 2019 // right shift could lose result mantissa bits. Adjust exponent instead 2020 // of performing excessive shift. 2021 if (shift < 0 && isFiniteNonZero()) { 2022 int exponentChange = significandMSB() + 1 - fromSemantics.precision; 2023 if (exponent + exponentChange < toSemantics.minExponent) 2024 exponentChange = toSemantics.minExponent - exponent; 2025 if (exponentChange < shift) 2026 exponentChange = shift; 2027 if (exponentChange < 0) { 2028 shift -= exponentChange; 2029 exponent += exponentChange; 2030 } 2031 } 2032 2033 // If this is a truncation, perform the shift before we narrow the storage. 2034 if (shift < 0 && (isFiniteNonZero() || category==fcNaN)) 2035 lostFraction = shiftRight(significandParts(), oldPartCount, -shift); 2036 2037 // Fix the storage so it can hold to new value. 2038 if (newPartCount > oldPartCount) { 2039 // The new type requires more storage; make it available. 2040 integerPart *newParts; 2041 newParts = new integerPart[newPartCount]; 2042 APInt::tcSet(newParts, 0, newPartCount); 2043 if (isFiniteNonZero() || category==fcNaN) 2044 APInt::tcAssign(newParts, significandParts(), oldPartCount); 2045 freeSignificand(); 2046 significand.parts = newParts; 2047 } else if (newPartCount == 1 && oldPartCount != 1) { 2048 // Switch to built-in storage for a single part. 2049 integerPart newPart = 0; 2050 if (isFiniteNonZero() || category==fcNaN) 2051 newPart = significandParts()[0]; 2052 freeSignificand(); 2053 significand.part = newPart; 2054 } 2055 2056 // Now that we have the right storage, switch the semantics. 2057 semantics = &toSemantics; 2058 2059 // If this is an extension, perform the shift now that the storage is 2060 // available. 2061 if (shift > 0 && (isFiniteNonZero() || category==fcNaN)) 2062 APInt::tcShiftLeft(significandParts(), newPartCount, shift); 2063 2064 if (isFiniteNonZero()) { 2065 fs = normalize(rounding_mode, lostFraction); 2066 *losesInfo = (fs != opOK); 2067 } else if (category == fcNaN) { 2068 *losesInfo = lostFraction != lfExactlyZero || X86SpecialNan; 2069 2070 // For x87 extended precision, we want to make a NaN, not a special NaN if 2071 // the input wasn't special either. 2072 if (!X86SpecialNan && semantics == &APFloat::x87DoubleExtended) 2073 APInt::tcSetBit(significandParts(), semantics->precision - 1); 2074 2075 // gcc forces the Quiet bit on, which means (float)(double)(float_sNan) 2076 // does not give you back the same bits. This is dubious, and we 2077 // don't currently do it. You're really supposed to get 2078 // an invalid operation signal at runtime, but nobody does that. 2079 fs = opOK; 2080 } else { 2081 *losesInfo = false; 2082 fs = opOK; 2083 } 2084 2085 return fs; 2086 } 2087 2088 /* Convert a floating point number to an integer according to the 2089 rounding mode. If the rounded integer value is out of range this 2090 returns an invalid operation exception and the contents of the 2091 destination parts are unspecified. If the rounded value is in 2092 range but the floating point number is not the exact integer, the C 2093 standard doesn't require an inexact exception to be raised. IEEE 2094 854 does require it so we do that. 2095 2096 Note that for conversions to integer type the C standard requires 2097 round-to-zero to always be used. */ 2098 APFloat::opStatus 2099 APFloat::convertToSignExtendedInteger(integerPart *parts, unsigned int width, 2100 bool isSigned, 2101 roundingMode rounding_mode, 2102 bool *isExact) const 2103 { 2104 lostFraction lost_fraction; 2105 const integerPart *src; 2106 unsigned int dstPartsCount, truncatedBits; 2107 2108 *isExact = false; 2109 2110 /* Handle the three special cases first. */ 2111 if (category == fcInfinity || category == fcNaN) 2112 return opInvalidOp; 2113 2114 dstPartsCount = partCountForBits(width); 2115 2116 if (category == fcZero) { 2117 APInt::tcSet(parts, 0, dstPartsCount); 2118 // Negative zero can't be represented as an int. 2119 *isExact = !sign; 2120 return opOK; 2121 } 2122 2123 src = significandParts(); 2124 2125 /* Step 1: place our absolute value, with any fraction truncated, in 2126 the destination. */ 2127 if (exponent < 0) { 2128 /* Our absolute value is less than one; truncate everything. */ 2129 APInt::tcSet(parts, 0, dstPartsCount); 2130 /* For exponent -1 the integer bit represents .5, look at that. 2131 For smaller exponents leftmost truncated bit is 0. */ 2132 truncatedBits = semantics->precision -1U - exponent; 2133 } else { 2134 /* We want the most significant (exponent + 1) bits; the rest are 2135 truncated. */ 2136 unsigned int bits = exponent + 1U; 2137 2138 /* Hopelessly large in magnitude? */ 2139 if (bits > width) 2140 return opInvalidOp; 2141 2142 if (bits < semantics->precision) { 2143 /* We truncate (semantics->precision - bits) bits. */ 2144 truncatedBits = semantics->precision - bits; 2145 APInt::tcExtract(parts, dstPartsCount, src, bits, truncatedBits); 2146 } else { 2147 /* We want at least as many bits as are available. */ 2148 APInt::tcExtract(parts, dstPartsCount, src, semantics->precision, 0); 2149 APInt::tcShiftLeft(parts, dstPartsCount, bits - semantics->precision); 2150 truncatedBits = 0; 2151 } 2152 } 2153 2154 /* Step 2: work out any lost fraction, and increment the absolute 2155 value if we would round away from zero. */ 2156 if (truncatedBits) { 2157 lost_fraction = lostFractionThroughTruncation(src, partCount(), 2158 truncatedBits); 2159 if (lost_fraction != lfExactlyZero && 2160 roundAwayFromZero(rounding_mode, lost_fraction, truncatedBits)) { 2161 if (APInt::tcIncrement(parts, dstPartsCount)) 2162 return opInvalidOp; /* Overflow. */ 2163 } 2164 } else { 2165 lost_fraction = lfExactlyZero; 2166 } 2167 2168 /* Step 3: check if we fit in the destination. */ 2169 unsigned int omsb = APInt::tcMSB(parts, dstPartsCount) + 1; 2170 2171 if (sign) { 2172 if (!isSigned) { 2173 /* Negative numbers cannot be represented as unsigned. */ 2174 if (omsb != 0) 2175 return opInvalidOp; 2176 } else { 2177 /* It takes omsb bits to represent the unsigned integer value. 2178 We lose a bit for the sign, but care is needed as the 2179 maximally negative integer is a special case. */ 2180 if (omsb == width && APInt::tcLSB(parts, dstPartsCount) + 1 != omsb) 2181 return opInvalidOp; 2182 2183 /* This case can happen because of rounding. */ 2184 if (omsb > width) 2185 return opInvalidOp; 2186 } 2187 2188 APInt::tcNegate (parts, dstPartsCount); 2189 } else { 2190 if (omsb >= width + !isSigned) 2191 return opInvalidOp; 2192 } 2193 2194 if (lost_fraction == lfExactlyZero) { 2195 *isExact = true; 2196 return opOK; 2197 } else 2198 return opInexact; 2199 } 2200 2201 /* Same as convertToSignExtendedInteger, except we provide 2202 deterministic values in case of an invalid operation exception, 2203 namely zero for NaNs and the minimal or maximal value respectively 2204 for underflow or overflow. 2205 The *isExact output tells whether the result is exact, in the sense 2206 that converting it back to the original floating point type produces 2207 the original value. This is almost equivalent to result==opOK, 2208 except for negative zeroes. 2209 */ 2210 APFloat::opStatus 2211 APFloat::convertToInteger(integerPart *parts, unsigned int width, 2212 bool isSigned, 2213 roundingMode rounding_mode, bool *isExact) const 2214 { 2215 opStatus fs; 2216 2217 fs = convertToSignExtendedInteger(parts, width, isSigned, rounding_mode, 2218 isExact); 2219 2220 if (fs == opInvalidOp) { 2221 unsigned int bits, dstPartsCount; 2222 2223 dstPartsCount = partCountForBits(width); 2224 2225 if (category == fcNaN) 2226 bits = 0; 2227 else if (sign) 2228 bits = isSigned; 2229 else 2230 bits = width - isSigned; 2231 2232 APInt::tcSetLeastSignificantBits(parts, dstPartsCount, bits); 2233 if (sign && isSigned) 2234 APInt::tcShiftLeft(parts, dstPartsCount, width - 1); 2235 } 2236 2237 return fs; 2238 } 2239 2240 /* Same as convertToInteger(integerPart*, ...), except the result is returned in 2241 an APSInt, whose initial bit-width and signed-ness are used to determine the 2242 precision of the conversion. 2243 */ 2244 APFloat::opStatus 2245 APFloat::convertToInteger(APSInt &result, 2246 roundingMode rounding_mode, bool *isExact) const 2247 { 2248 unsigned bitWidth = result.getBitWidth(); 2249 SmallVector<uint64_t, 4> parts(result.getNumWords()); 2250 opStatus status = convertToInteger( 2251 parts.data(), bitWidth, result.isSigned(), rounding_mode, isExact); 2252 // Keeps the original signed-ness. 2253 result = APInt(bitWidth, parts); 2254 return status; 2255 } 2256 2257 /* Convert an unsigned integer SRC to a floating point number, 2258 rounding according to ROUNDING_MODE. The sign of the floating 2259 point number is not modified. */ 2260 APFloat::opStatus 2261 APFloat::convertFromUnsignedParts(const integerPart *src, 2262 unsigned int srcCount, 2263 roundingMode rounding_mode) 2264 { 2265 unsigned int omsb, precision, dstCount; 2266 integerPart *dst; 2267 lostFraction lost_fraction; 2268 2269 category = fcNormal; 2270 omsb = APInt::tcMSB(src, srcCount) + 1; 2271 dst = significandParts(); 2272 dstCount = partCount(); 2273 precision = semantics->precision; 2274 2275 /* We want the most significant PRECISION bits of SRC. There may not 2276 be that many; extract what we can. */ 2277 if (precision <= omsb) { 2278 exponent = omsb - 1; 2279 lost_fraction = lostFractionThroughTruncation(src, srcCount, 2280 omsb - precision); 2281 APInt::tcExtract(dst, dstCount, src, precision, omsb - precision); 2282 } else { 2283 exponent = precision - 1; 2284 lost_fraction = lfExactlyZero; 2285 APInt::tcExtract(dst, dstCount, src, omsb, 0); 2286 } 2287 2288 return normalize(rounding_mode, lost_fraction); 2289 } 2290 2291 APFloat::opStatus 2292 APFloat::convertFromAPInt(const APInt &Val, 2293 bool isSigned, 2294 roundingMode rounding_mode) 2295 { 2296 unsigned int partCount = Val.getNumWords(); 2297 APInt api = Val; 2298 2299 sign = false; 2300 if (isSigned && api.isNegative()) { 2301 sign = true; 2302 api = -api; 2303 } 2304 2305 return convertFromUnsignedParts(api.getRawData(), partCount, rounding_mode); 2306 } 2307 2308 /* Convert a two's complement integer SRC to a floating point number, 2309 rounding according to ROUNDING_MODE. ISSIGNED is true if the 2310 integer is signed, in which case it must be sign-extended. */ 2311 APFloat::opStatus 2312 APFloat::convertFromSignExtendedInteger(const integerPart *src, 2313 unsigned int srcCount, 2314 bool isSigned, 2315 roundingMode rounding_mode) 2316 { 2317 opStatus status; 2318 2319 if (isSigned && 2320 APInt::tcExtractBit(src, srcCount * integerPartWidth - 1)) { 2321 integerPart *copy; 2322 2323 /* If we're signed and negative negate a copy. */ 2324 sign = true; 2325 copy = new integerPart[srcCount]; 2326 APInt::tcAssign(copy, src, srcCount); 2327 APInt::tcNegate(copy, srcCount); 2328 status = convertFromUnsignedParts(copy, srcCount, rounding_mode); 2329 delete [] copy; 2330 } else { 2331 sign = false; 2332 status = convertFromUnsignedParts(src, srcCount, rounding_mode); 2333 } 2334 2335 return status; 2336 } 2337 2338 /* FIXME: should this just take a const APInt reference? */ 2339 APFloat::opStatus 2340 APFloat::convertFromZeroExtendedInteger(const integerPart *parts, 2341 unsigned int width, bool isSigned, 2342 roundingMode rounding_mode) 2343 { 2344 unsigned int partCount = partCountForBits(width); 2345 APInt api = APInt(width, makeArrayRef(parts, partCount)); 2346 2347 sign = false; 2348 if (isSigned && APInt::tcExtractBit(parts, width - 1)) { 2349 sign = true; 2350 api = -api; 2351 } 2352 2353 return convertFromUnsignedParts(api.getRawData(), partCount, rounding_mode); 2354 } 2355 2356 APFloat::opStatus 2357 APFloat::convertFromHexadecimalString(StringRef s, roundingMode rounding_mode) 2358 { 2359 lostFraction lost_fraction = lfExactlyZero; 2360 2361 category = fcNormal; 2362 zeroSignificand(); 2363 exponent = 0; 2364 2365 integerPart *significand = significandParts(); 2366 unsigned partsCount = partCount(); 2367 unsigned bitPos = partsCount * integerPartWidth; 2368 bool computedTrailingFraction = false; 2369 2370 // Skip leading zeroes and any (hexa)decimal point. 2371 StringRef::iterator begin = s.begin(); 2372 StringRef::iterator end = s.end(); 2373 StringRef::iterator dot; 2374 StringRef::iterator p = skipLeadingZeroesAndAnyDot(begin, end, &dot); 2375 StringRef::iterator firstSignificantDigit = p; 2376 2377 while (p != end) { 2378 integerPart hex_value; 2379 2380 if (*p == '.') { 2381 assert(dot == end && "String contains multiple dots"); 2382 dot = p++; 2383 continue; 2384 } 2385 2386 hex_value = hexDigitValue(*p); 2387 if (hex_value == -1U) 2388 break; 2389 2390 p++; 2391 2392 // Store the number while we have space. 2393 if (bitPos) { 2394 bitPos -= 4; 2395 hex_value <<= bitPos % integerPartWidth; 2396 significand[bitPos / integerPartWidth] |= hex_value; 2397 } else if (!computedTrailingFraction) { 2398 lost_fraction = trailingHexadecimalFraction(p, end, hex_value); 2399 computedTrailingFraction = true; 2400 } 2401 } 2402 2403 /* Hex floats require an exponent but not a hexadecimal point. */ 2404 assert(p != end && "Hex strings require an exponent"); 2405 assert((*p == 'p' || *p == 'P') && "Invalid character in significand"); 2406 assert(p != begin && "Significand has no digits"); 2407 assert((dot == end || p - begin != 1) && "Significand has no digits"); 2408 2409 /* Ignore the exponent if we are zero. */ 2410 if (p != firstSignificantDigit) { 2411 int expAdjustment; 2412 2413 /* Implicit hexadecimal point? */ 2414 if (dot == end) 2415 dot = p; 2416 2417 /* Calculate the exponent adjustment implicit in the number of 2418 significant digits. */ 2419 expAdjustment = static_cast<int>(dot - firstSignificantDigit); 2420 if (expAdjustment < 0) 2421 expAdjustment++; 2422 expAdjustment = expAdjustment * 4 - 1; 2423 2424 /* Adjust for writing the significand starting at the most 2425 significant nibble. */ 2426 expAdjustment += semantics->precision; 2427 expAdjustment -= partsCount * integerPartWidth; 2428 2429 /* Adjust for the given exponent. */ 2430 exponent = totalExponent(p + 1, end, expAdjustment); 2431 } 2432 2433 return normalize(rounding_mode, lost_fraction); 2434 } 2435 2436 APFloat::opStatus 2437 APFloat::roundSignificandWithExponent(const integerPart *decSigParts, 2438 unsigned sigPartCount, int exp, 2439 roundingMode rounding_mode) 2440 { 2441 unsigned int parts, pow5PartCount; 2442 fltSemantics calcSemantics = { 32767, -32767, 0, 0 }; 2443 integerPart pow5Parts[maxPowerOfFiveParts]; 2444 bool isNearest; 2445 2446 isNearest = (rounding_mode == rmNearestTiesToEven || 2447 rounding_mode == rmNearestTiesToAway); 2448 2449 parts = partCountForBits(semantics->precision + 11); 2450 2451 /* Calculate pow(5, abs(exp)). */ 2452 pow5PartCount = powerOf5(pow5Parts, exp >= 0 ? exp: -exp); 2453 2454 for (;; parts *= 2) { 2455 opStatus sigStatus, powStatus; 2456 unsigned int excessPrecision, truncatedBits; 2457 2458 calcSemantics.precision = parts * integerPartWidth - 1; 2459 excessPrecision = calcSemantics.precision - semantics->precision; 2460 truncatedBits = excessPrecision; 2461 2462 APFloat decSig = APFloat::getZero(calcSemantics, sign); 2463 APFloat pow5(calcSemantics); 2464 2465 sigStatus = decSig.convertFromUnsignedParts(decSigParts, sigPartCount, 2466 rmNearestTiesToEven); 2467 powStatus = pow5.convertFromUnsignedParts(pow5Parts, pow5PartCount, 2468 rmNearestTiesToEven); 2469 /* Add exp, as 10^n = 5^n * 2^n. */ 2470 decSig.exponent += exp; 2471 2472 lostFraction calcLostFraction; 2473 integerPart HUerr, HUdistance; 2474 unsigned int powHUerr; 2475 2476 if (exp >= 0) { 2477 /* multiplySignificand leaves the precision-th bit set to 1. */ 2478 calcLostFraction = decSig.multiplySignificand(pow5, nullptr); 2479 powHUerr = powStatus != opOK; 2480 } else { 2481 calcLostFraction = decSig.divideSignificand(pow5); 2482 /* Denormal numbers have less precision. */ 2483 if (decSig.exponent < semantics->minExponent) { 2484 excessPrecision += (semantics->minExponent - decSig.exponent); 2485 truncatedBits = excessPrecision; 2486 if (excessPrecision > calcSemantics.precision) 2487 excessPrecision = calcSemantics.precision; 2488 } 2489 /* Extra half-ulp lost in reciprocal of exponent. */ 2490 powHUerr = (powStatus == opOK && calcLostFraction == lfExactlyZero) ? 0:2; 2491 } 2492 2493 /* Both multiplySignificand and divideSignificand return the 2494 result with the integer bit set. */ 2495 assert(APInt::tcExtractBit 2496 (decSig.significandParts(), calcSemantics.precision - 1) == 1); 2497 2498 HUerr = HUerrBound(calcLostFraction != lfExactlyZero, sigStatus != opOK, 2499 powHUerr); 2500 HUdistance = 2 * ulpsFromBoundary(decSig.significandParts(), 2501 excessPrecision, isNearest); 2502 2503 /* Are we guaranteed to round correctly if we truncate? */ 2504 if (HUdistance >= HUerr) { 2505 APInt::tcExtract(significandParts(), partCount(), decSig.significandParts(), 2506 calcSemantics.precision - excessPrecision, 2507 excessPrecision); 2508 /* Take the exponent of decSig. If we tcExtract-ed less bits 2509 above we must adjust our exponent to compensate for the 2510 implicit right shift. */ 2511 exponent = (decSig.exponent + semantics->precision 2512 - (calcSemantics.precision - excessPrecision)); 2513 calcLostFraction = lostFractionThroughTruncation(decSig.significandParts(), 2514 decSig.partCount(), 2515 truncatedBits); 2516 return normalize(rounding_mode, calcLostFraction); 2517 } 2518 } 2519 } 2520 2521 APFloat::opStatus 2522 APFloat::convertFromDecimalString(StringRef str, roundingMode rounding_mode) 2523 { 2524 decimalInfo D; 2525 opStatus fs; 2526 2527 /* Scan the text. */ 2528 StringRef::iterator p = str.begin(); 2529 interpretDecimal(p, str.end(), &D); 2530 2531 /* Handle the quick cases. First the case of no significant digits, 2532 i.e. zero, and then exponents that are obviously too large or too 2533 small. Writing L for log 10 / log 2, a number d.ddddd*10^exp 2534 definitely overflows if 2535 2536 (exp - 1) * L >= maxExponent 2537 2538 and definitely underflows to zero where 2539 2540 (exp + 1) * L <= minExponent - precision 2541 2542 With integer arithmetic the tightest bounds for L are 2543 2544 93/28 < L < 196/59 [ numerator <= 256 ] 2545 42039/12655 < L < 28738/8651 [ numerator <= 65536 ] 2546 */ 2547 2548 // Test if we have a zero number allowing for strings with no null terminators 2549 // and zero decimals with non-zero exponents. 2550 // 2551 // We computed firstSigDigit by ignoring all zeros and dots. Thus if 2552 // D->firstSigDigit equals str.end(), every digit must be a zero and there can 2553 // be at most one dot. On the other hand, if we have a zero with a non-zero 2554 // exponent, then we know that D.firstSigDigit will be non-numeric. 2555 if (D.firstSigDigit == str.end() || decDigitValue(*D.firstSigDigit) >= 10U) { 2556 category = fcZero; 2557 fs = opOK; 2558 2559 /* Check whether the normalized exponent is high enough to overflow 2560 max during the log-rebasing in the max-exponent check below. */ 2561 } else if (D.normalizedExponent - 1 > INT_MAX / 42039) { 2562 fs = handleOverflow(rounding_mode); 2563 2564 /* If it wasn't, then it also wasn't high enough to overflow max 2565 during the log-rebasing in the min-exponent check. Check that it 2566 won't overflow min in either check, then perform the min-exponent 2567 check. */ 2568 } else if (D.normalizedExponent - 1 < INT_MIN / 42039 || 2569 (D.normalizedExponent + 1) * 28738 <= 2570 8651 * (semantics->minExponent - (int) semantics->precision)) { 2571 /* Underflow to zero and round. */ 2572 category = fcNormal; 2573 zeroSignificand(); 2574 fs = normalize(rounding_mode, lfLessThanHalf); 2575 2576 /* We can finally safely perform the max-exponent check. */ 2577 } else if ((D.normalizedExponent - 1) * 42039 2578 >= 12655 * semantics->maxExponent) { 2579 /* Overflow and round. */ 2580 fs = handleOverflow(rounding_mode); 2581 } else { 2582 integerPart *decSignificand; 2583 unsigned int partCount; 2584 2585 /* A tight upper bound on number of bits required to hold an 2586 N-digit decimal integer is N * 196 / 59. Allocate enough space 2587 to hold the full significand, and an extra part required by 2588 tcMultiplyPart. */ 2589 partCount = static_cast<unsigned int>(D.lastSigDigit - D.firstSigDigit) + 1; 2590 partCount = partCountForBits(1 + 196 * partCount / 59); 2591 decSignificand = new integerPart[partCount + 1]; 2592 partCount = 0; 2593 2594 /* Convert to binary efficiently - we do almost all multiplication 2595 in an integerPart. When this would overflow do we do a single 2596 bignum multiplication, and then revert again to multiplication 2597 in an integerPart. */ 2598 do { 2599 integerPart decValue, val, multiplier; 2600 2601 val = 0; 2602 multiplier = 1; 2603 2604 do { 2605 if (*p == '.') { 2606 p++; 2607 if (p == str.end()) { 2608 break; 2609 } 2610 } 2611 decValue = decDigitValue(*p++); 2612 assert(decValue < 10U && "Invalid character in significand"); 2613 multiplier *= 10; 2614 val = val * 10 + decValue; 2615 /* The maximum number that can be multiplied by ten with any 2616 digit added without overflowing an integerPart. */ 2617 } while (p <= D.lastSigDigit && multiplier <= (~ (integerPart) 0 - 9) / 10); 2618 2619 /* Multiply out the current part. */ 2620 APInt::tcMultiplyPart(decSignificand, decSignificand, multiplier, val, 2621 partCount, partCount + 1, false); 2622 2623 /* If we used another part (likely but not guaranteed), increase 2624 the count. */ 2625 if (decSignificand[partCount]) 2626 partCount++; 2627 } while (p <= D.lastSigDigit); 2628 2629 category = fcNormal; 2630 fs = roundSignificandWithExponent(decSignificand, partCount, 2631 D.exponent, rounding_mode); 2632 2633 delete [] decSignificand; 2634 } 2635 2636 return fs; 2637 } 2638 2639 bool 2640 APFloat::convertFromStringSpecials(StringRef str) { 2641 if (str.equals("inf") || str.equals("INFINITY")) { 2642 makeInf(false); 2643 return true; 2644 } 2645 2646 if (str.equals("-inf") || str.equals("-INFINITY")) { 2647 makeInf(true); 2648 return true; 2649 } 2650 2651 if (str.equals("nan") || str.equals("NaN")) { 2652 makeNaN(false, false); 2653 return true; 2654 } 2655 2656 if (str.equals("-nan") || str.equals("-NaN")) { 2657 makeNaN(false, true); 2658 return true; 2659 } 2660 2661 return false; 2662 } 2663 2664 APFloat::opStatus 2665 APFloat::convertFromString(StringRef str, roundingMode rounding_mode) 2666 { 2667 assert(!str.empty() && "Invalid string length"); 2668 2669 // Handle special cases. 2670 if (convertFromStringSpecials(str)) 2671 return opOK; 2672 2673 /* Handle a leading minus sign. */ 2674 StringRef::iterator p = str.begin(); 2675 size_t slen = str.size(); 2676 sign = *p == '-' ? 1 : 0; 2677 if (*p == '-' || *p == '+') { 2678 p++; 2679 slen--; 2680 assert(slen && "String has no digits"); 2681 } 2682 2683 if (slen >= 2 && p[0] == '0' && (p[1] == 'x' || p[1] == 'X')) { 2684 assert(slen - 2 && "Invalid string"); 2685 return convertFromHexadecimalString(StringRef(p + 2, slen - 2), 2686 rounding_mode); 2687 } 2688 2689 return convertFromDecimalString(StringRef(p, slen), rounding_mode); 2690 } 2691 2692 /* Write out a hexadecimal representation of the floating point value 2693 to DST, which must be of sufficient size, in the C99 form 2694 [-]0xh.hhhhp[+-]d. Return the number of characters written, 2695 excluding the terminating NUL. 2696 2697 If UPPERCASE, the output is in upper case, otherwise in lower case. 2698 2699 HEXDIGITS digits appear altogether, rounding the value if 2700 necessary. If HEXDIGITS is 0, the minimal precision to display the 2701 number precisely is used instead. If nothing would appear after 2702 the decimal point it is suppressed. 2703 2704 The decimal exponent is always printed and has at least one digit. 2705 Zero values display an exponent of zero. Infinities and NaNs 2706 appear as "infinity" or "nan" respectively. 2707 2708 The above rules are as specified by C99. There is ambiguity about 2709 what the leading hexadecimal digit should be. This implementation 2710 uses whatever is necessary so that the exponent is displayed as 2711 stored. This implies the exponent will fall within the IEEE format 2712 range, and the leading hexadecimal digit will be 0 (for denormals), 2713 1 (normal numbers) or 2 (normal numbers rounded-away-from-zero with 2714 any other digits zero). 2715 */ 2716 unsigned int 2717 APFloat::convertToHexString(char *dst, unsigned int hexDigits, 2718 bool upperCase, roundingMode rounding_mode) const 2719 { 2720 char *p; 2721 2722 p = dst; 2723 if (sign) 2724 *dst++ = '-'; 2725 2726 switch (category) { 2727 case fcInfinity: 2728 memcpy (dst, upperCase ? infinityU: infinityL, sizeof infinityU - 1); 2729 dst += sizeof infinityL - 1; 2730 break; 2731 2732 case fcNaN: 2733 memcpy (dst, upperCase ? NaNU: NaNL, sizeof NaNU - 1); 2734 dst += sizeof NaNU - 1; 2735 break; 2736 2737 case fcZero: 2738 *dst++ = '0'; 2739 *dst++ = upperCase ? 'X': 'x'; 2740 *dst++ = '0'; 2741 if (hexDigits > 1) { 2742 *dst++ = '.'; 2743 memset (dst, '0', hexDigits - 1); 2744 dst += hexDigits - 1; 2745 } 2746 *dst++ = upperCase ? 'P': 'p'; 2747 *dst++ = '0'; 2748 break; 2749 2750 case fcNormal: 2751 dst = convertNormalToHexString (dst, hexDigits, upperCase, rounding_mode); 2752 break; 2753 } 2754 2755 *dst = 0; 2756 2757 return static_cast<unsigned int>(dst - p); 2758 } 2759 2760 /* Does the hard work of outputting the correctly rounded hexadecimal 2761 form of a normal floating point number with the specified number of 2762 hexadecimal digits. If HEXDIGITS is zero the minimum number of 2763 digits necessary to print the value precisely is output. */ 2764 char * 2765 APFloat::convertNormalToHexString(char *dst, unsigned int hexDigits, 2766 bool upperCase, 2767 roundingMode rounding_mode) const 2768 { 2769 unsigned int count, valueBits, shift, partsCount, outputDigits; 2770 const char *hexDigitChars; 2771 const integerPart *significand; 2772 char *p; 2773 bool roundUp; 2774 2775 *dst++ = '0'; 2776 *dst++ = upperCase ? 'X': 'x'; 2777 2778 roundUp = false; 2779 hexDigitChars = upperCase ? hexDigitsUpper: hexDigitsLower; 2780 2781 significand = significandParts(); 2782 partsCount = partCount(); 2783 2784 /* +3 because the first digit only uses the single integer bit, so 2785 we have 3 virtual zero most-significant-bits. */ 2786 valueBits = semantics->precision + 3; 2787 shift = integerPartWidth - valueBits % integerPartWidth; 2788 2789 /* The natural number of digits required ignoring trailing 2790 insignificant zeroes. */ 2791 outputDigits = (valueBits - significandLSB () + 3) / 4; 2792 2793 /* hexDigits of zero means use the required number for the 2794 precision. Otherwise, see if we are truncating. If we are, 2795 find out if we need to round away from zero. */ 2796 if (hexDigits) { 2797 if (hexDigits < outputDigits) { 2798 /* We are dropping non-zero bits, so need to check how to round. 2799 "bits" is the number of dropped bits. */ 2800 unsigned int bits; 2801 lostFraction fraction; 2802 2803 bits = valueBits - hexDigits * 4; 2804 fraction = lostFractionThroughTruncation (significand, partsCount, bits); 2805 roundUp = roundAwayFromZero(rounding_mode, fraction, bits); 2806 } 2807 outputDigits = hexDigits; 2808 } 2809 2810 /* Write the digits consecutively, and start writing in the location 2811 of the hexadecimal point. We move the most significant digit 2812 left and add the hexadecimal point later. */ 2813 p = ++dst; 2814 2815 count = (valueBits + integerPartWidth - 1) / integerPartWidth; 2816 2817 while (outputDigits && count) { 2818 integerPart part; 2819 2820 /* Put the most significant integerPartWidth bits in "part". */ 2821 if (--count == partsCount) 2822 part = 0; /* An imaginary higher zero part. */ 2823 else 2824 part = significand[count] << shift; 2825 2826 if (count && shift) 2827 part |= significand[count - 1] >> (integerPartWidth - shift); 2828 2829 /* Convert as much of "part" to hexdigits as we can. */ 2830 unsigned int curDigits = integerPartWidth / 4; 2831 2832 if (curDigits > outputDigits) 2833 curDigits = outputDigits; 2834 dst += partAsHex (dst, part, curDigits, hexDigitChars); 2835 outputDigits -= curDigits; 2836 } 2837 2838 if (roundUp) { 2839 char *q = dst; 2840 2841 /* Note that hexDigitChars has a trailing '0'. */ 2842 do { 2843 q--; 2844 *q = hexDigitChars[hexDigitValue (*q) + 1]; 2845 } while (*q == '0'); 2846 assert(q >= p); 2847 } else { 2848 /* Add trailing zeroes. */ 2849 memset (dst, '0', outputDigits); 2850 dst += outputDigits; 2851 } 2852 2853 /* Move the most significant digit to before the point, and if there 2854 is something after the decimal point add it. This must come 2855 after rounding above. */ 2856 p[-1] = p[0]; 2857 if (dst -1 == p) 2858 dst--; 2859 else 2860 p[0] = '.'; 2861 2862 /* Finally output the exponent. */ 2863 *dst++ = upperCase ? 'P': 'p'; 2864 2865 return writeSignedDecimal (dst, exponent); 2866 } 2867 2868 hash_code llvm::hash_value(const APFloat &Arg) { 2869 if (!Arg.isFiniteNonZero()) 2870 return hash_combine((uint8_t)Arg.category, 2871 // NaN has no sign, fix it at zero. 2872 Arg.isNaN() ? (uint8_t)0 : (uint8_t)Arg.sign, 2873 Arg.semantics->precision); 2874 2875 // Normal floats need their exponent and significand hashed. 2876 return hash_combine((uint8_t)Arg.category, (uint8_t)Arg.sign, 2877 Arg.semantics->precision, Arg.exponent, 2878 hash_combine_range( 2879 Arg.significandParts(), 2880 Arg.significandParts() + Arg.partCount())); 2881 } 2882 2883 // Conversion from APFloat to/from host float/double. It may eventually be 2884 // possible to eliminate these and have everybody deal with APFloats, but that 2885 // will take a while. This approach will not easily extend to long double. 2886 // Current implementation requires integerPartWidth==64, which is correct at 2887 // the moment but could be made more general. 2888 2889 // Denormals have exponent minExponent in APFloat, but minExponent-1 in 2890 // the actual IEEE respresentations. We compensate for that here. 2891 2892 APInt 2893 APFloat::convertF80LongDoubleAPFloatToAPInt() const 2894 { 2895 assert(semantics == (const llvm::fltSemantics*)&x87DoubleExtended); 2896 assert(partCount()==2); 2897 2898 uint64_t myexponent, mysignificand; 2899 2900 if (isFiniteNonZero()) { 2901 myexponent = exponent+16383; //bias 2902 mysignificand = significandParts()[0]; 2903 if (myexponent==1 && !(mysignificand & 0x8000000000000000ULL)) 2904 myexponent = 0; // denormal 2905 } else if (category==fcZero) { 2906 myexponent = 0; 2907 mysignificand = 0; 2908 } else if (category==fcInfinity) { 2909 myexponent = 0x7fff; 2910 mysignificand = 0x8000000000000000ULL; 2911 } else { 2912 assert(category == fcNaN && "Unknown category"); 2913 myexponent = 0x7fff; 2914 mysignificand = significandParts()[0]; 2915 } 2916 2917 uint64_t words[2]; 2918 words[0] = mysignificand; 2919 words[1] = ((uint64_t)(sign & 1) << 15) | 2920 (myexponent & 0x7fffLL); 2921 return APInt(80, words); 2922 } 2923 2924 APInt 2925 APFloat::convertPPCDoubleDoubleAPFloatToAPInt() const 2926 { 2927 assert(semantics == (const llvm::fltSemantics*)&PPCDoubleDouble); 2928 assert(partCount()==2); 2929 2930 uint64_t words[2]; 2931 opStatus fs; 2932 bool losesInfo; 2933 2934 // Convert number to double. To avoid spurious underflows, we re- 2935 // normalize against the "double" minExponent first, and only *then* 2936 // truncate the mantissa. The result of that second conversion 2937 // may be inexact, but should never underflow. 2938 // Declare fltSemantics before APFloat that uses it (and 2939 // saves pointer to it) to ensure correct destruction order. 2940 fltSemantics extendedSemantics = *semantics; 2941 extendedSemantics.minExponent = IEEEdouble.minExponent; 2942 APFloat extended(*this); 2943 fs = extended.convert(extendedSemantics, rmNearestTiesToEven, &losesInfo); 2944 assert(fs == opOK && !losesInfo); 2945 (void)fs; 2946 2947 APFloat u(extended); 2948 fs = u.convert(IEEEdouble, rmNearestTiesToEven, &losesInfo); 2949 assert(fs == opOK || fs == opInexact); 2950 (void)fs; 2951 words[0] = *u.convertDoubleAPFloatToAPInt().getRawData(); 2952 2953 // If conversion was exact or resulted in a special case, we're done; 2954 // just set the second double to zero. Otherwise, re-convert back to 2955 // the extended format and compute the difference. This now should 2956 // convert exactly to double. 2957 if (u.isFiniteNonZero() && losesInfo) { 2958 fs = u.convert(extendedSemantics, rmNearestTiesToEven, &losesInfo); 2959 assert(fs == opOK && !losesInfo); 2960 (void)fs; 2961 2962 APFloat v(extended); 2963 v.subtract(u, rmNearestTiesToEven); 2964 fs = v.convert(IEEEdouble, rmNearestTiesToEven, &losesInfo); 2965 assert(fs == opOK && !losesInfo); 2966 (void)fs; 2967 words[1] = *v.convertDoubleAPFloatToAPInt().getRawData(); 2968 } else { 2969 words[1] = 0; 2970 } 2971 2972 return APInt(128, words); 2973 } 2974 2975 APInt 2976 APFloat::convertQuadrupleAPFloatToAPInt() const 2977 { 2978 assert(semantics == (const llvm::fltSemantics*)&IEEEquad); 2979 assert(partCount()==2); 2980 2981 uint64_t myexponent, mysignificand, mysignificand2; 2982 2983 if (isFiniteNonZero()) { 2984 myexponent = exponent+16383; //bias 2985 mysignificand = significandParts()[0]; 2986 mysignificand2 = significandParts()[1]; 2987 if (myexponent==1 && !(mysignificand2 & 0x1000000000000LL)) 2988 myexponent = 0; // denormal 2989 } else if (category==fcZero) { 2990 myexponent = 0; 2991 mysignificand = mysignificand2 = 0; 2992 } else if (category==fcInfinity) { 2993 myexponent = 0x7fff; 2994 mysignificand = mysignificand2 = 0; 2995 } else { 2996 assert(category == fcNaN && "Unknown category!"); 2997 myexponent = 0x7fff; 2998 mysignificand = significandParts()[0]; 2999 mysignificand2 = significandParts()[1]; 3000 } 3001 3002 uint64_t words[2]; 3003 words[0] = mysignificand; 3004 words[1] = ((uint64_t)(sign & 1) << 63) | 3005 ((myexponent & 0x7fff) << 48) | 3006 (mysignificand2 & 0xffffffffffffLL); 3007 3008 return APInt(128, words); 3009 } 3010 3011 APInt 3012 APFloat::convertDoubleAPFloatToAPInt() const 3013 { 3014 assert(semantics == (const llvm::fltSemantics*)&IEEEdouble); 3015 assert(partCount()==1); 3016 3017 uint64_t myexponent, mysignificand; 3018 3019 if (isFiniteNonZero()) { 3020 myexponent = exponent+1023; //bias 3021 mysignificand = *significandParts(); 3022 if (myexponent==1 && !(mysignificand & 0x10000000000000LL)) 3023 myexponent = 0; // denormal 3024 } else if (category==fcZero) { 3025 myexponent = 0; 3026 mysignificand = 0; 3027 } else if (category==fcInfinity) { 3028 myexponent = 0x7ff; 3029 mysignificand = 0; 3030 } else { 3031 assert(category == fcNaN && "Unknown category!"); 3032 myexponent = 0x7ff; 3033 mysignificand = *significandParts(); 3034 } 3035 3036 return APInt(64, ((((uint64_t)(sign & 1) << 63) | 3037 ((myexponent & 0x7ff) << 52) | 3038 (mysignificand & 0xfffffffffffffLL)))); 3039 } 3040 3041 APInt 3042 APFloat::convertFloatAPFloatToAPInt() const 3043 { 3044 assert(semantics == (const llvm::fltSemantics*)&IEEEsingle); 3045 assert(partCount()==1); 3046 3047 uint32_t myexponent, mysignificand; 3048 3049 if (isFiniteNonZero()) { 3050 myexponent = exponent+127; //bias 3051 mysignificand = (uint32_t)*significandParts(); 3052 if (myexponent == 1 && !(mysignificand & 0x800000)) 3053 myexponent = 0; // denormal 3054 } else if (category==fcZero) { 3055 myexponent = 0; 3056 mysignificand = 0; 3057 } else if (category==fcInfinity) { 3058 myexponent = 0xff; 3059 mysignificand = 0; 3060 } else { 3061 assert(category == fcNaN && "Unknown category!"); 3062 myexponent = 0xff; 3063 mysignificand = (uint32_t)*significandParts(); 3064 } 3065 3066 return APInt(32, (((sign&1) << 31) | ((myexponent&0xff) << 23) | 3067 (mysignificand & 0x7fffff))); 3068 } 3069 3070 APInt 3071 APFloat::convertHalfAPFloatToAPInt() const 3072 { 3073 assert(semantics == (const llvm::fltSemantics*)&IEEEhalf); 3074 assert(partCount()==1); 3075 3076 uint32_t myexponent, mysignificand; 3077 3078 if (isFiniteNonZero()) { 3079 myexponent = exponent+15; //bias 3080 mysignificand = (uint32_t)*significandParts(); 3081 if (myexponent == 1 && !(mysignificand & 0x400)) 3082 myexponent = 0; // denormal 3083 } else if (category==fcZero) { 3084 myexponent = 0; 3085 mysignificand = 0; 3086 } else if (category==fcInfinity) { 3087 myexponent = 0x1f; 3088 mysignificand = 0; 3089 } else { 3090 assert(category == fcNaN && "Unknown category!"); 3091 myexponent = 0x1f; 3092 mysignificand = (uint32_t)*significandParts(); 3093 } 3094 3095 return APInt(16, (((sign&1) << 15) | ((myexponent&0x1f) << 10) | 3096 (mysignificand & 0x3ff))); 3097 } 3098 3099 // This function creates an APInt that is just a bit map of the floating 3100 // point constant as it would appear in memory. It is not a conversion, 3101 // and treating the result as a normal integer is unlikely to be useful. 3102 3103 APInt 3104 APFloat::bitcastToAPInt() const 3105 { 3106 if (semantics == (const llvm::fltSemantics*)&IEEEhalf) 3107 return convertHalfAPFloatToAPInt(); 3108 3109 if (semantics == (const llvm::fltSemantics*)&IEEEsingle) 3110 return convertFloatAPFloatToAPInt(); 3111 3112 if (semantics == (const llvm::fltSemantics*)&IEEEdouble) 3113 return convertDoubleAPFloatToAPInt(); 3114 3115 if (semantics == (const llvm::fltSemantics*)&IEEEquad) 3116 return convertQuadrupleAPFloatToAPInt(); 3117 3118 if (semantics == (const llvm::fltSemantics*)&PPCDoubleDouble) 3119 return convertPPCDoubleDoubleAPFloatToAPInt(); 3120 3121 assert(semantics == (const llvm::fltSemantics*)&x87DoubleExtended && 3122 "unknown format!"); 3123 return convertF80LongDoubleAPFloatToAPInt(); 3124 } 3125 3126 float 3127 APFloat::convertToFloat() const 3128 { 3129 assert(semantics == (const llvm::fltSemantics*)&IEEEsingle && 3130 "Float semantics are not IEEEsingle"); 3131 APInt api = bitcastToAPInt(); 3132 return api.bitsToFloat(); 3133 } 3134 3135 double 3136 APFloat::convertToDouble() const 3137 { 3138 assert(semantics == (const llvm::fltSemantics*)&IEEEdouble && 3139 "Float semantics are not IEEEdouble"); 3140 APInt api = bitcastToAPInt(); 3141 return api.bitsToDouble(); 3142 } 3143 3144 /// Integer bit is explicit in this format. Intel hardware (387 and later) 3145 /// does not support these bit patterns: 3146 /// exponent = all 1's, integer bit 0, significand 0 ("pseudoinfinity") 3147 /// exponent = all 1's, integer bit 0, significand nonzero ("pseudoNaN") 3148 /// exponent = 0, integer bit 1 ("pseudodenormal") 3149 /// exponent!=0 nor all 1's, integer bit 0 ("unnormal") 3150 /// At the moment, the first two are treated as NaNs, the second two as Normal. 3151 void 3152 APFloat::initFromF80LongDoubleAPInt(const APInt &api) 3153 { 3154 assert(api.getBitWidth()==80); 3155 uint64_t i1 = api.getRawData()[0]; 3156 uint64_t i2 = api.getRawData()[1]; 3157 uint64_t myexponent = (i2 & 0x7fff); 3158 uint64_t mysignificand = i1; 3159 3160 initialize(&APFloat::x87DoubleExtended); 3161 assert(partCount()==2); 3162 3163 sign = static_cast<unsigned int>(i2>>15); 3164 if (myexponent==0 && mysignificand==0) { 3165 // exponent, significand meaningless 3166 category = fcZero; 3167 } else if (myexponent==0x7fff && mysignificand==0x8000000000000000ULL) { 3168 // exponent, significand meaningless 3169 category = fcInfinity; 3170 } else if (myexponent==0x7fff && mysignificand!=0x8000000000000000ULL) { 3171 // exponent meaningless 3172 category = fcNaN; 3173 significandParts()[0] = mysignificand; 3174 significandParts()[1] = 0; 3175 } else { 3176 category = fcNormal; 3177 exponent = myexponent - 16383; 3178 significandParts()[0] = mysignificand; 3179 significandParts()[1] = 0; 3180 if (myexponent==0) // denormal 3181 exponent = -16382; 3182 } 3183 } 3184 3185 void 3186 APFloat::initFromPPCDoubleDoubleAPInt(const APInt &api) 3187 { 3188 assert(api.getBitWidth()==128); 3189 uint64_t i1 = api.getRawData()[0]; 3190 uint64_t i2 = api.getRawData()[1]; 3191 opStatus fs; 3192 bool losesInfo; 3193 3194 // Get the first double and convert to our format. 3195 initFromDoubleAPInt(APInt(64, i1)); 3196 fs = convert(PPCDoubleDouble, rmNearestTiesToEven, &losesInfo); 3197 assert(fs == opOK && !losesInfo); 3198 (void)fs; 3199 3200 // Unless we have a special case, add in second double. 3201 if (isFiniteNonZero()) { 3202 APFloat v(IEEEdouble, APInt(64, i2)); 3203 fs = v.convert(PPCDoubleDouble, rmNearestTiesToEven, &losesInfo); 3204 assert(fs == opOK && !losesInfo); 3205 (void)fs; 3206 3207 add(v, rmNearestTiesToEven); 3208 } 3209 } 3210 3211 void 3212 APFloat::initFromQuadrupleAPInt(const APInt &api) 3213 { 3214 assert(api.getBitWidth()==128); 3215 uint64_t i1 = api.getRawData()[0]; 3216 uint64_t i2 = api.getRawData()[1]; 3217 uint64_t myexponent = (i2 >> 48) & 0x7fff; 3218 uint64_t mysignificand = i1; 3219 uint64_t mysignificand2 = i2 & 0xffffffffffffLL; 3220 3221 initialize(&APFloat::IEEEquad); 3222 assert(partCount()==2); 3223 3224 sign = static_cast<unsigned int>(i2>>63); 3225 if (myexponent==0 && 3226 (mysignificand==0 && mysignificand2==0)) { 3227 // exponent, significand meaningless 3228 category = fcZero; 3229 } else if (myexponent==0x7fff && 3230 (mysignificand==0 && mysignificand2==0)) { 3231 // exponent, significand meaningless 3232 category = fcInfinity; 3233 } else if (myexponent==0x7fff && 3234 (mysignificand!=0 || mysignificand2 !=0)) { 3235 // exponent meaningless 3236 category = fcNaN; 3237 significandParts()[0] = mysignificand; 3238 significandParts()[1] = mysignificand2; 3239 } else { 3240 category = fcNormal; 3241 exponent = myexponent - 16383; 3242 significandParts()[0] = mysignificand; 3243 significandParts()[1] = mysignificand2; 3244 if (myexponent==0) // denormal 3245 exponent = -16382; 3246 else 3247 significandParts()[1] |= 0x1000000000000LL; // integer bit 3248 } 3249 } 3250 3251 void 3252 APFloat::initFromDoubleAPInt(const APInt &api) 3253 { 3254 assert(api.getBitWidth()==64); 3255 uint64_t i = *api.getRawData(); 3256 uint64_t myexponent = (i >> 52) & 0x7ff; 3257 uint64_t mysignificand = i & 0xfffffffffffffLL; 3258 3259 initialize(&APFloat::IEEEdouble); 3260 assert(partCount()==1); 3261 3262 sign = static_cast<unsigned int>(i>>63); 3263 if (myexponent==0 && mysignificand==0) { 3264 // exponent, significand meaningless 3265 category = fcZero; 3266 } else if (myexponent==0x7ff && mysignificand==0) { 3267 // exponent, significand meaningless 3268 category = fcInfinity; 3269 } else if (myexponent==0x7ff && mysignificand!=0) { 3270 // exponent meaningless 3271 category = fcNaN; 3272 *significandParts() = mysignificand; 3273 } else { 3274 category = fcNormal; 3275 exponent = myexponent - 1023; 3276 *significandParts() = mysignificand; 3277 if (myexponent==0) // denormal 3278 exponent = -1022; 3279 else 3280 *significandParts() |= 0x10000000000000LL; // integer bit 3281 } 3282 } 3283 3284 void 3285 APFloat::initFromFloatAPInt(const APInt & api) 3286 { 3287 assert(api.getBitWidth()==32); 3288 uint32_t i = (uint32_t)*api.getRawData(); 3289 uint32_t myexponent = (i >> 23) & 0xff; 3290 uint32_t mysignificand = i & 0x7fffff; 3291 3292 initialize(&APFloat::IEEEsingle); 3293 assert(partCount()==1); 3294 3295 sign = i >> 31; 3296 if (myexponent==0 && mysignificand==0) { 3297 // exponent, significand meaningless 3298 category = fcZero; 3299 } else if (myexponent==0xff && mysignificand==0) { 3300 // exponent, significand meaningless 3301 category = fcInfinity; 3302 } else if (myexponent==0xff && mysignificand!=0) { 3303 // sign, exponent, significand meaningless 3304 category = fcNaN; 3305 *significandParts() = mysignificand; 3306 } else { 3307 category = fcNormal; 3308 exponent = myexponent - 127; //bias 3309 *significandParts() = mysignificand; 3310 if (myexponent==0) // denormal 3311 exponent = -126; 3312 else 3313 *significandParts() |= 0x800000; // integer bit 3314 } 3315 } 3316 3317 void 3318 APFloat::initFromHalfAPInt(const APInt & api) 3319 { 3320 assert(api.getBitWidth()==16); 3321 uint32_t i = (uint32_t)*api.getRawData(); 3322 uint32_t myexponent = (i >> 10) & 0x1f; 3323 uint32_t mysignificand = i & 0x3ff; 3324 3325 initialize(&APFloat::IEEEhalf); 3326 assert(partCount()==1); 3327 3328 sign = i >> 15; 3329 if (myexponent==0 && mysignificand==0) { 3330 // exponent, significand meaningless 3331 category = fcZero; 3332 } else if (myexponent==0x1f && mysignificand==0) { 3333 // exponent, significand meaningless 3334 category = fcInfinity; 3335 } else if (myexponent==0x1f && mysignificand!=0) { 3336 // sign, exponent, significand meaningless 3337 category = fcNaN; 3338 *significandParts() = mysignificand; 3339 } else { 3340 category = fcNormal; 3341 exponent = myexponent - 15; //bias 3342 *significandParts() = mysignificand; 3343 if (myexponent==0) // denormal 3344 exponent = -14; 3345 else 3346 *significandParts() |= 0x400; // integer bit 3347 } 3348 } 3349 3350 /// Treat api as containing the bits of a floating point number. Currently 3351 /// we infer the floating point type from the size of the APInt. The 3352 /// isIEEE argument distinguishes between PPC128 and IEEE128 (not meaningful 3353 /// when the size is anything else). 3354 void 3355 APFloat::initFromAPInt(const fltSemantics* Sem, const APInt& api) 3356 { 3357 if (Sem == &IEEEhalf) 3358 return initFromHalfAPInt(api); 3359 if (Sem == &IEEEsingle) 3360 return initFromFloatAPInt(api); 3361 if (Sem == &IEEEdouble) 3362 return initFromDoubleAPInt(api); 3363 if (Sem == &x87DoubleExtended) 3364 return initFromF80LongDoubleAPInt(api); 3365 if (Sem == &IEEEquad) 3366 return initFromQuadrupleAPInt(api); 3367 if (Sem == &PPCDoubleDouble) 3368 return initFromPPCDoubleDoubleAPInt(api); 3369 3370 llvm_unreachable(nullptr); 3371 } 3372 3373 APFloat 3374 APFloat::getAllOnesValue(unsigned BitWidth, bool isIEEE) 3375 { 3376 switch (BitWidth) { 3377 case 16: 3378 return APFloat(IEEEhalf, APInt::getAllOnesValue(BitWidth)); 3379 case 32: 3380 return APFloat(IEEEsingle, APInt::getAllOnesValue(BitWidth)); 3381 case 64: 3382 return APFloat(IEEEdouble, APInt::getAllOnesValue(BitWidth)); 3383 case 80: 3384 return APFloat(x87DoubleExtended, APInt::getAllOnesValue(BitWidth)); 3385 case 128: 3386 if (isIEEE) 3387 return APFloat(IEEEquad, APInt::getAllOnesValue(BitWidth)); 3388 return APFloat(PPCDoubleDouble, APInt::getAllOnesValue(BitWidth)); 3389 default: 3390 llvm_unreachable("Unknown floating bit width"); 3391 } 3392 } 3393 3394 unsigned APFloat::getSizeInBits(const fltSemantics &Sem) { 3395 return Sem.sizeInBits; 3396 } 3397 3398 /// Make this number the largest magnitude normal number in the given 3399 /// semantics. 3400 void APFloat::makeLargest(bool Negative) { 3401 // We want (in interchange format): 3402 // sign = {Negative} 3403 // exponent = 1..10 3404 // significand = 1..1 3405 category = fcNormal; 3406 sign = Negative; 3407 exponent = semantics->maxExponent; 3408 3409 // Use memset to set all but the highest integerPart to all ones. 3410 integerPart *significand = significandParts(); 3411 unsigned PartCount = partCount(); 3412 memset(significand, 0xFF, sizeof(integerPart)*(PartCount - 1)); 3413 3414 // Set the high integerPart especially setting all unused top bits for 3415 // internal consistency. 3416 const unsigned NumUnusedHighBits = 3417 PartCount*integerPartWidth - semantics->precision; 3418 significand[PartCount - 1] = (NumUnusedHighBits < integerPartWidth) 3419 ? (~integerPart(0) >> NumUnusedHighBits) 3420 : 0; 3421 } 3422 3423 /// Make this number the smallest magnitude denormal number in the given 3424 /// semantics. 3425 void APFloat::makeSmallest(bool Negative) { 3426 // We want (in interchange format): 3427 // sign = {Negative} 3428 // exponent = 0..0 3429 // significand = 0..01 3430 category = fcNormal; 3431 sign = Negative; 3432 exponent = semantics->minExponent; 3433 APInt::tcSet(significandParts(), 1, partCount()); 3434 } 3435 3436 3437 APFloat APFloat::getLargest(const fltSemantics &Sem, bool Negative) { 3438 // We want (in interchange format): 3439 // sign = {Negative} 3440 // exponent = 1..10 3441 // significand = 1..1 3442 APFloat Val(Sem, uninitialized); 3443 Val.makeLargest(Negative); 3444 return Val; 3445 } 3446 3447 APFloat APFloat::getSmallest(const fltSemantics &Sem, bool Negative) { 3448 // We want (in interchange format): 3449 // sign = {Negative} 3450 // exponent = 0..0 3451 // significand = 0..01 3452 APFloat Val(Sem, uninitialized); 3453 Val.makeSmallest(Negative); 3454 return Val; 3455 } 3456 3457 APFloat APFloat::getSmallestNormalized(const fltSemantics &Sem, bool Negative) { 3458 APFloat Val(Sem, uninitialized); 3459 3460 // We want (in interchange format): 3461 // sign = {Negative} 3462 // exponent = 0..0 3463 // significand = 10..0 3464 3465 Val.category = fcNormal; 3466 Val.zeroSignificand(); 3467 Val.sign = Negative; 3468 Val.exponent = Sem.minExponent; 3469 Val.significandParts()[partCountForBits(Sem.precision)-1] |= 3470 (((integerPart) 1) << ((Sem.precision - 1) % integerPartWidth)); 3471 3472 return Val; 3473 } 3474 3475 APFloat::APFloat(const fltSemantics &Sem, const APInt &API) { 3476 initFromAPInt(&Sem, API); 3477 } 3478 3479 APFloat::APFloat(float f) { 3480 initFromAPInt(&IEEEsingle, APInt::floatToBits(f)); 3481 } 3482 3483 APFloat::APFloat(double d) { 3484 initFromAPInt(&IEEEdouble, APInt::doubleToBits(d)); 3485 } 3486 3487 namespace { 3488 void append(SmallVectorImpl<char> &Buffer, StringRef Str) { 3489 Buffer.append(Str.begin(), Str.end()); 3490 } 3491 3492 /// Removes data from the given significand until it is no more 3493 /// precise than is required for the desired precision. 3494 void AdjustToPrecision(APInt &significand, 3495 int &exp, unsigned FormatPrecision) { 3496 unsigned bits = significand.getActiveBits(); 3497 3498 // 196/59 is a very slight overestimate of lg_2(10). 3499 unsigned bitsRequired = (FormatPrecision * 196 + 58) / 59; 3500 3501 if (bits <= bitsRequired) return; 3502 3503 unsigned tensRemovable = (bits - bitsRequired) * 59 / 196; 3504 if (!tensRemovable) return; 3505 3506 exp += tensRemovable; 3507 3508 APInt divisor(significand.getBitWidth(), 1); 3509 APInt powten(significand.getBitWidth(), 10); 3510 while (true) { 3511 if (tensRemovable & 1) 3512 divisor *= powten; 3513 tensRemovable >>= 1; 3514 if (!tensRemovable) break; 3515 powten *= powten; 3516 } 3517 3518 significand = significand.udiv(divisor); 3519 3520 // Truncate the significand down to its active bit count. 3521 significand = significand.trunc(significand.getActiveBits()); 3522 } 3523 3524 3525 void AdjustToPrecision(SmallVectorImpl<char> &buffer, 3526 int &exp, unsigned FormatPrecision) { 3527 unsigned N = buffer.size(); 3528 if (N <= FormatPrecision) return; 3529 3530 // The most significant figures are the last ones in the buffer. 3531 unsigned FirstSignificant = N - FormatPrecision; 3532 3533 // Round. 3534 // FIXME: this probably shouldn't use 'round half up'. 3535 3536 // Rounding down is just a truncation, except we also want to drop 3537 // trailing zeros from the new result. 3538 if (buffer[FirstSignificant - 1] < '5') { 3539 while (FirstSignificant < N && buffer[FirstSignificant] == '0') 3540 FirstSignificant++; 3541 3542 exp += FirstSignificant; 3543 buffer.erase(&buffer[0], &buffer[FirstSignificant]); 3544 return; 3545 } 3546 3547 // Rounding up requires a decimal add-with-carry. If we continue 3548 // the carry, the newly-introduced zeros will just be truncated. 3549 for (unsigned I = FirstSignificant; I != N; ++I) { 3550 if (buffer[I] == '9') { 3551 FirstSignificant++; 3552 } else { 3553 buffer[I]++; 3554 break; 3555 } 3556 } 3557 3558 // If we carried through, we have exactly one digit of precision. 3559 if (FirstSignificant == N) { 3560 exp += FirstSignificant; 3561 buffer.clear(); 3562 buffer.push_back('1'); 3563 return; 3564 } 3565 3566 exp += FirstSignificant; 3567 buffer.erase(&buffer[0], &buffer[FirstSignificant]); 3568 } 3569 } 3570 3571 void APFloat::toString(SmallVectorImpl<char> &Str, 3572 unsigned FormatPrecision, 3573 unsigned FormatMaxPadding) const { 3574 switch (category) { 3575 case fcInfinity: 3576 if (isNegative()) 3577 return append(Str, "-Inf"); 3578 else 3579 return append(Str, "+Inf"); 3580 3581 case fcNaN: return append(Str, "NaN"); 3582 3583 case fcZero: 3584 if (isNegative()) 3585 Str.push_back('-'); 3586 3587 if (!FormatMaxPadding) 3588 append(Str, "0.0E+0"); 3589 else 3590 Str.push_back('0'); 3591 return; 3592 3593 case fcNormal: 3594 break; 3595 } 3596 3597 if (isNegative()) 3598 Str.push_back('-'); 3599 3600 // Decompose the number into an APInt and an exponent. 3601 int exp = exponent - ((int) semantics->precision - 1); 3602 APInt significand(semantics->precision, 3603 makeArrayRef(significandParts(), 3604 partCountForBits(semantics->precision))); 3605 3606 // Set FormatPrecision if zero. We want to do this before we 3607 // truncate trailing zeros, as those are part of the precision. 3608 if (!FormatPrecision) { 3609 // We use enough digits so the number can be round-tripped back to an 3610 // APFloat. The formula comes from "How to Print Floating-Point Numbers 3611 // Accurately" by Steele and White. 3612 // FIXME: Using a formula based purely on the precision is conservative; 3613 // we can print fewer digits depending on the actual value being printed. 3614 3615 // FormatPrecision = 2 + floor(significandBits / lg_2(10)) 3616 FormatPrecision = 2 + semantics->precision * 59 / 196; 3617 } 3618 3619 // Ignore trailing binary zeros. 3620 int trailingZeros = significand.countTrailingZeros(); 3621 exp += trailingZeros; 3622 significand = significand.lshr(trailingZeros); 3623 3624 // Change the exponent from 2^e to 10^e. 3625 if (exp == 0) { 3626 // Nothing to do. 3627 } else if (exp > 0) { 3628 // Just shift left. 3629 significand = significand.zext(semantics->precision + exp); 3630 significand <<= exp; 3631 exp = 0; 3632 } else { /* exp < 0 */ 3633 int texp = -exp; 3634 3635 // We transform this using the identity: 3636 // (N)(2^-e) == (N)(5^e)(10^-e) 3637 // This means we have to multiply N (the significand) by 5^e. 3638 // To avoid overflow, we have to operate on numbers large 3639 // enough to store N * 5^e: 3640 // log2(N * 5^e) == log2(N) + e * log2(5) 3641 // <= semantics->precision + e * 137 / 59 3642 // (log_2(5) ~ 2.321928 < 2.322034 ~ 137/59) 3643 3644 unsigned precision = semantics->precision + (137 * texp + 136) / 59; 3645 3646 // Multiply significand by 5^e. 3647 // N * 5^0101 == N * 5^(1*1) * 5^(0*2) * 5^(1*4) * 5^(0*8) 3648 significand = significand.zext(precision); 3649 APInt five_to_the_i(precision, 5); 3650 while (true) { 3651 if (texp & 1) significand *= five_to_the_i; 3652 3653 texp >>= 1; 3654 if (!texp) break; 3655 five_to_the_i *= five_to_the_i; 3656 } 3657 } 3658 3659 AdjustToPrecision(significand, exp, FormatPrecision); 3660 3661 SmallVector<char, 256> buffer; 3662 3663 // Fill the buffer. 3664 unsigned precision = significand.getBitWidth(); 3665 APInt ten(precision, 10); 3666 APInt digit(precision, 0); 3667 3668 bool inTrail = true; 3669 while (significand != 0) { 3670 // digit <- significand % 10 3671 // significand <- significand / 10 3672 APInt::udivrem(significand, ten, significand, digit); 3673 3674 unsigned d = digit.getZExtValue(); 3675 3676 // Drop trailing zeros. 3677 if (inTrail && !d) exp++; 3678 else { 3679 buffer.push_back((char) ('0' + d)); 3680 inTrail = false; 3681 } 3682 } 3683 3684 assert(!buffer.empty() && "no characters in buffer!"); 3685 3686 // Drop down to FormatPrecision. 3687 // TODO: don't do more precise calculations above than are required. 3688 AdjustToPrecision(buffer, exp, FormatPrecision); 3689 3690 unsigned NDigits = buffer.size(); 3691 3692 // Check whether we should use scientific notation. 3693 bool FormatScientific; 3694 if (!FormatMaxPadding) 3695 FormatScientific = true; 3696 else { 3697 if (exp >= 0) { 3698 // 765e3 --> 765000 3699 // ^^^ 3700 // But we shouldn't make the number look more precise than it is. 3701 FormatScientific = ((unsigned) exp > FormatMaxPadding || 3702 NDigits + (unsigned) exp > FormatPrecision); 3703 } else { 3704 // Power of the most significant digit. 3705 int MSD = exp + (int) (NDigits - 1); 3706 if (MSD >= 0) { 3707 // 765e-2 == 7.65 3708 FormatScientific = false; 3709 } else { 3710 // 765e-5 == 0.00765 3711 // ^ ^^ 3712 FormatScientific = ((unsigned) -MSD) > FormatMaxPadding; 3713 } 3714 } 3715 } 3716 3717 // Scientific formatting is pretty straightforward. 3718 if (FormatScientific) { 3719 exp += (NDigits - 1); 3720 3721 Str.push_back(buffer[NDigits-1]); 3722 Str.push_back('.'); 3723 if (NDigits == 1) 3724 Str.push_back('0'); 3725 else 3726 for (unsigned I = 1; I != NDigits; ++I) 3727 Str.push_back(buffer[NDigits-1-I]); 3728 Str.push_back('E'); 3729 3730 Str.push_back(exp >= 0 ? '+' : '-'); 3731 if (exp < 0) exp = -exp; 3732 SmallVector<char, 6> expbuf; 3733 do { 3734 expbuf.push_back((char) ('0' + (exp % 10))); 3735 exp /= 10; 3736 } while (exp); 3737 for (unsigned I = 0, E = expbuf.size(); I != E; ++I) 3738 Str.push_back(expbuf[E-1-I]); 3739 return; 3740 } 3741 3742 // Non-scientific, positive exponents. 3743 if (exp >= 0) { 3744 for (unsigned I = 0; I != NDigits; ++I) 3745 Str.push_back(buffer[NDigits-1-I]); 3746 for (unsigned I = 0; I != (unsigned) exp; ++I) 3747 Str.push_back('0'); 3748 return; 3749 } 3750 3751 // Non-scientific, negative exponents. 3752 3753 // The number of digits to the left of the decimal point. 3754 int NWholeDigits = exp + (int) NDigits; 3755 3756 unsigned I = 0; 3757 if (NWholeDigits > 0) { 3758 for (; I != (unsigned) NWholeDigits; ++I) 3759 Str.push_back(buffer[NDigits-I-1]); 3760 Str.push_back('.'); 3761 } else { 3762 unsigned NZeros = 1 + (unsigned) -NWholeDigits; 3763 3764 Str.push_back('0'); 3765 Str.push_back('.'); 3766 for (unsigned Z = 1; Z != NZeros; ++Z) 3767 Str.push_back('0'); 3768 } 3769 3770 for (; I != NDigits; ++I) 3771 Str.push_back(buffer[NDigits-I-1]); 3772 } 3773 3774 bool APFloat::getExactInverse(APFloat *inv) const { 3775 // Special floats and denormals have no exact inverse. 3776 if (!isFiniteNonZero()) 3777 return false; 3778 3779 // Check that the number is a power of two by making sure that only the 3780 // integer bit is set in the significand. 3781 if (significandLSB() != semantics->precision - 1) 3782 return false; 3783 3784 // Get the inverse. 3785 APFloat reciprocal(*semantics, 1ULL); 3786 if (reciprocal.divide(*this, rmNearestTiesToEven) != opOK) 3787 return false; 3788 3789 // Avoid multiplication with a denormal, it is not safe on all platforms and 3790 // may be slower than a normal division. 3791 if (reciprocal.isDenormal()) 3792 return false; 3793 3794 assert(reciprocal.isFiniteNonZero() && 3795 reciprocal.significandLSB() == reciprocal.semantics->precision - 1); 3796 3797 if (inv) 3798 *inv = reciprocal; 3799 3800 return true; 3801 } 3802 3803 bool APFloat::isSignaling() const { 3804 if (!isNaN()) 3805 return false; 3806 3807 // IEEE-754R 2008 6.2.1: A signaling NaN bit string should be encoded with the 3808 // first bit of the trailing significand being 0. 3809 return !APInt::tcExtractBit(significandParts(), semantics->precision - 2); 3810 } 3811 3812 /// IEEE-754R 2008 5.3.1: nextUp/nextDown. 3813 /// 3814 /// *NOTE* since nextDown(x) = -nextUp(-x), we only implement nextUp with 3815 /// appropriate sign switching before/after the computation. 3816 APFloat::opStatus APFloat::next(bool nextDown) { 3817 // If we are performing nextDown, swap sign so we have -x. 3818 if (nextDown) 3819 changeSign(); 3820 3821 // Compute nextUp(x) 3822 opStatus result = opOK; 3823 3824 // Handle each float category separately. 3825 switch (category) { 3826 case fcInfinity: 3827 // nextUp(+inf) = +inf 3828 if (!isNegative()) 3829 break; 3830 // nextUp(-inf) = -getLargest() 3831 makeLargest(true); 3832 break; 3833 case fcNaN: 3834 // IEEE-754R 2008 6.2 Par 2: nextUp(sNaN) = qNaN. Set Invalid flag. 3835 // IEEE-754R 2008 6.2: nextUp(qNaN) = qNaN. Must be identity so we do not 3836 // change the payload. 3837 if (isSignaling()) { 3838 result = opInvalidOp; 3839 // For consistency, propagate the sign of the sNaN to the qNaN. 3840 makeNaN(false, isNegative(), nullptr); 3841 } 3842 break; 3843 case fcZero: 3844 // nextUp(pm 0) = +getSmallest() 3845 makeSmallest(false); 3846 break; 3847 case fcNormal: 3848 // nextUp(-getSmallest()) = -0 3849 if (isSmallest() && isNegative()) { 3850 APInt::tcSet(significandParts(), 0, partCount()); 3851 category = fcZero; 3852 exponent = 0; 3853 break; 3854 } 3855 3856 // nextUp(getLargest()) == INFINITY 3857 if (isLargest() && !isNegative()) { 3858 APInt::tcSet(significandParts(), 0, partCount()); 3859 category = fcInfinity; 3860 exponent = semantics->maxExponent + 1; 3861 break; 3862 } 3863 3864 // nextUp(normal) == normal + inc. 3865 if (isNegative()) { 3866 // If we are negative, we need to decrement the significand. 3867 3868 // We only cross a binade boundary that requires adjusting the exponent 3869 // if: 3870 // 1. exponent != semantics->minExponent. This implies we are not in the 3871 // smallest binade or are dealing with denormals. 3872 // 2. Our significand excluding the integral bit is all zeros. 3873 bool WillCrossBinadeBoundary = 3874 exponent != semantics->minExponent && isSignificandAllZeros(); 3875 3876 // Decrement the significand. 3877 // 3878 // We always do this since: 3879 // 1. If we are dealing with a non-binade decrement, by definition we 3880 // just decrement the significand. 3881 // 2. If we are dealing with a normal -> normal binade decrement, since 3882 // we have an explicit integral bit the fact that all bits but the 3883 // integral bit are zero implies that subtracting one will yield a 3884 // significand with 0 integral bit and 1 in all other spots. Thus we 3885 // must just adjust the exponent and set the integral bit to 1. 3886 // 3. If we are dealing with a normal -> denormal binade decrement, 3887 // since we set the integral bit to 0 when we represent denormals, we 3888 // just decrement the significand. 3889 integerPart *Parts = significandParts(); 3890 APInt::tcDecrement(Parts, partCount()); 3891 3892 if (WillCrossBinadeBoundary) { 3893 // Our result is a normal number. Do the following: 3894 // 1. Set the integral bit to 1. 3895 // 2. Decrement the exponent. 3896 APInt::tcSetBit(Parts, semantics->precision - 1); 3897 exponent--; 3898 } 3899 } else { 3900 // If we are positive, we need to increment the significand. 3901 3902 // We only cross a binade boundary that requires adjusting the exponent if 3903 // the input is not a denormal and all of said input's significand bits 3904 // are set. If all of said conditions are true: clear the significand, set 3905 // the integral bit to 1, and increment the exponent. If we have a 3906 // denormal always increment since moving denormals and the numbers in the 3907 // smallest normal binade have the same exponent in our representation. 3908 bool WillCrossBinadeBoundary = !isDenormal() && isSignificandAllOnes(); 3909 3910 if (WillCrossBinadeBoundary) { 3911 integerPart *Parts = significandParts(); 3912 APInt::tcSet(Parts, 0, partCount()); 3913 APInt::tcSetBit(Parts, semantics->precision - 1); 3914 assert(exponent != semantics->maxExponent && 3915 "We can not increment an exponent beyond the maxExponent allowed" 3916 " by the given floating point semantics."); 3917 exponent++; 3918 } else { 3919 incrementSignificand(); 3920 } 3921 } 3922 break; 3923 } 3924 3925 // If we are performing nextDown, swap sign so we have -nextUp(-x) 3926 if (nextDown) 3927 changeSign(); 3928 3929 return result; 3930 } 3931 3932 void 3933 APFloat::makeInf(bool Negative) { 3934 category = fcInfinity; 3935 sign = Negative; 3936 exponent = semantics->maxExponent + 1; 3937 APInt::tcSet(significandParts(), 0, partCount()); 3938 } 3939 3940 void 3941 APFloat::makeZero(bool Negative) { 3942 category = fcZero; 3943 sign = Negative; 3944 exponent = semantics->minExponent-1; 3945 APInt::tcSet(significandParts(), 0, partCount()); 3946 } 3947 3948 APFloat llvm::scalbn(APFloat X, int Exp) { 3949 if (X.isInfinity() || X.isZero() || X.isNaN()) 3950 return X; 3951 3952 auto MaxExp = X.getSemantics().maxExponent; 3953 auto MinExp = X.getSemantics().minExponent; 3954 if (Exp > (MaxExp - X.exponent)) 3955 // Overflow saturates to infinity. 3956 return APFloat::getInf(X.getSemantics(), X.isNegative()); 3957 if (Exp < (MinExp - X.exponent)) 3958 // Underflow saturates to zero. 3959 return APFloat::getZero(X.getSemantics(), X.isNegative()); 3960 3961 X.exponent += Exp; 3962 return X; 3963 } 3964