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