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