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