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, 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       APInt::tcShiftLeft(fullSignificand, newPartsCount,
 910                          extendedPrecision - omsb);
 911       exponent -= extendedPrecision - omsb;
 912     }
 913
 914     /* Create new semantics.  */
 915     extendedSemantics = *semantics;
 916     extendedSemantics.precision = extendedPrecision;
 917
 918     if (newPartsCount == 1)
 919       significand.part = fullSignificand[0];
 920     else
 921       significand.parts = fullSignificand;
 922     semantics = &extendedSemantics;
 923
 924     APFloat extendedAddend(*addend);
 925     status = extendedAddend.convert(extendedSemantics, rmTowardZero, &ignored);
 926     assert(status == opOK);
 927     lost_fraction = addOrSubtractSignificand(extendedAddend, false);
 928
 929     /* Restore our state.  */
 930     if (newPartsCount == 1)
 931       fullSignificand[0] = significand.part;
 932     significand = savedSignificand;
 933     semantics = savedSemantics;
 934
 935     omsb = APInt::tcMSB(fullSignificand, newPartsCount) + 1;
 936   }
 937
 938   exponent -= (precision - 1);
 939
 940   if (omsb > precision) {
 941     unsigned int bits, significantParts;
 942     lostFraction lf;
 943
 944     bits = omsb - precision;
 945     significantParts = partCountForBits(omsb);
 946     lf = shiftRight(fullSignificand, significantParts, bits);
 947     lost_fraction = combineLostFractions(lf, lost_fraction);
 948     exponent += bits;
 949   }
 950
 951   APInt::tcAssign(lhsSignificand, fullSignificand, partsCount);
 952
 953   if (newPartsCount > 4)
 954     delete [] fullSignificand;
 955
 956   return lost_fraction;
 957 }
 958
 959 /* Multiply the significands of LHS and RHS to DST.  */
 960 lostFraction
 961 APFloat::divideSignificand(const APFloat &rhs)
 962 {
 963   unsigned int bit, i, partsCount;
 964   const integerPart *rhsSignificand;
 965   integerPart *lhsSignificand, *dividend, *divisor;
 966   integerPart scratch[4];
 967   lostFraction lost_fraction;
 968
 969   assert(semantics == rhs.semantics);
 970
 971   lhsSignificand = significandParts();
 972   rhsSignificand = rhs.significandParts();
 973   partsCount = partCount();
 974
 975   if (partsCount > 2)
 976     dividend = new integerPart[partsCount * 2];
 977   else
 978     dividend = scratch;
 979
 980   divisor = dividend + partsCount;
 981
 982   /* Copy the dividend and divisor as they will be modified in-place.  */
 983   for (i = 0; i < partsCount; i++) {
 984     dividend[i] = lhsSignificand[i];
 985     divisor[i] = rhsSignificand[i];
 986     lhsSignificand[i] = 0;
 987   }
 988
 989   exponent -= rhs.exponent;
 990
 991   unsigned int precision = semantics->precision;
 992
 993   /* Normalize the divisor.  */
 994   bit = precision - APInt::tcMSB(divisor, partsCount) - 1;
 995   if (bit) {
 996     exponent += bit;
 997     APInt::tcShiftLeft(divisor, partsCount, bit);
 998   }
 999
1000   /* Normalize the dividend.  */
1001   bit = precision - APInt::tcMSB(dividend, partsCount) - 1;
1002   if (bit) {
1003     exponent -= bit;
1004     APInt::tcShiftLeft(dividend, partsCount, bit);
1005   }
1006
1007   /* Ensure the dividend >= divisor initially for the loop below.
1008      Incidentally, this means that the division loop below is
1009      guaranteed to set the integer bit to one.  */
1010   if (APInt::tcCompare(dividend, divisor, partsCount) < 0) {
1011     exponent--;
1012     APInt::tcShiftLeft(dividend, partsCount, 1);
1013     assert(APInt::tcCompare(dividend, divisor, partsCount) >= 0);
1014   }
1015
1016   /* Long division.  */
1017   for (bit = precision; bit; bit -= 1) {
1018     if (APInt::tcCompare(dividend, divisor, partsCount) >= 0) {
1019       APInt::tcSubtract(dividend, divisor, 0, partsCount);
1020       APInt::tcSetBit(lhsSignificand, bit - 1);
1021     }
1022
1023     APInt::tcShiftLeft(dividend, partsCount, 1);
1024   }
1025
1026   /* Figure out the lost fraction.  */
1027   int cmp = APInt::tcCompare(dividend, divisor, partsCount);
1028
1029   if (cmp > 0)
1030     lost_fraction = lfMoreThanHalf;
1031   else if (cmp == 0)
1032     lost_fraction = lfExactlyHalf;
1033   else if (APInt::tcIsZero(dividend, partsCount))
1034     lost_fraction = lfExactlyZero;
1035   else
1036     lost_fraction = lfLessThanHalf;
1037
1038   if (partsCount > 2)
1039     delete [] dividend;
1040
1041   return lost_fraction;
1042 }
1043
1044 unsigned int
1045 APFloat::significandMSB() const
1046 {
1047   return APInt::tcMSB(significandParts(), partCount());
1048 }
1049
1050 unsigned int
1051 APFloat::significandLSB() const
1052 {
1053   return APInt::tcLSB(significandParts(), partCount());
1054 }
1055
1056 /* Note that a zero result is NOT normalized to fcZero.  */
1057 lostFraction
1058 APFloat::shiftSignificandRight(unsigned int bits)
1059 {
1060   /* Our exponent should not overflow.  */
1061   assert((exponent_t) (exponent + bits) >= exponent);
1062
1063   exponent += bits;
1064
1065   return shiftRight(significandParts(), partCount(), bits);
1066 }
1067
1068 /* Shift the significand left BITS bits, subtract BITS from its exponent.  */
1069 void
1070 APFloat::shiftSignificandLeft(unsigned int bits)
1071 {
1072   assert(bits < semantics->precision);
1073
1074   if (bits) {
1075     unsigned int partsCount = partCount();
1076
1077     APInt::tcShiftLeft(significandParts(), partsCount, bits);
1078     exponent -= bits;
1079
1080     assert(!APInt::tcIsZero(significandParts(), partsCount));
1081   }
1082 }
1083
1084 APFloat::cmpResult
1085 APFloat::compareAbsoluteValue(const APFloat &rhs) const
1086 {
1087   int compare;
1088
1089   assert(semantics == rhs.semantics);
1090   assert(category == fcNormal);
1091   assert(rhs.category == fcNormal);
1092
1093   compare = exponent - rhs.exponent;
1094
1095   /* If exponents are equal, do an unsigned bignum comparison of the
1096      significands.  */
1097   if (compare == 0)
1098     compare = APInt::tcCompare(significandParts(), rhs.significandParts(),
1099                                partCount());
1100
1101   if (compare > 0)
1102     return cmpGreaterThan;
1103   else if (compare < 0)
1104     return cmpLessThan;
1105   else
1106     return cmpEqual;
1107 }
1108
1109 /* Handle overflow.  Sign is preserved.  We either become infinity or
1110    the largest finite number.  */
1111 APFloat::opStatus
1112 APFloat::handleOverflow(roundingMode rounding_mode)
1113 {
1114   /* Infinity?  */
1115   if (rounding_mode == rmNearestTiesToEven ||
1116       rounding_mode == rmNearestTiesToAway ||
1117       (rounding_mode == rmTowardPositive && !sign) ||
1118       (rounding_mode == rmTowardNegative && sign)) {
1119     category = fcInfinity;
1120     return (opStatus) (opOverflow | opInexact);
1121   }
1122
1123   /* Otherwise we become the largest finite number.  */
1124   category = fcNormal;
1125   exponent = semantics->maxExponent;
1126   APInt::tcSetLeastSignificantBits(significandParts(), partCount(),
1127                                    semantics->precision);
1128
1129   return opInexact;
1130 }
1131
1132 /* Returns TRUE if, when truncating the current number, with BIT the
1133    new LSB, with the given lost fraction and rounding mode, the result
1134    would need to be rounded away from zero (i.e., by increasing the
1135    signficand).  This routine must work for fcZero of both signs, and
1136    fcNormal numbers.  */
1137 bool
1138 APFloat::roundAwayFromZero(roundingMode rounding_mode,
1139                            lostFraction lost_fraction,
1140                            unsigned int bit) const
1141 {
1142   /* NaNs and infinities should not have lost fractions.  */
1143   assert(category == fcNormal || category == fcZero);
1144
1145   /* Current callers never pass this so we don't handle it.  */
1146   assert(lost_fraction != lfExactlyZero);
1147
1148   switch (rounding_mode) {
1149   default:
1150     llvm_unreachable(0);
1151
1152   case rmNearestTiesToAway:
1153     return lost_fraction == lfExactlyHalf || lost_fraction == lfMoreThanHalf;
1154
1155   case rmNearestTiesToEven:
1156     if (lost_fraction == lfMoreThanHalf)
1157       return true;
1158
1159     /* Our zeroes don't have a significand to test.  */
1160     if (lost_fraction == lfExactlyHalf && category != fcZero)
1161       return APInt::tcExtractBit(significandParts(), bit);
1162
1163     return false;
1164
1165   case rmTowardZero:
1166     return false;
1167
1168   case rmTowardPositive:
1169     return sign == false;
1170
1171   case rmTowardNegative:
1172     return sign == true;
1173   }
1174 }
1175
1176 APFloat::opStatus
1177 APFloat::normalize(roundingMode rounding_mode,
1178                    lostFraction lost_fraction)
1179 {
1180   unsigned int omsb;                /* One, not zero, based MSB.  */
1181   int exponentChange;
1182
1183   if (category != fcNormal)
1184     return opOK;
1185
1186   /* Before rounding normalize the exponent of fcNormal numbers.  */
1187   omsb = significandMSB() + 1;
1188
1189   if (omsb) {
1190     /* OMSB is numbered from 1.  We want to place it in the integer
1191        bit numbered PRECISON if possible, with a compensating change in
1192        the exponent.  */
1193     exponentChange = omsb - semantics->precision;
1194
1195     /* If the resulting exponent is too high, overflow according to
1196        the rounding mode.  */
1197     if (exponent + exponentChange > semantics->maxExponent)
1198       return handleOverflow(rounding_mode);
1199
1200     /* Subnormal numbers have exponent minExponent, and their MSB
1201        is forced based on that.  */
1202     if (exponent + exponentChange < semantics->minExponent)
1203       exponentChange = semantics->minExponent - exponent;
1204
1205     /* Shifting left is easy as we don't lose precision.  */
1206     if (exponentChange < 0) {
1207       assert(lost_fraction == lfExactlyZero);
1208
1209       shiftSignificandLeft(-exponentChange);
1210
1211       return opOK;
1212     }
1213
1214     if (exponentChange > 0) {
1215       lostFraction lf;
1216
1217       /* Shift right and capture any new lost fraction.  */
1218       lf = shiftSignificandRight(exponentChange);
1219
1220       lost_fraction = combineLostFractions(lf, lost_fraction);
1221
1222       /* Keep OMSB up-to-date.  */
1223       if (omsb > (unsigned) exponentChange)
1224         omsb -= exponentChange;
1225       else
1226         omsb = 0;
1227     }
1228   }
1229
1230   /* Now round the number according to rounding_mode given the lost
1231      fraction.  */
1232
1233   /* As specified in IEEE 754, since we do not trap we do not report
1234      underflow for exact results.  */
1235   if (lost_fraction == lfExactlyZero) {
1236     /* Canonicalize zeroes.  */
1237     if (omsb == 0)
1238       category = fcZero;
1239
1240     return opOK;
1241   }
1242
1243   /* Increment the significand if we're rounding away from zero.  */
1244   if (roundAwayFromZero(rounding_mode, lost_fraction, 0)) {
1245     if (omsb == 0)
1246       exponent = semantics->minExponent;
1247
1248     incrementSignificand();
1249     omsb = significandMSB() + 1;
1250
1251     /* Did the significand increment overflow?  */
1252     if (omsb == (unsigned) semantics->precision + 1) {
1253       /* Renormalize by incrementing the exponent and shifting our
1254          significand right one.  However if we already have the
1255          maximum exponent we overflow to infinity.  */
1256       if (exponent == semantics->maxExponent) {
1257         category = fcInfinity;
1258
1259         return (opStatus) (opOverflow | opInexact);
1260       }
1261
1262       shiftSignificandRight(1);
1263
1264       return opInexact;
1265     }
1266   }
1267
1268   /* The normal case - we were and are not denormal, and any
1269      significand increment above didn't overflow.  */
1270   if (omsb == semantics->precision)
1271     return opInexact;
1272
1273   /* We have a non-zero denormal.  */
1274   assert(omsb < semantics->precision);
1275
1276   /* Canonicalize zeroes.  */
1277   if (omsb == 0)
1278     category = fcZero;
1279
1280   /* The fcZero case is a denormal that underflowed to zero.  */
1281   return (opStatus) (opUnderflow | opInexact);
1282 }
1283
1284 APFloat::opStatus
1285 APFloat::addOrSubtractSpecials(const APFloat &rhs, bool subtract)
1286 {
1287   switch (convolve(category, rhs.category)) {
1288   default:
1289     llvm_unreachable(0);
1290
1291   case convolve(fcNaN, fcZero):
1292   case convolve(fcNaN, fcNormal):
1293   case convolve(fcNaN, fcInfinity):
1294   case convolve(fcNaN, fcNaN):
1295   case convolve(fcNormal, fcZero):
1296   case convolve(fcInfinity, fcNormal):
1297   case convolve(fcInfinity, fcZero):
1298     return opOK;
1299
1300   case convolve(fcZero, fcNaN):
1301   case convolve(fcNormal, fcNaN):
1302   case convolve(fcInfinity, fcNaN):
1303     category = fcNaN;
1304     copySignificand(rhs);
1305     return opOK;
1306
1307   case convolve(fcNormal, fcInfinity):
1308   case convolve(fcZero, fcInfinity):
1309     category = fcInfinity;
1310     sign = rhs.sign ^ subtract;
1311     return opOK;
1312
1313   case convolve(fcZero, fcNormal):
1314     assign(rhs);
1315     sign = rhs.sign ^ subtract;
1316     return opOK;
1317
1318   case convolve(fcZero, fcZero):
1319     /* Sign depends on rounding mode; handled by caller.  */
1320     return opOK;
1321
1322   case convolve(fcInfinity, fcInfinity):
1323     /* Differently signed infinities can only be validly
1324        subtracted.  */
1325     if (((sign ^ rhs.sign)!=0) != subtract) {
1326       makeNaN();
1327       return opInvalidOp;
1328     }
1329
1330     return opOK;
1331
1332   case convolve(fcNormal, fcNormal):
1333     return opDivByZero;
1334   }
1335 }
1336
1337 /* Add or subtract two normal numbers.  */
1338 lostFraction
1339 APFloat::addOrSubtractSignificand(const APFloat &rhs, bool subtract)
1340 {
1341   integerPart carry;
1342   lostFraction lost_fraction;
1343   int bits;
1344
1345   /* Determine if the operation on the absolute values is effectively
1346      an addition or subtraction.  */
1347   subtract ^= (sign ^ rhs.sign) ? true : false;
1348
1349   /* Are we bigger exponent-wise than the RHS?  */
1350   bits = exponent - rhs.exponent;
1351
1352   /* Subtraction is more subtle than one might naively expect.  */
1353   if (subtract) {
1354     APFloat temp_rhs(rhs);
1355     bool reverse;
1356
1357     if (bits == 0) {
1358       reverse = compareAbsoluteValue(temp_rhs) == cmpLessThan;
1359       lost_fraction = lfExactlyZero;
1360     } else if (bits > 0) {
1361       lost_fraction = temp_rhs.shiftSignificandRight(bits - 1);
1362       shiftSignificandLeft(1);
1363       reverse = false;
1364     } else {
1365       lost_fraction = shiftSignificandRight(-bits - 1);
1366       temp_rhs.shiftSignificandLeft(1);
1367       reverse = true;
1368     }
1369
1370     if (reverse) {
1371       carry = temp_rhs.subtractSignificand
1372         (*this, lost_fraction != lfExactlyZero);
1373       copySignificand(temp_rhs);
1374       sign = !sign;
1375     } else {
1376       carry = subtractSignificand
1377         (temp_rhs, lost_fraction != lfExactlyZero);
1378     }
1379
1380     /* Invert the lost fraction - it was on the RHS and
1381        subtracted.  */
1382     if (lost_fraction == lfLessThanHalf)
1383       lost_fraction = lfMoreThanHalf;
1384     else if (lost_fraction == lfMoreThanHalf)
1385       lost_fraction = lfLessThanHalf;
1386
1387     /* The code above is intended to ensure that no borrow is
1388        necessary.  */
1389     assert(!carry);
1390   } else {
1391     if (bits > 0) {
1392       APFloat temp_rhs(rhs);
1393
1394       lost_fraction = temp_rhs.shiftSignificandRight(bits);
1395       carry = addSignificand(temp_rhs);
1396     } else {
1397       lost_fraction = shiftSignificandRight(-bits);
1398       carry = addSignificand(rhs);
1399     }
1400
1401     /* We have a guard bit; generating a carry cannot happen.  */
1402     assert(!carry);
1403   }
1404
1405   return lost_fraction;
1406 }
1407
1408 APFloat::opStatus
1409 APFloat::multiplySpecials(const APFloat &rhs)
1410 {
1411   switch (convolve(category, rhs.category)) {
1412   default:
1413     llvm_unreachable(0);
1414
1415   case convolve(fcNaN, fcZero):
1416   case convolve(fcNaN, fcNormal):
1417   case convolve(fcNaN, fcInfinity):
1418   case convolve(fcNaN, fcNaN):
1419     return opOK;
1420
1421   case convolve(fcZero, fcNaN):
1422   case convolve(fcNormal, fcNaN):
1423   case convolve(fcInfinity, fcNaN):
1424     category = fcNaN;
1425     copySignificand(rhs);
1426     return opOK;
1427
1428   case convolve(fcNormal, fcInfinity):
1429   case convolve(fcInfinity, fcNormal):
1430   case convolve(fcInfinity, fcInfinity):
1431     category = fcInfinity;
1432     return opOK;
1433
1434   case convolve(fcZero, fcNormal):
1435   case convolve(fcNormal, fcZero):
1436   case convolve(fcZero, fcZero):
1437     category = fcZero;
1438     return opOK;
1439
1440   case convolve(fcZero, fcInfinity):
1441   case convolve(fcInfinity, fcZero):
1442     makeNaN();
1443     return opInvalidOp;
1444
1445   case convolve(fcNormal, fcNormal):
1446     return opOK;
1447   }
1448 }
1449
1450 APFloat::opStatus
1451 APFloat::divideSpecials(const APFloat &rhs)
1452 {
1453   switch (convolve(category, rhs.category)) {
1454   default:
1455     llvm_unreachable(0);
1456
1457   case convolve(fcNaN, fcZero):
1458   case convolve(fcNaN, fcNormal):
1459   case convolve(fcNaN, fcInfinity):
1460   case convolve(fcNaN, fcNaN):
1461   case convolve(fcInfinity, fcZero):
1462   case convolve(fcInfinity, fcNormal):
1463   case convolve(fcZero, fcInfinity):
1464   case convolve(fcZero, fcNormal):
1465     return opOK;
1466
1467   case convolve(fcZero, fcNaN):
1468   case convolve(fcNormal, fcNaN):
1469   case convolve(fcInfinity, fcNaN):
1470     category = fcNaN;
1471     copySignificand(rhs);
1472     return opOK;
1473
1474   case convolve(fcNormal, fcInfinity):
1475     category = fcZero;
1476     return opOK;
1477
1478   case convolve(fcNormal, fcZero):
1479     category = fcInfinity;
1480     return opDivByZero;
1481
1482   case convolve(fcInfinity, fcInfinity):
1483   case convolve(fcZero, fcZero):
1484     makeNaN();
1485     return opInvalidOp;
1486
1487   case convolve(fcNormal, fcNormal):
1488     return opOK;
1489   }
1490 }
1491
1492 APFloat::opStatus
1493 APFloat::modSpecials(const APFloat &rhs)
1494 {
1495   switch (convolve(category, rhs.category)) {
1496   default:
1497     llvm_unreachable(0);
1498
1499   case convolve(fcNaN, fcZero):
1500   case convolve(fcNaN, fcNormal):
1501   case convolve(fcNaN, fcInfinity):
1502   case convolve(fcNaN, fcNaN):
1503   case convolve(fcZero, fcInfinity):
1504   case convolve(fcZero, fcNormal):
1505   case convolve(fcNormal, fcInfinity):
1506     return opOK;
1507
1508   case convolve(fcZero, fcNaN):
1509   case convolve(fcNormal, fcNaN):
1510   case convolve(fcInfinity, fcNaN):
1511     category = fcNaN;
1512     copySignificand(rhs);
1513     return opOK;
1514
1515   case convolve(fcNormal, fcZero):
1516   case convolve(fcInfinity, fcZero):
1517   case convolve(fcInfinity, fcNormal):
1518   case convolve(fcInfinity, fcInfinity):
1519   case convolve(fcZero, fcZero):
1520     makeNaN();
1521     return opInvalidOp;
1522
1523   case convolve(fcNormal, fcNormal):
1524     return opOK;
1525   }
1526 }
1527
1528 /* Change sign.  */
1529 void
1530 APFloat::changeSign()
1531 {
1532   /* Look mummy, this one's easy.  */
1533   sign = !sign;
1534 }
1535
1536 void
1537 APFloat::clearSign()
1538 {
1539   /* So is this one. */
1540   sign = 0;
1541 }
1542
1543 void
1544 APFloat::copySign(const APFloat &rhs)
1545 {
1546   /* And this one. */
1547   sign = rhs.sign;
1548 }
1549
1550 /* Normalized addition or subtraction.  */
1551 APFloat::opStatus
1552 APFloat::addOrSubtract(const APFloat &rhs, roundingMode rounding_mode,
1553                        bool subtract)
1554 {
1555   opStatus fs;
1556
1557   assertArithmeticOK(*semantics);
1558
1559   fs = addOrSubtractSpecials(rhs, subtract);
1560
1561   /* This return code means it was not a simple case.  */
1562   if (fs == opDivByZero) {
1563     lostFraction lost_fraction;
1564
1565     lost_fraction = addOrSubtractSignificand(rhs, subtract);
1566     fs = normalize(rounding_mode, lost_fraction);
1567
1568     /* Can only be zero if we lost no fraction.  */
1569     assert(category != fcZero || lost_fraction == lfExactlyZero);
1570   }
1571
1572   /* If two numbers add (exactly) to zero, IEEE 754 decrees it is a
1573      positive zero unless rounding to minus infinity, except that
1574      adding two like-signed zeroes gives that zero.  */
1575   if (category == fcZero) {
1576     if (rhs.category != fcZero || (sign == rhs.sign) == subtract)
1577       sign = (rounding_mode == rmTowardNegative);
1578   }
1579
1580   return fs;
1581 }
1582
1583 /* Normalized addition.  */
1584 APFloat::opStatus
1585 APFloat::add(const APFloat &rhs, roundingMode rounding_mode)
1586 {
1587   return addOrSubtract(rhs, rounding_mode, false);
1588 }
1589
1590 /* Normalized subtraction.  */
1591 APFloat::opStatus
1592 APFloat::subtract(const APFloat &rhs, roundingMode rounding_mode)
1593 {
1594   return addOrSubtract(rhs, rounding_mode, true);
1595 }
1596
1597 /* Normalized multiply.  */
1598 APFloat::opStatus
1599 APFloat::multiply(const APFloat &rhs, roundingMode rounding_mode)
1600 {
1601   opStatus fs;
1602
1603   assertArithmeticOK(*semantics);
1604   sign ^= rhs.sign;
1605   fs = multiplySpecials(rhs);
1606
1607   if (category == fcNormal) {
1608     lostFraction lost_fraction = multiplySignificand(rhs, 0);
1609     fs = normalize(rounding_mode, lost_fraction);
1610     if (lost_fraction != lfExactlyZero)
1611       fs = (opStatus) (fs | opInexact);
1612   }
1613
1614   return fs;
1615 }
1616
1617 /* Normalized divide.  */
1618 APFloat::opStatus
1619 APFloat::divide(const APFloat &rhs, roundingMode rounding_mode)
1620 {
1621   opStatus fs;
1622
1623   assertArithmeticOK(*semantics);
1624   sign ^= rhs.sign;
1625   fs = divideSpecials(rhs);
1626
1627   if (category == fcNormal) {
1628     lostFraction lost_fraction = divideSignificand(rhs);
1629     fs = normalize(rounding_mode, lost_fraction);
1630     if (lost_fraction != lfExactlyZero)
1631       fs = (opStatus) (fs | opInexact);
1632   }
1633
1634   return fs;
1635 }
1636
1637 /* Normalized remainder.  This is not currently correct in all cases.  */
1638 APFloat::opStatus
1639 APFloat::remainder(const APFloat &rhs)
1640 {
1641   opStatus fs;
1642   APFloat V = *this;
1643   unsigned int origSign = sign;
1644
1645   assertArithmeticOK(*semantics);
1646   fs = V.divide(rhs, rmNearestTiesToEven);
1647   if (fs == opDivByZero)
1648     return fs;
1649
1650   int parts = partCount();
1651   integerPart *x = new integerPart[parts];
1652   bool ignored;
1653   fs = V.convertToInteger(x, parts * integerPartWidth, true,
1654                           rmNearestTiesToEven, &ignored);
1655   if (fs==opInvalidOp)
1656     return fs;
1657
1658   fs = V.convertFromZeroExtendedInteger(x, parts * integerPartWidth, true,
1659                                         rmNearestTiesToEven);
1660   assert(fs==opOK);   // should always work
1661
1662   fs = V.multiply(rhs, rmNearestTiesToEven);
1663   assert(fs==opOK || fs==opInexact);   // should not overflow or underflow
1664
1665   fs = subtract(V, rmNearestTiesToEven);
1666   assert(fs==opOK || fs==opInexact);   // likewise
1667
1668   if (isZero())
1669     sign = origSign;    // IEEE754 requires this
1670   delete[] x;
1671   return fs;
1672 }
1673
1674 /* Normalized llvm frem (C fmod).
1675    This is not currently correct in all cases.  */
1676 APFloat::opStatus
1677 APFloat::mod(const APFloat &rhs, roundingMode rounding_mode)
1678 {
1679   opStatus fs;
1680   assertArithmeticOK(*semantics);
1681   fs = modSpecials(rhs);
1682
1683   if (category == fcNormal && rhs.category == fcNormal) {
1684     APFloat V = *this;
1685     unsigned int origSign = sign;
1686
1687     fs = V.divide(rhs, rmNearestTiesToEven);
1688     if (fs == opDivByZero)
1689       return fs;
1690
1691     int parts = partCount();
1692     integerPart *x = new integerPart[parts];
1693     bool ignored;
1694     fs = V.convertToInteger(x, parts * integerPartWidth, true,
1695                             rmTowardZero, &ignored);
1696     if (fs==opInvalidOp)
1697       return fs;
1698
1699     fs = V.convertFromZeroExtendedInteger(x, parts * integerPartWidth, true,
1700                                           rmNearestTiesToEven);
1701     assert(fs==opOK);   // should always work
1702
1703     fs = V.multiply(rhs, rounding_mode);
1704     assert(fs==opOK || fs==opInexact);   // should not overflow or underflow
1705
1706     fs = subtract(V, rounding_mode);
1707     assert(fs==opOK || fs==opInexact);   // likewise
1708
1709     if (isZero())
1710       sign = origSign;    // IEEE754 requires this
1711     delete[] x;
1712   }
1713   return fs;
1714 }
1715
1716 /* Normalized fused-multiply-add.  */
1717 APFloat::opStatus
1718 APFloat::fusedMultiplyAdd(const APFloat &multiplicand,
1719                           const APFloat &addend,
1720                           roundingMode rounding_mode)
1721 {
1722   opStatus fs;
1723
1724   assertArithmeticOK(*semantics);
1725
1726   /* Post-multiplication sign, before addition.  */
1727   sign ^= multiplicand.sign;
1728
1729   /* If and only if all arguments are normal do we need to do an
1730      extended-precision calculation.  */
1731   if (category == fcNormal &&
1732       multiplicand.category == fcNormal &&
1733       addend.category == fcNormal) {
1734     lostFraction lost_fraction;
1735
1736     lost_fraction = multiplySignificand(multiplicand, &addend);
1737     fs = normalize(rounding_mode, lost_fraction);
1738     if (lost_fraction != lfExactlyZero)
1739       fs = (opStatus) (fs | opInexact);
1740
1741     /* If two numbers add (exactly) to zero, IEEE 754 decrees it is a
1742        positive zero unless rounding to minus infinity, except that
1743        adding two like-signed zeroes gives that zero.  */
1744     if (category == fcZero && sign != addend.sign)
1745       sign = (rounding_mode == rmTowardNegative);
1746   } else {
1747     fs = multiplySpecials(multiplicand);
1748
1749     /* FS can only be opOK or opInvalidOp.  There is no more work
1750        to do in the latter case.  The IEEE-754R standard says it is
1751        implementation-defined in this case whether, if ADDEND is a
1752        quiet NaN, we raise invalid op; this implementation does so.
1753
1754        If we need to do the addition we can do so with normal
1755        precision.  */
1756     if (fs == opOK)
1757       fs = addOrSubtract(addend, rounding_mode, false);
1758   }
1759
1760   return fs;
1761 }
1762
1763 /* Comparison requires normalized numbers.  */
1764 APFloat::cmpResult
1765 APFloat::compare(const APFloat &rhs) const
1766 {
1767   cmpResult result;
1768
1769   assertArithmeticOK(*semantics);
1770   assert(semantics == rhs.semantics);
1771
1772   switch (convolve(category, rhs.category)) {
1773   default:
1774     llvm_unreachable(0);
1775
1776   case convolve(fcNaN, fcZero):
1777   case convolve(fcNaN, fcNormal):
1778   case convolve(fcNaN, fcInfinity):
1779   case convolve(fcNaN, fcNaN):
1780   case convolve(fcZero, fcNaN):
1781   case convolve(fcNormal, fcNaN):
1782   case convolve(fcInfinity, fcNaN):
1783     return cmpUnordered;
1784
1785   case convolve(fcInfinity, fcNormal):
1786   case convolve(fcInfinity, fcZero):
1787   case convolve(fcNormal, fcZero):
1788     if (sign)
1789       return cmpLessThan;
1790     else
1791       return cmpGreaterThan;
1792
1793   case convolve(fcNormal, fcInfinity):
1794   case convolve(fcZero, fcInfinity):
1795   case convolve(fcZero, fcNormal):
1796     if (rhs.sign)
1797       return cmpGreaterThan;
1798     else
1799       return cmpLessThan;
1800
1801   case convolve(fcInfinity, fcInfinity):
1802     if (sign == rhs.sign)
1803       return cmpEqual;
1804     else if (sign)
1805       return cmpLessThan;
1806     else
1807       return cmpGreaterThan;
1808
1809   case convolve(fcZero, fcZero):
1810     return cmpEqual;
1811
1812   case convolve(fcNormal, fcNormal):
1813     break;
1814   }
1815
1816   /* Two normal numbers.  Do they have the same sign?  */
1817   if (sign != rhs.sign) {
1818     if (sign)
1819       result = cmpLessThan;
1820     else
1821       result = cmpGreaterThan;
1822   } else {
1823     /* Compare absolute values; invert result if negative.  */
1824     result = compareAbsoluteValue(rhs);
1825
1826     if (sign) {
1827       if (result == cmpLessThan)
1828         result = cmpGreaterThan;
1829       else if (result == cmpGreaterThan)
1830         result = cmpLessThan;
1831     }
1832   }
1833
1834   return result;
1835 }
1836
1837 /// APFloat::convert - convert a value of one floating point type to another.
1838 /// The return value corresponds to the IEEE754 exceptions.  *losesInfo
1839 /// records whether the transformation lost information, i.e. whether
1840 /// converting the result back to the original type will produce the
1841 /// original value (this is almost the same as return value==fsOK, but there
1842 /// are edge cases where this is not so).
1843
1844 APFloat::opStatus
1845 APFloat::convert(const fltSemantics &toSemantics,
1846                  roundingMode rounding_mode, bool *losesInfo)
1847 {
1848   lostFraction lostFraction;
1849   unsigned int newPartCount, oldPartCount;
1850   opStatus fs;
1851
1852   assertArithmeticOK(*semantics);
1853   assertArithmeticOK(toSemantics);
1854   lostFraction = lfExactlyZero;
1855   newPartCount = partCountForBits(toSemantics.precision + 1);
1856   oldPartCount = partCount();
1857
1858   /* Handle storage complications.  If our new form is wider,
1859      re-allocate our bit pattern into wider storage.  If it is
1860      narrower, we ignore the excess parts, but if narrowing to a
1861      single part we need to free the old storage.
1862      Be careful not to reference significandParts for zeroes
1863      and infinities, since it aborts.  */
1864   if (newPartCount > oldPartCount) {
1865     integerPart *newParts;
1866     newParts = new integerPart[newPartCount];
1867     APInt::tcSet(newParts, 0, newPartCount);
1868     if (category==fcNormal || category==fcNaN)
1869       APInt::tcAssign(newParts, significandParts(), oldPartCount);
1870     freeSignificand();
1871     significand.parts = newParts;
1872   } else if (newPartCount < oldPartCount) {
1873     /* Capture any lost fraction through truncation of parts so we get
1874        correct rounding whilst normalizing.  */
1875     if (category==fcNormal)
1876       lostFraction = lostFractionThroughTruncation
1877         (significandParts(), oldPartCount, toSemantics.precision);
1878     if (newPartCount == 1) {
1879         integerPart newPart = 0;
1880         if (category==fcNormal || category==fcNaN)
1881           newPart = significandParts()[0];
1882         freeSignificand();
1883         significand.part = newPart;
1884     }
1885   }
1886
1887   if (category == fcNormal) {
1888     /* Re-interpret our bit-pattern.  */
1889     exponent += toSemantics.precision - semantics->precision;
1890     semantics = &toSemantics;
1891     fs = normalize(rounding_mode, lostFraction);
1892     *losesInfo = (fs != opOK);
1893   } else if (category == fcNaN) {
1894     int shift = toSemantics.precision - semantics->precision;
1895     // Do this now so significandParts gets the right answer
1896     const fltSemantics *oldSemantics = semantics;
1897     semantics = &toSemantics;
1898     *losesInfo = false;
1899     // No normalization here, just truncate
1900     if (shift>0)
1901       APInt::tcShiftLeft(significandParts(), newPartCount, shift);
1902     else if (shift < 0) {
1903       unsigned ushift = -shift;
1904       // Figure out if we are losing information.  This happens
1905       // if are shifting out something other than 0s, or if the x87 long
1906       // double input did not have its integer bit set (pseudo-NaN), or if the
1907       // x87 long double input did not have its QNan bit set (because the x87
1908       // hardware sets this bit when converting a lower-precision NaN to
1909       // x87 long double).
1910       if (APInt::tcLSB(significandParts(), newPartCount) < ushift)
1911         *losesInfo = true;
1912       if (oldSemantics == &APFloat::x87DoubleExtended &&
1913           (!(*significandParts() & 0x8000000000000000ULL) ||
1914            !(*significandParts() & 0x4000000000000000ULL)))
1915         *losesInfo = true;
1916       APInt::tcShiftRight(significandParts(), newPartCount, ushift);
1917     }
1918     // gcc forces the Quiet bit on, which means (float)(double)(float_sNan)
1919     // does not give you back the same bits.  This is dubious, and we
1920     // don't currently do it.  You're really supposed to get
1921     // an invalid operation signal at runtime, but nobody does that.
1922     fs = opOK;
1923   } else {
1924     semantics = &toSemantics;
1925     fs = opOK;
1926     *losesInfo = false;
1927   }
1928
1929   return fs;
1930 }
1931
1932 /* Convert a floating point number to an integer according to the
1933    rounding mode.  If the rounded integer value is out of range this
1934    returns an invalid operation exception and the contents of the
1935    destination parts are unspecified.  If the rounded value is in
1936    range but the floating point number is not the exact integer, the C
1937    standard doesn't require an inexact exception to be raised.  IEEE
1938    854 does require it so we do that.
1939
1940    Note that for conversions to integer type the C standard requires
1941    round-to-zero to always be used.  */
1942 APFloat::opStatus
1943 APFloat::convertToSignExtendedInteger(integerPart *parts, unsigned int width,
1944                                       bool isSigned,
1945                                       roundingMode rounding_mode,
1946                                       bool *isExact) const
1947 {
1948   lostFraction lost_fraction;
1949   const integerPart *src;
1950   unsigned int dstPartsCount, truncatedBits;
1951
1952   assertArithmeticOK(*semantics);
1953
1954   *isExact = false;
1955
1956   /* Handle the three special cases first.  */
1957   if (category == fcInfinity || category == fcNaN)
1958     return opInvalidOp;
1959
1960   dstPartsCount = partCountForBits(width);
1961
1962   if (category == fcZero) {
1963     APInt::tcSet(parts, 0, dstPartsCount);
1964     // Negative zero can't be represented as an int.
1965     *isExact = !sign;
1966     return opOK;
1967   }
1968
1969   src = significandParts();
1970
1971   /* Step 1: place our absolute value, with any fraction truncated, in
1972      the destination.  */
1973   if (exponent < 0) {
1974     /* Our absolute value is less than one; truncate everything.  */
1975     APInt::tcSet(parts, 0, dstPartsCount);
1976     /* For exponent -1 the integer bit represents .5, look at that.
1977        For smaller exponents leftmost truncated bit is 0. */
1978     truncatedBits = semantics->precision -1U - exponent;
1979   } else {
1980     /* We want the most significant (exponent + 1) bits; the rest are
1981        truncated.  */
1982     unsigned int bits = exponent + 1U;
1983
1984     /* Hopelessly large in magnitude?  */
1985     if (bits > width)
1986       return opInvalidOp;
1987
1988     if (bits < semantics->precision) {
1989       /* We truncate (semantics->precision - bits) bits.  */
1990       truncatedBits = semantics->precision - bits;
1991       APInt::tcExtract(parts, dstPartsCount, src, bits, truncatedBits);
1992     } else {
1993       /* We want at least as many bits as are available.  */
1994       APInt::tcExtract(parts, dstPartsCount, src, semantics->precision, 0);
1995       APInt::tcShiftLeft(parts, dstPartsCount, bits - semantics->precision);
1996       truncatedBits = 0;
1997     }
1998   }
1999
2000   /* Step 2: work out any lost fraction, and increment the absolute
2001      value if we would round away from zero.  */
2002   if (truncatedBits) {
2003     lost_fraction = lostFractionThroughTruncation(src, partCount(),
2004                                                   truncatedBits);
2005     if (lost_fraction != lfExactlyZero &&
2006         roundAwayFromZero(rounding_mode, lost_fraction, truncatedBits)) {
2007       if (APInt::tcIncrement(parts, dstPartsCount))
2008         return opInvalidOp;     /* Overflow.  */
2009     }
2010   } else {
2011     lost_fraction = lfExactlyZero;
2012   }
2013
2014   /* Step 3: check if we fit in the destination.  */
2015   unsigned int omsb = APInt::tcMSB(parts, dstPartsCount) + 1;
2016
2017   if (sign) {
2018     if (!isSigned) {
2019       /* Negative numbers cannot be represented as unsigned.  */
2020       if (omsb != 0)
2021         return opInvalidOp;
2022     } else {
2023       /* It takes omsb bits to represent the unsigned integer value.
2024          We lose a bit for the sign, but care is needed as the
2025          maximally negative integer is a special case.  */
2026       if (omsb == width && APInt::tcLSB(parts, dstPartsCount) + 1 != omsb)
2027         return opInvalidOp;
2028
2029       /* This case can happen because of rounding.  */
2030       if (omsb > width)
2031         return opInvalidOp;
2032     }
2033
2034     APInt::tcNegate (parts, dstPartsCount);
2035   } else {
2036     if (omsb >= width + !isSigned)
2037       return opInvalidOp;
2038   }
2039
2040   if (lost_fraction == lfExactlyZero) {
2041     *isExact = true;
2042     return opOK;
2043   } else
2044     return opInexact;
2045 }
2046
2047 /* Same as convertToSignExtendedInteger, except we provide
2048    deterministic values in case of an invalid operation exception,
2049    namely zero for NaNs and the minimal or maximal value respectively
2050    for underflow or overflow.
2051    The *isExact output tells whether the result is exact, in the sense
2052    that converting it back to the original floating point type produces
2053    the original value.  This is almost equivalent to result==opOK,
2054    except for negative zeroes.
2055 */
2056 APFloat::opStatus
2057 APFloat::convertToInteger(integerPart *parts, unsigned int width,
2058                           bool isSigned,
2059                           roundingMode rounding_mode, bool *isExact) const
2060 {
2061   opStatus fs;
2062
2063   fs = convertToSignExtendedInteger(parts, width, isSigned, rounding_mode,
2064                                     isExact);
2065
2066   if (fs == opInvalidOp) {
2067     unsigned int bits, dstPartsCount;
2068
2069     dstPartsCount = partCountForBits(width);
2070
2071     if (category == fcNaN)
2072       bits = 0;
2073     else if (sign)
2074       bits = isSigned;
2075     else
2076       bits = width - isSigned;
2077
2078     APInt::tcSetLeastSignificantBits(parts, dstPartsCount, bits);
2079     if (sign && isSigned)
2080       APInt::tcShiftLeft(parts, dstPartsCount, width - 1);
2081   }
2082
2083   return fs;
2084 }
2085
2086 /* Convert an unsigned integer SRC to a floating point number,
2087    rounding according to ROUNDING_MODE.  The sign of the floating
2088    point number is not modified.  */
2089 APFloat::opStatus
2090 APFloat::convertFromUnsignedParts(const integerPart *src,
2091                                   unsigned int srcCount,
2092                                   roundingMode rounding_mode)
2093 {
2094   unsigned int omsb, precision, dstCount;
2095   integerPart *dst;
2096   lostFraction lost_fraction;
2097
2098   assertArithmeticOK(*semantics);
2099   category = fcNormal;
2100   omsb = APInt::tcMSB(src, srcCount) + 1;
2101   dst = significandParts();
2102   dstCount = partCount();
2103   precision = semantics->precision;
2104
2105   /* We want the most significant PRECISON bits of SRC.  There may not
2106      be that many; extract what we can.  */
2107   if (precision <= omsb) {
2108     exponent = omsb - 1;
2109     lost_fraction = lostFractionThroughTruncation(src, srcCount,
2110                                                   omsb - precision);
2111     APInt::tcExtract(dst, dstCount, src, precision, omsb - precision);
2112   } else {
2113     exponent = precision - 1;
2114     lost_fraction = lfExactlyZero;
2115     APInt::tcExtract(dst, dstCount, src, omsb, 0);
2116   }
2117
2118   return normalize(rounding_mode, lost_fraction);
2119 }
2120
2121 APFloat::opStatus
2122 APFloat::convertFromAPInt(const APInt &Val,
2123                           bool isSigned,
2124                           roundingMode rounding_mode)
2125 {
2126   unsigned int partCount = Val.getNumWords();
2127   APInt api = Val;
2128
2129   sign = false;
2130   if (isSigned && api.isNegative()) {
2131     sign = true;
2132     api = -api;
2133   }
2134
2135   return convertFromUnsignedParts(api.getRawData(), partCount, rounding_mode);
2136 }
2137
2138 /* Convert a two's complement integer SRC to a floating point number,
2139    rounding according to ROUNDING_MODE.  ISSIGNED is true if the
2140    integer is signed, in which case it must be sign-extended.  */
2141 APFloat::opStatus
2142 APFloat::convertFromSignExtendedInteger(const integerPart *src,
2143                                         unsigned int srcCount,
2144                                         bool isSigned,
2145                                         roundingMode rounding_mode)
2146 {
2147   opStatus status;
2148
2149   assertArithmeticOK(*semantics);
2150   if (isSigned &&
2151       APInt::tcExtractBit(src, srcCount * integerPartWidth - 1)) {
2152     integerPart *copy;
2153
2154     /* If we're signed and negative negate a copy.  */
2155     sign = true;
2156     copy = new integerPart[srcCount];
2157     APInt::tcAssign(copy, src, srcCount);
2158     APInt::tcNegate(copy, srcCount);
2159     status = convertFromUnsignedParts(copy, srcCount, rounding_mode);
2160     delete [] copy;
2161   } else {
2162     sign = false;
2163     status = convertFromUnsignedParts(src, srcCount, rounding_mode);
2164   }
2165
2166   return status;
2167 }
2168
2169 /* FIXME: should this just take a const APInt reference?  */
2170 APFloat::opStatus
2171 APFloat::convertFromZeroExtendedInteger(const integerPart *parts,
2172                                         unsigned int width, bool isSigned,
2173                                         roundingMode rounding_mode)
2174 {
2175   unsigned int partCount = partCountForBits(width);
2176   APInt api = APInt(width, partCount, parts);
2177
2178   sign = false;
2179   if (isSigned && APInt::tcExtractBit(parts, width - 1)) {
2180     sign = true;
2181     api = -api;
2182   }
2183
2184   return convertFromUnsignedParts(api.getRawData(), partCount, rounding_mode);
2185 }
2186
2187 APFloat::opStatus
2188 APFloat::convertFromHexadecimalString(StringRef s, roundingMode rounding_mode)
2189 {
2190   lostFraction lost_fraction = lfExactlyZero;
2191   integerPart *significand;
2192   unsigned int bitPos, partsCount;
2193   StringRef::iterator dot, firstSignificantDigit;
2194
2195   zeroSignificand();
2196   exponent = 0;
2197   category = fcNormal;
2198
2199   significand = significandParts();
2200   partsCount = partCount();
2201   bitPos = partsCount * integerPartWidth;
2202
2203   /* Skip leading zeroes and any (hexa)decimal point.  */
2204   StringRef::iterator begin = s.begin();
2205   StringRef::iterator end = s.end();
2206   StringRef::iterator p = skipLeadingZeroesAndAnyDot(begin, end, &dot);
2207   firstSignificantDigit = p;
2208
2209   for (; p != end;) {
2210     integerPart hex_value;
2211
2212     if (*p == '.') {
2213       assert(dot == end && "String contains multiple dots");
2214       dot = p++;
2215       if (p == end) {
2216         break;
2217       }
2218     }
2219
2220     hex_value = hexDigitValue(*p);
2221     if (hex_value == -1U) {
2222       break;
2223     }
2224
2225     p++;
2226
2227     if (p == end) {
2228       break;
2229     } else {
2230       /* Store the number whilst 4-bit nibbles remain.  */
2231       if (bitPos) {
2232         bitPos -= 4;
2233         hex_value <<= bitPos % integerPartWidth;
2234         significand[bitPos / integerPartWidth] |= hex_value;
2235       } else {
2236         lost_fraction = trailingHexadecimalFraction(p, end, hex_value);
2237         while (p != end && hexDigitValue(*p) != -1U)
2238           p++;
2239         break;
2240       }
2241     }
2242   }
2243
2244   /* Hex floats require an exponent but not a hexadecimal point.  */
2245   assert(p != end && "Hex strings require an exponent");
2246   assert((*p == 'p' || *p == 'P') && "Invalid character in significand");
2247   assert(p != begin && "Significand has no digits");
2248   assert((dot == end || p - begin != 1) && "Significand has no digits");
2249
2250   /* Ignore the exponent if we are zero.  */
2251   if (p != firstSignificantDigit) {
2252     int expAdjustment;
2253
2254     /* Implicit hexadecimal point?  */
2255     if (dot == end)
2256       dot = p;
2257
2258     /* Calculate the exponent adjustment implicit in the number of
2259        significant digits.  */
2260     expAdjustment = static_cast<int>(dot - firstSignificantDigit);
2261     if (expAdjustment < 0)
2262       expAdjustment++;
2263     expAdjustment = expAdjustment * 4 - 1;
2264
2265     /* Adjust for writing the significand starting at the most
2266        significant nibble.  */
2267     expAdjustment += semantics->precision;
2268     expAdjustment -= partsCount * integerPartWidth;
2269
2270     /* Adjust for the given exponent.  */
2271     exponent = totalExponent(p + 1, end, expAdjustment);
2272   }
2273
2274   return normalize(rounding_mode, lost_fraction);
2275 }
2276
2277 APFloat::opStatus
2278 APFloat::roundSignificandWithExponent(const integerPart *decSigParts,
2279                                       unsigned sigPartCount, int exp,
2280                                       roundingMode rounding_mode)
2281 {
2282   unsigned int parts, pow5PartCount;
2283   fltSemantics calcSemantics = { 32767, -32767, 0, true };
2284   integerPart pow5Parts[maxPowerOfFiveParts];
2285   bool isNearest;
2286
2287   isNearest = (rounding_mode == rmNearestTiesToEven ||
2288                rounding_mode == rmNearestTiesToAway);
2289
2290   parts = partCountForBits(semantics->precision + 11);
2291
2292   /* Calculate pow(5, abs(exp)).  */
2293   pow5PartCount = powerOf5(pow5Parts, exp >= 0 ? exp: -exp);
2294
2295   for (;; parts *= 2) {
2296     opStatus sigStatus, powStatus;
2297     unsigned int excessPrecision, truncatedBits;
2298
2299     calcSemantics.precision = parts * integerPartWidth - 1;
2300     excessPrecision = calcSemantics.precision - semantics->precision;
2301     truncatedBits = excessPrecision;
2302
2303     APFloat decSig(calcSemantics, fcZero, sign);
2304     APFloat pow5(calcSemantics, fcZero, false);
2305
2306     sigStatus = decSig.convertFromUnsignedParts(decSigParts, sigPartCount,
2307                                                 rmNearestTiesToEven);
2308     powStatus = pow5.convertFromUnsignedParts(pow5Parts, pow5PartCount,
2309                                               rmNearestTiesToEven);
2310     /* Add exp, as 10^n = 5^n * 2^n.  */
2311     decSig.exponent += exp;
2312
2313     lostFraction calcLostFraction;
2314     integerPart HUerr, HUdistance;
2315     unsigned int powHUerr;
2316
2317     if (exp >= 0) {
2318       /* multiplySignificand leaves the precision-th bit set to 1.  */
2319       calcLostFraction = decSig.multiplySignificand(pow5, NULL);
2320       powHUerr = powStatus != opOK;
2321     } else {
2322       calcLostFraction = decSig.divideSignificand(pow5);
2323       /* Denormal numbers have less precision.  */
2324       if (decSig.exponent < semantics->minExponent) {
2325         excessPrecision += (semantics->minExponent - decSig.exponent);
2326         truncatedBits = excessPrecision;
2327         if (excessPrecision > calcSemantics.precision)
2328           excessPrecision = calcSemantics.precision;
2329       }
2330       /* Extra half-ulp lost in reciprocal of exponent.  */
2331       powHUerr = (powStatus == opOK && calcLostFraction == lfExactlyZero) ? 0:2;
2332     }
2333
2334     /* Both multiplySignificand and divideSignificand return the
2335        result with the integer bit set.  */
2336     assert(APInt::tcExtractBit
2337            (decSig.significandParts(), calcSemantics.precision - 1) == 1);
2338
2339     HUerr = HUerrBound(calcLostFraction != lfExactlyZero, sigStatus != opOK,
2340                        powHUerr);
2341     HUdistance = 2 * ulpsFromBoundary(decSig.significandParts(),
2342                                       excessPrecision, isNearest);
2343
2344     /* Are we guaranteed to round correctly if we truncate?  */
2345     if (HUdistance >= HUerr) {
2346       APInt::tcExtract(significandParts(), partCount(), decSig.significandParts(),
2347                        calcSemantics.precision - excessPrecision,
2348                        excessPrecision);
2349       /* Take the exponent of decSig.  If we tcExtract-ed less bits
2350          above we must adjust our exponent to compensate for the
2351          implicit right shift.  */
2352       exponent = (decSig.exponent + semantics->precision
2353                   - (calcSemantics.precision - excessPrecision));
2354       calcLostFraction = lostFractionThroughTruncation(decSig.significandParts(),
2355                                                        decSig.partCount(),
2356                                                        truncatedBits);
2357       return normalize(rounding_mode, calcLostFraction);
2358     }
2359   }
2360 }
2361
2362 APFloat::opStatus
2363 APFloat::convertFromDecimalString(StringRef str, roundingMode rounding_mode)
2364 {
2365   decimalInfo D;
2366   opStatus fs;
2367
2368   /* Scan the text.  */
2369   StringRef::iterator p = str.begin();
2370   interpretDecimal(p, str.end(), &D);
2371
2372   /* Handle the quick cases.  First the case of no significant digits,
2373      i.e. zero, and then exponents that are obviously too large or too
2374      small.  Writing L for log 10 / log 2, a number d.ddddd*10^exp
2375      definitely overflows if
2376
2377            (exp - 1) * L >= maxExponent
2378
2379      and definitely underflows to zero where
2380
2381            (exp + 1) * L <= minExponent - precision
2382
2383      With integer arithmetic the tightest bounds for L are
2384
2385            93/28 < L < 196/59            [ numerator <= 256 ]
2386            42039/12655 < L < 28738/8651  [ numerator <= 65536 ]
2387   */
2388
2389   if (decDigitValue(*D.firstSigDigit) >= 10U) {
2390     category = fcZero;
2391     fs = opOK;
2392
2393   /* Check whether the normalized exponent is high enough to overflow
2394      max during the log-rebasing in the max-exponent check below. */
2395   } else if (D.normalizedExponent - 1 > INT_MAX / 42039) {
2396     fs = handleOverflow(rounding_mode);
2397
2398   /* If it wasn't, then it also wasn't high enough to overflow max
2399      during the log-rebasing in the min-exponent check.  Check that it
2400      won't overflow min in either check, then perform the min-exponent
2401      check. */
2402   } else if (D.normalizedExponent - 1 < INT_MIN / 42039 ||
2403              (D.normalizedExponent + 1) * 28738 <=
2404                8651 * (semantics->minExponent - (int) semantics->precision)) {
2405     /* Underflow to zero and round.  */
2406     zeroSignificand();
2407     fs = normalize(rounding_mode, lfLessThanHalf);
2408
2409   /* We can finally safely perform the max-exponent check. */
2410   } else if ((D.normalizedExponent - 1) * 42039
2411              >= 12655 * semantics->maxExponent) {
2412     /* Overflow and round.  */
2413     fs = handleOverflow(rounding_mode);
2414   } else {
2415     integerPart *decSignificand;
2416     unsigned int partCount;
2417
2418     /* A tight upper bound on number of bits required to hold an
2419        N-digit decimal integer is N * 196 / 59.  Allocate enough space
2420        to hold the full significand, and an extra part required by
2421        tcMultiplyPart.  */
2422     partCount = static_cast<unsigned int>(D.lastSigDigit - D.firstSigDigit) + 1;
2423     partCount = partCountForBits(1 + 196 * partCount / 59);
2424     decSignificand = new integerPart[partCount + 1];
2425     partCount = 0;
2426
2427     /* Convert to binary efficiently - we do almost all multiplication
2428        in an integerPart.  When this would overflow do we do a single
2429        bignum multiplication, and then revert again to multiplication
2430        in an integerPart.  */
2431     do {
2432       integerPart decValue, val, multiplier;
2433
2434       val = 0;
2435       multiplier = 1;
2436
2437       do {
2438         if (*p == '.') {
2439           p++;
2440           if (p == str.end()) {
2441             break;
2442           }
2443         }
2444         decValue = decDigitValue(*p++);
2445         assert(decValue < 10U && "Invalid character in significand");
2446         multiplier *= 10;
2447         val = val * 10 + decValue;
2448         /* The maximum number that can be multiplied by ten with any
2449            digit added without overflowing an integerPart.  */
2450       } while (p <= D.lastSigDigit && multiplier <= (~ (integerPart) 0 - 9) / 10);
2451
2452       /* Multiply out the current part.  */
2453       APInt::tcMultiplyPart(decSignificand, decSignificand, multiplier, val,
2454                             partCount, partCount + 1, false);
2455
2456       /* If we used another part (likely but not guaranteed), increase
2457          the count.  */
2458       if (decSignificand[partCount])
2459         partCount++;
2460     } while (p <= D.lastSigDigit);
2461
2462     category = fcNormal;
2463     fs = roundSignificandWithExponent(decSignificand, partCount,
2464                                       D.exponent, rounding_mode);
2465
2466     delete [] decSignificand;
2467   }
2468
2469   return fs;
2470 }
2471
2472 APFloat::opStatus
2473 APFloat::convertFromString(StringRef str, roundingMode rounding_mode)
2474 {
2475   assertArithmeticOK(*semantics);
2476   assert(!str.empty() && "Invalid string length");
2477
2478   /* Handle a leading minus sign.  */
2479   StringRef::iterator p = str.begin();
2480   size_t slen = str.size();
2481   sign = *p == '-' ? 1 : 0;
2482   if (*p == '-' || *p == '+') {
2483     p++;
2484     slen--;
2485     assert(slen && "String has no digits");
2486   }
2487
2488   if (slen >= 2 && p[0] == '0' && (p[1] == 'x' || p[1] == 'X')) {
2489     assert(slen - 2 && "Invalid string");
2490     return convertFromHexadecimalString(StringRef(p + 2, slen - 2),
2491                                         rounding_mode);
2492   }
2493
2494   return convertFromDecimalString(StringRef(p, slen), rounding_mode);
2495 }
2496
2497 /* Write out a hexadecimal representation of the floating point value
2498    to DST, which must be of sufficient size, in the C99 form
2499    [-]0xh.hhhhp[+-]d.  Return the number of characters written,
2500    excluding the terminating NUL.
2501
2502    If UPPERCASE, the output is in upper case, otherwise in lower case.
2503
2504    HEXDIGITS digits appear altogether, rounding the value if
2505    necessary.  If HEXDIGITS is 0, the minimal precision to display the
2506    number precisely is used instead.  If nothing would appear after
2507    the decimal point it is suppressed.
2508
2509    The decimal exponent is always printed and has at least one digit.
2510    Zero values display an exponent of zero.  Infinities and NaNs
2511    appear as "infinity" or "nan" respectively.
2512
2513    The above rules are as specified by C99.  There is ambiguity about
2514    what the leading hexadecimal digit should be.  This implementation
2515    uses whatever is necessary so that the exponent is displayed as
2516    stored.  This implies the exponent will fall within the IEEE format
2517    range, and the leading hexadecimal digit will be 0 (for denormals),
2518    1 (normal numbers) or 2 (normal numbers rounded-away-from-zero with
2519    any other digits zero).
2520 */
2521 unsigned int
2522 APFloat::convertToHexString(char *dst, unsigned int hexDigits,
2523                             bool upperCase, roundingMode rounding_mode) const
2524 {
2525   char *p;
2526
2527   assertArithmeticOK(*semantics);
2528
2529   p = dst;
2530   if (sign)
2531     *dst++ = '-';
2532
2533   switch (category) {
2534   case fcInfinity:
2535     memcpy (dst, upperCase ? infinityU: infinityL, sizeof infinityU - 1);
2536     dst += sizeof infinityL - 1;
2537     break;
2538
2539   case fcNaN:
2540     memcpy (dst, upperCase ? NaNU: NaNL, sizeof NaNU - 1);
2541     dst += sizeof NaNU - 1;
2542     break;
2543
2544   case fcZero:
2545     *dst++ = '0';
2546     *dst++ = upperCase ? 'X': 'x';
2547     *dst++ = '0';
2548     if (hexDigits > 1) {
2549       *dst++ = '.';
2550       memset (dst, '0', hexDigits - 1);
2551       dst += hexDigits - 1;
2552     }
2553     *dst++ = upperCase ? 'P': 'p';
2554     *dst++ = '0';
2555     break;
2556
2557   case fcNormal:
2558     dst = convertNormalToHexString (dst, hexDigits, upperCase, rounding_mode);
2559     break;
2560   }
2561
2562   *dst = 0;
2563
2564   return static_cast<unsigned int>(dst - p);
2565 }
2566
2567 /* Does the hard work of outputting the correctly rounded hexadecimal
2568    form of a normal floating point number with the specified number of
2569    hexadecimal digits.  If HEXDIGITS is zero the minimum number of
2570    digits necessary to print the value precisely is output.  */
2571 char *
2572 APFloat::convertNormalToHexString(char *dst, unsigned int hexDigits,
2573                                   bool upperCase,
2574                                   roundingMode rounding_mode) const
2575 {
2576   unsigned int count, valueBits, shift, partsCount, outputDigits;
2577   const char *hexDigitChars;
2578   const integerPart *significand;
2579   char *p;
2580   bool roundUp;
2581
2582   *dst++ = '0';
2583   *dst++ = upperCase ? 'X': 'x';
2584
2585   roundUp = false;
2586   hexDigitChars = upperCase ? hexDigitsUpper: hexDigitsLower;
2587
2588   significand = significandParts();
2589   partsCount = partCount();
2590
2591   /* +3 because the first digit only uses the single integer bit, so
2592      we have 3 virtual zero most-significant-bits.  */
2593   valueBits = semantics->precision + 3;
2594   shift = integerPartWidth - valueBits % integerPartWidth;
2595
2596   /* The natural number of digits required ignoring trailing
2597      insignificant zeroes.  */
2598   outputDigits = (valueBits - significandLSB () + 3) / 4;
2599
2600   /* hexDigits of zero means use the required number for the
2601      precision.  Otherwise, see if we are truncating.  If we are,
2602      find out if we need to round away from zero.  */
2603   if (hexDigits) {
2604     if (hexDigits < outputDigits) {
2605       /* We are dropping non-zero bits, so need to check how to round.
2606          "bits" is the number of dropped bits.  */
2607       unsigned int bits;
2608       lostFraction fraction;
2609
2610       bits = valueBits - hexDigits * 4;
2611       fraction = lostFractionThroughTruncation (significand, partsCount, bits);
2612       roundUp = roundAwayFromZero(rounding_mode, fraction, bits);
2613     }
2614     outputDigits = hexDigits;
2615   }
2616
2617   /* Write the digits consecutively, and start writing in the location
2618      of the hexadecimal point.  We move the most significant digit
2619      left and add the hexadecimal point later.  */
2620   p = ++dst;
2621
2622   count = (valueBits + integerPartWidth - 1) / integerPartWidth;
2623
2624   while (outputDigits && count) {
2625     integerPart part;
2626
2627     /* Put the most significant integerPartWidth bits in "part".  */
2628     if (--count == partsCount)
2629       part = 0;  /* An imaginary higher zero part.  */
2630     else
2631       part = significand[count] << shift;
2632
2633     if (count && shift)
2634       part |= significand[count - 1] >> (integerPartWidth - shift);
2635
2636     /* Convert as much of "part" to hexdigits as we can.  */
2637     unsigned int curDigits = integerPartWidth / 4;
2638
2639     if (curDigits > outputDigits)
2640       curDigits = outputDigits;
2641     dst += partAsHex (dst, part, curDigits, hexDigitChars);
2642     outputDigits -= curDigits;
2643   }
2644
2645   if (roundUp) {
2646     char *q = dst;
2647
2648     /* Note that hexDigitChars has a trailing '0'.  */
2649     do {
2650       q--;
2651       *q = hexDigitChars[hexDigitValue (*q) + 1];
2652     } while (*q == '0');
2653     assert(q >= p);
2654   } else {
2655     /* Add trailing zeroes.  */
2656     memset (dst, '0', outputDigits);
2657     dst += outputDigits;
2658   }
2659
2660   /* Move the most significant digit to before the point, and if there
2661      is something after the decimal point add it.  This must come
2662      after rounding above.  */
2663   p[-1] = p[0];
2664   if (dst -1 == p)
2665     dst--;
2666   else
2667     p[0] = '.';
2668
2669   /* Finally output the exponent.  */
2670   *dst++ = upperCase ? 'P': 'p';
2671
2672   return writeSignedDecimal (dst, exponent);
2673 }
2674
2675 // For good performance it is desirable for different APFloats
2676 // to produce different integers.
2677 uint32_t
2678 APFloat::getHashValue() const
2679 {
2680   if (category==fcZero) return sign<<8 | semantics->precision ;
2681   else if (category==fcInfinity) return sign<<9 | semantics->precision;
2682   else if (category==fcNaN) return 1<<10 | semantics->precision;
2683   else {
2684     uint32_t hash = sign<<11 | semantics->precision | exponent<<12;
2685     const integerPart* p = significandParts();
2686     for (int i=partCount(); i>0; i--, p++)
2687       hash ^= ((uint32_t)*p) ^ (uint32_t)((*p)>>32);
2688     return hash;
2689   }
2690 }
2691
2692 // Conversion from APFloat to/from host float/double.  It may eventually be
2693 // possible to eliminate these and have everybody deal with APFloats, but that
2694 // will take a while.  This approach will not easily extend to long double.
2695 // Current implementation requires integerPartWidth==64, which is correct at
2696 // the moment but could be made more general.
2697
2698 // Denormals have exponent minExponent in APFloat, but minExponent-1 in
2699 // the actual IEEE respresentations.  We compensate for that here.
2700
2701 APInt
2702 APFloat::convertF80LongDoubleAPFloatToAPInt() const
2703 {
2704   assert(semantics == (const llvm::fltSemantics*)&x87DoubleExtended);
2705   assert(partCount()==2);
2706
2707   uint64_t myexponent, mysignificand;
2708
2709   if (category==fcNormal) {
2710     myexponent = exponent+16383; //bias
2711     mysignificand = significandParts()[0];
2712     if (myexponent==1 && !(mysignificand & 0x8000000000000000ULL))
2713       myexponent = 0;   // denormal
2714   } else if (category==fcZero) {
2715     myexponent = 0;
2716     mysignificand = 0;
2717   } else if (category==fcInfinity) {
2718     myexponent = 0x7fff;
2719     mysignificand = 0x8000000000000000ULL;
2720   } else {
2721     assert(category == fcNaN && "Unknown category");
2722     myexponent = 0x7fff;
2723     mysignificand = significandParts()[0];
2724   }
2725
2726   uint64_t words[2];
2727   words[0] = mysignificand;
2728   words[1] =  ((uint64_t)(sign & 1) << 15) |
2729               (myexponent & 0x7fffLL);
2730   return APInt(80, 2, words);
2731 }
2732
2733 APInt
2734 APFloat::convertPPCDoubleDoubleAPFloatToAPInt() const
2735 {
2736   assert(semantics == (const llvm::fltSemantics*)&PPCDoubleDouble);
2737   assert(partCount()==2);
2738
2739   uint64_t myexponent, mysignificand, myexponent2, mysignificand2;
2740
2741   if (category==fcNormal) {
2742     myexponent = exponent + 1023; //bias
2743     myexponent2 = exponent2 + 1023;
2744     mysignificand = significandParts()[0];
2745     mysignificand2 = significandParts()[1];
2746     if (myexponent==1 && !(mysignificand & 0x10000000000000LL))
2747       myexponent = 0;   // denormal
2748     if (myexponent2==1 && !(mysignificand2 & 0x10000000000000LL))
2749       myexponent2 = 0;   // denormal
2750   } else if (category==fcZero) {
2751     myexponent = 0;
2752     mysignificand = 0;
2753     myexponent2 = 0;
2754     mysignificand2 = 0;
2755   } else if (category==fcInfinity) {
2756     myexponent = 0x7ff;
2757     myexponent2 = 0;
2758     mysignificand = 0;
2759     mysignificand2 = 0;
2760   } else {
2761     assert(category == fcNaN && "Unknown category");
2762     myexponent = 0x7ff;
2763     mysignificand = significandParts()[0];
2764     myexponent2 = exponent2;
2765     mysignificand2 = significandParts()[1];
2766   }
2767
2768   uint64_t words[2];
2769   words[0] =  ((uint64_t)(sign & 1) << 63) |
2770               ((myexponent & 0x7ff) <<  52) |
2771               (mysignificand & 0xfffffffffffffLL);
2772   words[1] =  ((uint64_t)(sign2 & 1) << 63) |
2773               ((myexponent2 & 0x7ff) <<  52) |
2774               (mysignificand2 & 0xfffffffffffffLL);
2775   return APInt(128, 2, words);
2776 }
2777
2778 APInt
2779 APFloat::convertQuadrupleAPFloatToAPInt() const
2780 {
2781   assert(semantics == (const llvm::fltSemantics*)&IEEEquad);
2782   assert(partCount()==2);
2783
2784   uint64_t myexponent, mysignificand, mysignificand2;
2785
2786   if (category==fcNormal) {
2787     myexponent = exponent+16383; //bias
2788     mysignificand = significandParts()[0];
2789     mysignificand2 = significandParts()[1];
2790     if (myexponent==1 && !(mysignificand2 & 0x1000000000000LL))
2791       myexponent = 0;   // denormal
2792   } else if (category==fcZero) {
2793     myexponent = 0;
2794     mysignificand = mysignificand2 = 0;
2795   } else if (category==fcInfinity) {
2796     myexponent = 0x7fff;
2797     mysignificand = mysignificand2 = 0;
2798   } else {
2799     assert(category == fcNaN && "Unknown category!");
2800     myexponent = 0x7fff;
2801     mysignificand = significandParts()[0];
2802     mysignificand2 = significandParts()[1];
2803   }
2804
2805   uint64_t words[2];
2806   words[0] = mysignificand;
2807   words[1] = ((uint64_t)(sign & 1) << 63) |
2808              ((myexponent & 0x7fff) << 48) |
2809              (mysignificand2 & 0xffffffffffffLL);
2810
2811   return APInt(128, 2, words);
2812 }
2813
2814 APInt
2815 APFloat::convertDoubleAPFloatToAPInt() const
2816 {
2817   assert(semantics == (const llvm::fltSemantics*)&IEEEdouble);
2818   assert(partCount()==1);
2819
2820   uint64_t myexponent, mysignificand;
2821
2822   if (category==fcNormal) {
2823     myexponent = exponent+1023; //bias
2824     mysignificand = *significandParts();
2825     if (myexponent==1 && !(mysignificand & 0x10000000000000LL))
2826       myexponent = 0;   // denormal
2827   } else if (category==fcZero) {
2828     myexponent = 0;
2829     mysignificand = 0;
2830   } else if (category==fcInfinity) {
2831     myexponent = 0x7ff;
2832     mysignificand = 0;
2833   } else {
2834     assert(category == fcNaN && "Unknown category!");
2835     myexponent = 0x7ff;
2836     mysignificand = *significandParts();
2837   }
2838
2839   return APInt(64, ((((uint64_t)(sign & 1) << 63) |
2840                      ((myexponent & 0x7ff) <<  52) |
2841                      (mysignificand & 0xfffffffffffffLL))));
2842 }
2843
2844 APInt
2845 APFloat::convertFloatAPFloatToAPInt() const
2846 {
2847   assert(semantics == (const llvm::fltSemantics*)&IEEEsingle);
2848   assert(partCount()==1);
2849
2850   uint32_t myexponent, mysignificand;
2851
2852   if (category==fcNormal) {
2853     myexponent = exponent+127; //bias
2854     mysignificand = (uint32_t)*significandParts();
2855     if (myexponent == 1 && !(mysignificand & 0x800000))
2856       myexponent = 0;   // denormal
2857   } else if (category==fcZero) {
2858     myexponent = 0;
2859     mysignificand = 0;
2860   } else if (category==fcInfinity) {
2861     myexponent = 0xff;
2862     mysignificand = 0;
2863   } else {
2864     assert(category == fcNaN && "Unknown category!");
2865     myexponent = 0xff;
2866     mysignificand = (uint32_t)*significandParts();
2867   }
2868
2869   return APInt(32, (((sign&1) << 31) | ((myexponent&0xff) << 23) |
2870                     (mysignificand & 0x7fffff)));
2871 }
2872
2873 APInt
2874 APFloat::convertHalfAPFloatToAPInt() const
2875 {
2876   assert(semantics == (const llvm::fltSemantics*)&IEEEhalf);
2877   assert(partCount()==1);
2878
2879   uint32_t myexponent, mysignificand;
2880
2881   if (category==fcNormal) {
2882     myexponent = exponent+15; //bias
2883     mysignificand = (uint32_t)*significandParts();
2884     if (myexponent == 1 && !(mysignificand & 0x400))
2885       myexponent = 0;   // denormal
2886   } else if (category==fcZero) {
2887     myexponent = 0;
2888     mysignificand = 0;
2889   } else if (category==fcInfinity) {
2890     myexponent = 0x1f;
2891     mysignificand = 0;
2892   } else {
2893     assert(category == fcNaN && "Unknown category!");
2894     myexponent = 0x1f;
2895     mysignificand = (uint32_t)*significandParts();
2896   }
2897
2898   return APInt(16, (((sign&1) << 15) | ((myexponent&0x1f) << 10) |
2899                     (mysignificand & 0x3ff)));
2900 }
2901
2902 // This function creates an APInt that is just a bit map of the floating
2903 // point constant as it would appear in memory.  It is not a conversion,
2904 // and treating the result as a normal integer is unlikely to be useful.
2905
2906 APInt
2907 APFloat::bitcastToAPInt() const
2908 {
2909   if (semantics == (const llvm::fltSemantics*)&IEEEhalf)
2910     return convertHalfAPFloatToAPInt();
2911
2912   if (semantics == (const llvm::fltSemantics*)&IEEEsingle)
2913     return convertFloatAPFloatToAPInt();
2914
2915   if (semantics == (const llvm::fltSemantics*)&IEEEdouble)
2916     return convertDoubleAPFloatToAPInt();
2917
2918   if (semantics == (const llvm::fltSemantics*)&IEEEquad)
2919     return convertQuadrupleAPFloatToAPInt();
2920
2921   if (semantics == (const llvm::fltSemantics*)&PPCDoubleDouble)
2922     return convertPPCDoubleDoubleAPFloatToAPInt();
2923
2924   assert(semantics == (const llvm::fltSemantics*)&x87DoubleExtended &&
2925          "unknown format!");
2926   return convertF80LongDoubleAPFloatToAPInt();
2927 }
2928
2929 float
2930 APFloat::convertToFloat() const
2931 {
2932   assert(semantics == (const llvm::fltSemantics*)&IEEEsingle &&
2933          "Float semantics are not IEEEsingle");
2934   APInt api = bitcastToAPInt();
2935   return api.bitsToFloat();
2936 }
2937
2938 double
2939 APFloat::convertToDouble() const
2940 {
2941   assert(semantics == (const llvm::fltSemantics*)&IEEEdouble &&
2942          "Float semantics are not IEEEdouble");
2943   APInt api = bitcastToAPInt();
2944   return api.bitsToDouble();
2945 }
2946
2947 /// Integer bit is explicit in this format.  Intel hardware (387 and later)
2948 /// does not support these bit patterns:
2949 ///  exponent = all 1's, integer bit 0, significand 0 ("pseudoinfinity")
2950 ///  exponent = all 1's, integer bit 0, significand nonzero ("pseudoNaN")
2951 ///  exponent = 0, integer bit 1 ("pseudodenormal")
2952 ///  exponent!=0 nor all 1's, integer bit 0 ("unnormal")
2953 /// At the moment, the first two are treated as NaNs, the second two as Normal.
2954 void
2955 APFloat::initFromF80LongDoubleAPInt(const APInt &api)
2956 {
2957   assert(api.getBitWidth()==80);
2958   uint64_t i1 = api.getRawData()[0];
2959   uint64_t i2 = api.getRawData()[1];
2960   uint64_t myexponent = (i2 & 0x7fff);
2961   uint64_t mysignificand = i1;
2962
2963   initialize(&APFloat::x87DoubleExtended);
2964   assert(partCount()==2);
2965
2966   sign = static_cast<unsigned int>(i2>>15);
2967   if (myexponent==0 && mysignificand==0) {
2968     // exponent, significand meaningless
2969     category = fcZero;
2970   } else if (myexponent==0x7fff && mysignificand==0x8000000000000000ULL) {
2971     // exponent, significand meaningless
2972     category = fcInfinity;
2973   } else if (myexponent==0x7fff && mysignificand!=0x8000000000000000ULL) {
2974     // exponent meaningless
2975     category = fcNaN;
2976     significandParts()[0] = mysignificand;
2977     significandParts()[1] = 0;
2978   } else {
2979     category = fcNormal;
2980     exponent = myexponent - 16383;
2981     significandParts()[0] = mysignificand;
2982     significandParts()[1] = 0;
2983     if (myexponent==0)          // denormal
2984       exponent = -16382;
2985   }
2986 }
2987
2988 void
2989 APFloat::initFromPPCDoubleDoubleAPInt(const APInt &api)
2990 {
2991   assert(api.getBitWidth()==128);
2992   uint64_t i1 = api.getRawData()[0];
2993   uint64_t i2 = api.getRawData()[1];
2994   uint64_t myexponent = (i1 >> 52) & 0x7ff;
2995   uint64_t mysignificand = i1 & 0xfffffffffffffLL;
2996   uint64_t myexponent2 = (i2 >> 52) & 0x7ff;
2997   uint64_t mysignificand2 = i2 & 0xfffffffffffffLL;
2998
2999   initialize(&APFloat::PPCDoubleDouble);
3000   assert(partCount()==2);
3001
3002   sign = static_cast<unsigned int>(i1>>63);
3003   sign2 = static_cast<unsigned int>(i2>>63);
3004   if (myexponent==0 && mysignificand==0) {
3005     // exponent, significand meaningless
3006     // exponent2 and significand2 are required to be 0; we don't check
3007     category = fcZero;
3008   } else if (myexponent==0x7ff && mysignificand==0) {
3009     // exponent, significand meaningless
3010     // exponent2 and significand2 are required to be 0; we don't check
3011     category = fcInfinity;
3012   } else if (myexponent==0x7ff && mysignificand!=0) {
3013     // exponent meaningless.  So is the whole second word, but keep it
3014     // for determinism.
3015     category = fcNaN;
3016     exponent2 = myexponent2;
3017     significandParts()[0] = mysignificand;
3018     significandParts()[1] = mysignificand2;
3019   } else {
3020     category = fcNormal;
3021     // Note there is no category2; the second word is treated as if it is
3022     // fcNormal, although it might be something else considered by itself.
3023     exponent = myexponent - 1023;
3024     exponent2 = myexponent2 - 1023;
3025     significandParts()[0] = mysignificand;
3026     significandParts()[1] = mysignificand2;
3027     if (myexponent==0)          // denormal
3028       exponent = -1022;
3029     else
3030       significandParts()[0] |= 0x10000000000000LL;  // integer bit
3031     if (myexponent2==0)
3032       exponent2 = -1022;
3033     else
3034       significandParts()[1] |= 0x10000000000000LL;  // integer bit
3035   }
3036 }
3037
3038 void
3039 APFloat::initFromQuadrupleAPInt(const APInt &api)
3040 {
3041   assert(api.getBitWidth()==128);
3042   uint64_t i1 = api.getRawData()[0];
3043   uint64_t i2 = api.getRawData()[1];
3044   uint64_t myexponent = (i2 >> 48) & 0x7fff;
3045   uint64_t mysignificand  = i1;
3046   uint64_t mysignificand2 = i2 & 0xffffffffffffLL;
3047
3048   initialize(&APFloat::IEEEquad);
3049   assert(partCount()==2);
3050
3051   sign = static_cast<unsigned int>(i2>>63);
3052   if (myexponent==0 &&
3053       (mysignificand==0 && mysignificand2==0)) {
3054     // exponent, significand meaningless
3055     category = fcZero;
3056   } else if (myexponent==0x7fff &&
3057              (mysignificand==0 && mysignificand2==0)) {
3058     // exponent, significand meaningless
3059     category = fcInfinity;
3060   } else if (myexponent==0x7fff &&
3061              (mysignificand!=0 || mysignificand2 !=0)) {
3062     // exponent meaningless
3063     category = fcNaN;
3064     significandParts()[0] = mysignificand;
3065     significandParts()[1] = mysignificand2;
3066   } else {
3067     category = fcNormal;
3068     exponent = myexponent - 16383;
3069     significandParts()[0] = mysignificand;
3070     significandParts()[1] = mysignificand2;
3071     if (myexponent==0)          // denormal
3072       exponent = -16382;
3073     else
3074       significandParts()[1] |= 0x1000000000000LL;  // integer bit
3075   }
3076 }
3077
3078 void
3079 APFloat::initFromDoubleAPInt(const APInt &api)
3080 {
3081   assert(api.getBitWidth()==64);
3082   uint64_t i = *api.getRawData();
3083   uint64_t myexponent = (i >> 52) & 0x7ff;
3084   uint64_t mysignificand = i & 0xfffffffffffffLL;
3085
3086   initialize(&APFloat::IEEEdouble);
3087   assert(partCount()==1);
3088
3089   sign = static_cast<unsigned int>(i>>63);
3090   if (myexponent==0 && mysignificand==0) {
3091     // exponent, significand meaningless
3092     category = fcZero;
3093   } else if (myexponent==0x7ff && mysignificand==0) {
3094     // exponent, significand meaningless
3095     category = fcInfinity;
3096   } else if (myexponent==0x7ff && mysignificand!=0) {
3097     // exponent meaningless
3098     category = fcNaN;
3099     *significandParts() = mysignificand;
3100   } else {
3101     category = fcNormal;
3102     exponent = myexponent - 1023;
3103     *significandParts() = mysignificand;
3104     if (myexponent==0)          // denormal
3105       exponent = -1022;
3106     else
3107       *significandParts() |= 0x10000000000000LL;  // integer bit
3108   }
3109 }
3110
3111 void
3112 APFloat::initFromFloatAPInt(const APInt & api)
3113 {
3114   assert(api.getBitWidth()==32);
3115   uint32_t i = (uint32_t)*api.getRawData();
3116   uint32_t myexponent = (i >> 23) & 0xff;
3117   uint32_t mysignificand = i & 0x7fffff;
3118
3119   initialize(&APFloat::IEEEsingle);
3120   assert(partCount()==1);
3121
3122   sign = i >> 31;
3123   if (myexponent==0 && mysignificand==0) {
3124     // exponent, significand meaningless
3125     category = fcZero;
3126   } else if (myexponent==0xff && mysignificand==0) {
3127     // exponent, significand meaningless
3128     category = fcInfinity;
3129   } else if (myexponent==0xff && mysignificand!=0) {
3130     // sign, exponent, significand meaningless
3131     category = fcNaN;
3132     *significandParts() = mysignificand;
3133   } else {
3134     category = fcNormal;
3135     exponent = myexponent - 127;  //bias
3136     *significandParts() = mysignificand;
3137     if (myexponent==0)    // denormal
3138       exponent = -126;
3139     else
3140       *significandParts() |= 0x800000; // integer bit
3141   }
3142 }
3143
3144 void
3145 APFloat::initFromHalfAPInt(const APInt & api)
3146 {
3147   assert(api.getBitWidth()==16);
3148   uint32_t i = (uint32_t)*api.getRawData();
3149   uint32_t myexponent = (i >> 10) & 0x1f;
3150   uint32_t mysignificand = i & 0x3ff;
3151
3152   initialize(&APFloat::IEEEhalf);
3153   assert(partCount()==1);
3154
3155   sign = i >> 15;
3156   if (myexponent==0 && mysignificand==0) {
3157     // exponent, significand meaningless
3158     category = fcZero;
3159   } else if (myexponent==0x1f && mysignificand==0) {
3160     // exponent, significand meaningless
3161     category = fcInfinity;
3162   } else if (myexponent==0x1f && mysignificand!=0) {
3163     // sign, exponent, significand meaningless
3164     category = fcNaN;
3165     *significandParts() = mysignificand;
3166   } else {
3167     category = fcNormal;
3168     exponent = myexponent - 15;  //bias
3169     *significandParts() = mysignificand;
3170     if (myexponent==0)    // denormal
3171       exponent = -14;
3172     else
3173       *significandParts() |= 0x400; // integer bit
3174   }
3175 }
3176
3177 /// Treat api as containing the bits of a floating point number.  Currently
3178 /// we infer the floating point type from the size of the APInt.  The
3179 /// isIEEE argument distinguishes between PPC128 and IEEE128 (not meaningful
3180 /// when the size is anything else).
3181 void
3182 APFloat::initFromAPInt(const APInt& api, bool isIEEE)
3183 {
3184   if (api.getBitWidth() == 16)
3185     return initFromHalfAPInt(api);
3186   else if (api.getBitWidth() == 32)
3187     return initFromFloatAPInt(api);
3188   else if (api.getBitWidth()==64)
3189     return initFromDoubleAPInt(api);
3190   else if (api.getBitWidth()==80)
3191     return initFromF80LongDoubleAPInt(api);
3192   else if (api.getBitWidth()==128)
3193     return (isIEEE ?
3194             initFromQuadrupleAPInt(api) : initFromPPCDoubleDoubleAPInt(api));
3195   else
3196     llvm_unreachable(0);
3197 }
3198
3199 APFloat APFloat::getLargest(const fltSemantics &Sem, bool Negative) {
3200   APFloat Val(Sem, fcNormal, Negative);
3201
3202   // We want (in interchange format):
3203   //   sign = {Negative}
3204   //   exponent = 1..10
3205   //   significand = 1..1
3206
3207   Val.exponent = Sem.maxExponent; // unbiased
3208
3209   // 1-initialize all bits....
3210   Val.zeroSignificand();
3211   integerPart *significand = Val.significandParts();
3212   unsigned N = partCountForBits(Sem.precision);
3213   for (unsigned i = 0; i != N; ++i)
3214     significand[i] = ~((integerPart) 0);
3215
3216   // ...and then clear the top bits for internal consistency.
3217   significand[N-1] &=
3218     (((integerPart) 1) << ((Sem.precision % integerPartWidth) - 1)) - 1;
3219
3220   return Val;
3221 }
3222
3223 APFloat APFloat::getSmallest(const fltSemantics &Sem, bool Negative) {
3224   APFloat Val(Sem, fcNormal, Negative);
3225
3226   // We want (in interchange format):
3227   //   sign = {Negative}
3228   //   exponent = 0..0
3229   //   significand = 0..01
3230
3231   Val.exponent = Sem.minExponent; // unbiased
3232   Val.zeroSignificand();
3233   Val.significandParts()[0] = 1;
3234   return Val;
3235 }
3236
3237 APFloat APFloat::getSmallestNormalized(const fltSemantics &Sem, bool Negative) {
3238   APFloat Val(Sem, fcNormal, Negative);
3239
3240   // We want (in interchange format):
3241   //   sign = {Negative}
3242   //   exponent = 0..0
3243   //   significand = 10..0
3244
3245   Val.exponent = Sem.minExponent;
3246   Val.zeroSignificand();
3247   Val.significandParts()[partCountForBits(Sem.precision)-1] |=
3248     (((integerPart) 1) << ((Sem.precision % integerPartWidth) - 1));
3249
3250   return Val;
3251 }
3252
3253 APFloat::APFloat(const APInt& api, bool isIEEE)
3254 {
3255   initFromAPInt(api, isIEEE);
3256 }
3257
3258 APFloat::APFloat(float f)
3259 {
3260   APInt api = APInt(32, 0);
3261   initFromAPInt(api.floatToBits(f));
3262 }
3263
3264 APFloat::APFloat(double d)
3265 {
3266   APInt api = APInt(64, 0);
3267   initFromAPInt(api.doubleToBits(d));
3268 }
3269
3270 namespace {
3271   static void append(SmallVectorImpl<char> &Buffer,
3272                      unsigned N, const char *Str) {
3273     unsigned Start = Buffer.size();
3274     Buffer.set_size(Start + N);
3275     memcpy(&Buffer[Start], Str, N);
3276   }
3277
3278   template <unsigned N>
3279   void append(SmallVectorImpl<char> &Buffer, const char (&Str)[N]) {
3280     append(Buffer, N, Str);
3281   }
3282
3283   /// Removes data from the given significand until it is no more
3284   /// precise than is required for the desired precision.
3285   void AdjustToPrecision(APInt &significand,
3286                          int &exp, unsigned FormatPrecision) {
3287     unsigned bits = significand.getActiveBits();
3288
3289     // 196/59 is a very slight overestimate of lg_2(10).
3290     unsigned bitsRequired = (FormatPrecision * 196 + 58) / 59;
3291
3292     if (bits <= bitsRequired) return;
3293
3294     unsigned tensRemovable = (bits - bitsRequired) * 59 / 196;
3295     if (!tensRemovable) return;
3296
3297     exp += tensRemovable;
3298
3299     APInt divisor(significand.getBitWidth(), 1);
3300     APInt powten(significand.getBitWidth(), 10);
3301     while (true) {
3302       if (tensRemovable & 1)
3303         divisor *= powten;
3304       tensRemovable >>= 1;
3305       if (!tensRemovable) break;
3306       powten *= powten;
3307     }
3308
3309     significand = significand.udiv(divisor);
3310
3311     // Truncate the significand down to its active bit count, but
3312     // don't try to drop below 32.
3313     unsigned newPrecision = std::max(32U, significand.getActiveBits());
3314     significand.trunc(newPrecision);
3315   }
3316
3317
3318   void AdjustToPrecision(SmallVectorImpl<char> &buffer,
3319                          int &exp, unsigned FormatPrecision) {
3320     unsigned N = buffer.size();
3321     if (N <= FormatPrecision) return;
3322
3323     // The most significant figures are the last ones in the buffer.
3324     unsigned FirstSignificant = N - FormatPrecision;
3325
3326     // Round.
3327     // FIXME: this probably shouldn't use 'round half up'.
3328
3329     // Rounding down is just a truncation, except we also want to drop
3330     // trailing zeros from the new result.
3331     if (buffer[FirstSignificant - 1] < '5') {
3332       while (buffer[FirstSignificant] == '0')
3333         FirstSignificant++;
3334
3335       exp += FirstSignificant;
3336       buffer.erase(&buffer[0], &buffer[FirstSignificant]);
3337       return;
3338     }
3339
3340     // Rounding up requires a decimal add-with-carry.  If we continue
3341     // the carry, the newly-introduced zeros will just be truncated.
3342     for (unsigned I = FirstSignificant; I != N; ++I) {
3343       if (buffer[I] == '9') {
3344         FirstSignificant++;
3345       } else {
3346         buffer[I]++;
3347         break;
3348       }
3349     }
3350
3351     // If we carried through, we have exactly one digit of precision.
3352     if (FirstSignificant == N) {
3353       exp += FirstSignificant;
3354       buffer.clear();
3355       buffer.push_back('1');
3356       return;
3357     }
3358
3359     exp += FirstSignificant;
3360     buffer.erase(&buffer[0], &buffer[FirstSignificant]);
3361   }
3362 }
3363
3364 void APFloat::toString(SmallVectorImpl<char> &Str,
3365                        unsigned FormatPrecision,
3366                        unsigned FormatMaxPadding) const {
3367   switch (category) {
3368   case fcInfinity:
3369     if (isNegative())
3370       return append(Str, "-Inf");
3371     else
3372       return append(Str, "+Inf");
3373
3374   case fcNaN: return append(Str, "NaN");
3375
3376   case fcZero:
3377     if (isNegative())
3378       Str.push_back('-');
3379
3380     if (!FormatMaxPadding)
3381       append(Str, "0.0E+0");
3382     else
3383       Str.push_back('0');
3384     return;
3385
3386   case fcNormal:
3387     break;
3388   }
3389
3390   if (isNegative())
3391     Str.push_back('-');
3392
3393   // Decompose the number into an APInt and an exponent.
3394   int exp = exponent - ((int) semantics->precision - 1);
3395   APInt significand(semantics->precision,
3396                     partCountForBits(semantics->precision),
3397                     significandParts());
3398
3399   // Set FormatPrecision if zero.  We want to do this before we
3400   // truncate trailing zeros, as those are part of the precision.
3401   if (!FormatPrecision) {
3402     // It's an interesting question whether to use the nominal
3403     // precision or the active precision here for denormals.
3404
3405     // FormatPrecision = ceil(significandBits / lg_2(10))
3406     FormatPrecision = (semantics->precision * 59 + 195) / 196;
3407   }
3408
3409   // Ignore trailing binary zeros.
3410   int trailingZeros = significand.countTrailingZeros();
3411   exp += trailingZeros;
3412   significand = significand.lshr(trailingZeros);
3413
3414   // Change the exponent from 2^e to 10^e.
3415   if (exp == 0) {
3416     // Nothing to do.
3417   } else if (exp > 0) {
3418     // Just shift left.
3419     significand.zext(semantics->precision + exp);
3420     significand <<= exp;
3421     exp = 0;
3422   } else { /* exp < 0 */
3423     int texp = -exp;
3424
3425     // We transform this using the identity:
3426     //   (N)(2^-e) == (N)(5^e)(10^-e)
3427     // This means we have to multiply N (the significand) by 5^e.
3428     // To avoid overflow, we have to operate on numbers large
3429     // enough to store N * 5^e:
3430     //   log2(N * 5^e) == log2(N) + e * log2(5)
3431     //                 <= semantics->precision + e * 137 / 59
3432     //   (log_2(5) ~ 2.321928 < 2.322034 ~ 137/59)
3433
3434     unsigned precision = semantics->precision + 137 * texp / 59;
3435
3436     // Multiply significand by 5^e.
3437     //   N * 5^0101 == N * 5^(1*1) * 5^(0*2) * 5^(1*4) * 5^(0*8)
3438     significand.zext(precision);
3439     APInt five_to_the_i(precision, 5);
3440     while (true) {
3441       if (texp & 1) significand *= five_to_the_i;
3442
3443       texp >>= 1;
3444       if (!texp) break;
3445       five_to_the_i *= five_to_the_i;
3446     }
3447   }
3448
3449   AdjustToPrecision(significand, exp, FormatPrecision);
3450
3451   llvm::SmallVector<char, 256> buffer;
3452
3453   // Fill the buffer.
3454   unsigned precision = significand.getBitWidth();
3455   APInt ten(precision, 10);
3456   APInt digit(precision, 0);
3457
3458   bool inTrail = true;
3459   while (significand != 0) {
3460     // digit <- significand % 10
3461     // significand <- significand / 10
3462     APInt::udivrem(significand, ten, significand, digit);
3463
3464     unsigned d = digit.getZExtValue();
3465
3466     // Drop trailing zeros.
3467     if (inTrail && !d) exp++;
3468     else {
3469       buffer.push_back((char) ('0' + d));
3470       inTrail = false;
3471     }
3472   }
3473
3474   assert(!buffer.empty() && "no characters in buffer!");
3475
3476   // Drop down to FormatPrecision.
3477   // TODO: don't do more precise calculations above than are required.
3478   AdjustToPrecision(buffer, exp, FormatPrecision);
3479
3480   unsigned NDigits = buffer.size();
3481
3482   // Check whether we should use scientific notation.
3483   bool FormatScientific;
3484   if (!FormatMaxPadding)
3485     FormatScientific = true;
3486   else {
3487     if (exp >= 0) {
3488       // 765e3 --> 765000
3489       //              ^^^
3490       // But we shouldn't make the number look more precise than it is.
3491       FormatScientific = ((unsigned) exp > FormatMaxPadding ||
3492                           NDigits + (unsigned) exp > FormatPrecision);
3493     } else {
3494       // Power of the most significant digit.
3495       int MSD = exp + (int) (NDigits - 1);
3496       if (MSD >= 0) {
3497         // 765e-2 == 7.65
3498         FormatScientific = false;
3499       } else {
3500         // 765e-5 == 0.00765
3501         //           ^ ^^
3502         FormatScientific = ((unsigned) -MSD) > FormatMaxPadding;
3503       }
3504     }
3505   }
3506
3507   // Scientific formatting is pretty straightforward.
3508   if (FormatScientific) {
3509     exp += (NDigits - 1);
3510
3511     Str.push_back(buffer[NDigits-1]);
3512     Str.push_back('.');
3513     if (NDigits == 1)
3514       Str.push_back('0');
3515     else
3516       for (unsigned I = 1; I != NDigits; ++I)
3517         Str.push_back(buffer[NDigits-1-I]);
3518     Str.push_back('E');
3519
3520     Str.push_back(exp >= 0 ? '+' : '-');
3521     if (exp < 0) exp = -exp;
3522     SmallVector<char, 6> expbuf;
3523     do {
3524       expbuf.push_back((char) ('0' + (exp % 10)));
3525       exp /= 10;
3526     } while (exp);
3527     for (unsigned I = 0, E = expbuf.size(); I != E; ++I)
3528       Str.push_back(expbuf[E-1-I]);
3529     return;
3530   }
3531
3532   // Non-scientific, positive exponents.
3533   if (exp >= 0) {
3534     for (unsigned I = 0; I != NDigits; ++I)
3535       Str.push_back(buffer[NDigits-1-I]);
3536     for (unsigned I = 0; I != (unsigned) exp; ++I)
3537       Str.push_back('0');
3538     return;
3539   }
3540
3541   // Non-scientific, negative exponents.
3542
3543   // The number of digits to the left of the decimal point.
3544   int NWholeDigits = exp + (int) NDigits;
3545
3546   unsigned I = 0;
3547   if (NWholeDigits > 0) {
3548     for (; I != (unsigned) NWholeDigits; ++I)
3549       Str.push_back(buffer[NDigits-I-1]);
3550     Str.push_back('.');
3551   } else {
3552     unsigned NZeros = 1 + (unsigned) -NWholeDigits;
3553
3554     Str.push_back('0');
3555     Str.push_back('.');
3556     for (unsigned Z = 1; Z != NZeros; ++Z)
3557       Str.push_back('0');
3558   }
3559
3560   for (; I != NDigits; ++I)
3561     Str.push_back(buffer[NDigits-I-1]);
3562 }