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