lib/Support/APFloat.cpp

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