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