lib/Support/APFloat.cpp

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