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