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