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