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