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