lib/Support/APFloat.cpp

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