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