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