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