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