APFloat.cpp revision a471c2ecda37cd1bae0d94e832f002caa7b63216 
1//===-- APFloat.cpp - Implement APFloat class -----------------------------===//
2//
3//                     The LLVM Compiler Infrastructure
4//
5// This file was developed by Neil Booth and is distributed under the
6// University of Illinois Open Source License. See LICENSE.TXT for details.
7//
8//===----------------------------------------------------------------------===//
9//
10// This file implements a class to represent arbitrary precision floating
11// point values and provide a variety of arithmetic operations on them.
12//
13//===----------------------------------------------------------------------===//
14
15#include <cassert>
16#include <cstring>
17#include "llvm/ADT/APFloat.h"
18#include "llvm/Support/MathExtras.h"
19
20using namespace llvm;
21
22#define convolve(lhs, rhs) ((lhs) * 4 + (rhs))
23
24/* Assumed in hexadecimal significand parsing, and conversion to
25   hexadecimal strings.  */
26COMPILE_TIME_ASSERT(integerPartWidth % 4 == 0);
27
28namespace llvm {
29
30  /* Represents floating point arithmetic semantics.  */
31  struct fltSemantics {
32    /* The largest E such that 2^E is representable; this matches the
33       definition of IEEE 754.  */
34    exponent_t maxExponent;
35
36    /* The smallest E such that 2^E is a normalized number; this
37       matches the definition of IEEE 754.  */
38    exponent_t minExponent;
39
40    /* Number of bits in the significand.  This includes the integer
41       bit.  */
42    unsigned char precision;
43
44    /* If the target format has an implicit integer bit.  */
45    bool implicitIntegerBit;
46  };
47
48  const fltSemantics APFloat::IEEEsingle = { 127, -126, 24, true };
49  const fltSemantics APFloat::IEEEdouble = { 1023, -1022, 53, true };
50  const fltSemantics APFloat::IEEEquad = { 16383, -16382, 113, true };
51  const fltSemantics APFloat::x87DoubleExtended = { 16383, -16382, 64, false };
52  const fltSemantics APFloat::Bogus = { 0, 0, 0, false };
53
54  // The PowerPC format consists of two doubles.  It does not map cleanly
55  // onto the usual format above.  For now only storage of constants of
56  // this type is supported, no arithmetic.
57  const fltSemantics APFloat::PPCDoubleDouble = { 1023, -1022, 106, true };
58}
59
60/* Put a bunch of private, handy routines in an anonymous namespace.  */
61namespace {
62
63  inline unsigned int
64  partCountForBits(unsigned int bits)
65  {
66    return ((bits) + integerPartWidth - 1) / integerPartWidth;
67  }
68
69  unsigned int
70  digitValue(unsigned int c)
71  {
72    unsigned int r;
73
74    r = c - '0';
75    if(r <= 9)
76      return r;
77
78    return -1U;
79  }
80
81  unsigned int
82  hexDigitValue (unsigned int c)
83  {
84    unsigned int r;
85
86    r = c - '0';
87    if(r <= 9)
88      return r;
89
90    r = c - 'A';
91    if(r <= 5)
92      return r + 10;
93
94    r = c - 'a';
95    if(r <= 5)
96      return r + 10;
97
98    return -1U;
99  }
100
101  /* This is ugly and needs cleaning up, but I don't immediately see
102     how whilst remaining safe.  */
103  static int
104  totalExponent(const char *p, int exponentAdjustment)
105  {
106    integerPart unsignedExponent;
107    bool negative, overflow;
108    long exponent;
109
110    /* Move past the exponent letter and sign to the digits.  */
111    p++;
112    negative = *p == '-';
113    if(*p == '-' || *p == '+')
114      p++;
115
116    unsignedExponent = 0;
117    overflow = false;
118    for(;;) {
119      unsigned int value;
120
121      value = digitValue(*p);
122      if(value == -1U)
123        break;
124
125      p++;
126      unsignedExponent = unsignedExponent * 10 + value;
127      if(unsignedExponent > 65535)
128        overflow = true;
129    }
130
131    if(exponentAdjustment > 65535 || exponentAdjustment < -65536)
132      overflow = true;
133
134    if(!overflow) {
135      exponent = unsignedExponent;
136      if(negative)
137        exponent = -exponent;
138      exponent += exponentAdjustment;
139      if(exponent > 65535 || exponent < -65536)
140        overflow = true;
141    }
142
143    if(overflow)
144      exponent = negative ? -65536: 65535;
145
146    return exponent;
147  }
148
149  const char *
150  skipLeadingZeroesAndAnyDot(const char *p, const char **dot)
151  {
152    *dot = 0;
153    while(*p == '0')
154      p++;
155
156    if(*p == '.') {
157      *dot = p++;
158      while(*p == '0')
159        p++;
160    }
161
162    return p;
163  }
164
165  /* Return the trailing fraction of a hexadecimal number.
166     DIGITVALUE is the first hex digit of the fraction, P points to
167     the next digit.  */
168  lostFraction
169  trailingHexadecimalFraction(const char *p, unsigned int digitValue)
170  {
171    unsigned int hexDigit;
172
173    /* If the first trailing digit isn't 0 or 8 we can work out the
174       fraction immediately.  */
175    if(digitValue > 8)
176      return lfMoreThanHalf;
177    else if(digitValue < 8 && digitValue > 0)
178      return lfLessThanHalf;
179
180    /* Otherwise we need to find the first non-zero digit.  */
181    while(*p == '0')
182      p++;
183
184    hexDigit = hexDigitValue(*p);
185
186    /* If we ran off the end it is exactly zero or one-half, otherwise
187       a little more.  */
188    if(hexDigit == -1U)
189      return digitValue == 0 ? lfExactlyZero: lfExactlyHalf;
190    else
191      return digitValue == 0 ? lfLessThanHalf: lfMoreThanHalf;
192  }
193
194  /* Return the fraction lost were a bignum truncated losing the least
195     significant BITS bits.  */
196  lostFraction
197  lostFractionThroughTruncation(const integerPart *parts,
198                                unsigned int partCount,
199                                unsigned int bits)
200  {
201    unsigned int lsb;
202
203    lsb = APInt::tcLSB(parts, partCount);
204
205    /* Note this is guaranteed true if bits == 0, or LSB == -1U.  */
206    if(bits <= lsb)
207      return lfExactlyZero;
208    if(bits == lsb + 1)
209      return lfExactlyHalf;
210    if(bits <= partCount * integerPartWidth
211       && APInt::tcExtractBit(parts, bits - 1))
212      return lfMoreThanHalf;
213
214    return lfLessThanHalf;
215  }
216
217  /* Shift DST right BITS bits noting lost fraction.  */
218  lostFraction
219  shiftRight(integerPart *dst, unsigned int parts, unsigned int bits)
220  {
221    lostFraction lost_fraction;
222
223    lost_fraction = lostFractionThroughTruncation(dst, parts, bits);
224
225    APInt::tcShiftRight(dst, parts, bits);
226
227    return lost_fraction;
228  }
229
230  /* Combine the effect of two lost fractions.  */
231  lostFraction
232  combineLostFractions(lostFraction moreSignificant,
233                       lostFraction lessSignificant)
234  {
235    if(lessSignificant != lfExactlyZero) {
236      if(moreSignificant == lfExactlyZero)
237        moreSignificant = lfLessThanHalf;
238      else if(moreSignificant == lfExactlyHalf)
239        moreSignificant = lfMoreThanHalf;
240    }
241
242    return moreSignificant;
243  }
244
245  /* Zero at the end to avoid modular arithmetic when adding one; used
246     when rounding up during hexadecimal output.  */
247  static const char hexDigitsLower[] = "0123456789abcdef0";
248  static const char hexDigitsUpper[] = "0123456789ABCDEF0";
249  static const char infinityL[] = "infinity";
250  static const char infinityU[] = "INFINITY";
251  static const char NaNL[] = "nan";
252  static const char NaNU[] = "NAN";
253
254  /* Write out an integerPart in hexadecimal, starting with the most
255     significant nibble.  Write out exactly COUNT hexdigits, return
256     COUNT.  */
257  static unsigned int
258  partAsHex (char *dst, integerPart part, unsigned int count,
259             const char *hexDigitChars)
260  {
261    unsigned int result = count;
262
263    assert (count != 0 && count <= integerPartWidth / 4);
264
265    part >>= (integerPartWidth - 4 * count);
266    while (count--) {
267      dst[count] = hexDigitChars[part & 0xf];
268      part >>= 4;
269    }
270
271    return result;
272  }
273
274  /* Write out an unsigned decimal integer.  */
275  static char *
276  writeUnsignedDecimal (char *dst, unsigned int n)
277  {
278    char buff[40], *p;
279
280    p = buff;
281    do
282      *p++ = '0' + n % 10;
283    while (n /= 10);
284
285    do
286      *dst++ = *--p;
287    while (p != buff);
288
289    return dst;
290  }
291
292  /* Write out a signed decimal integer.  */
293  static char *
294  writeSignedDecimal (char *dst, int value)
295  {
296    if (value < 0) {
297      *dst++ = '-';
298      dst = writeUnsignedDecimal(dst, -(unsigned) value);
299    } else
300      dst = writeUnsignedDecimal(dst, value);
301
302    return dst;
303  }
304}
305
306/* Constructors.  */
307void
308APFloat::initialize(const fltSemantics *ourSemantics)
309{
310  unsigned int count;
311
312  semantics = ourSemantics;
313  count = partCount();
314  if(count > 1)
315    significand.parts = new integerPart[count];
316}
317
318void
319APFloat::freeSignificand()
320{
321  if(partCount() > 1)
322    delete [] significand.parts;
323}
324
325void
326APFloat::assign(const APFloat &rhs)
327{
328  assert(semantics == rhs.semantics);
329
330  sign = rhs.sign;
331  category = rhs.category;
332  exponent = rhs.exponent;
333  sign2 = rhs.sign2;
334  exponent2 = rhs.exponent2;
335  if(category == fcNormal || category == fcNaN)
336    copySignificand(rhs);
337}
338
339void
340APFloat::copySignificand(const APFloat &rhs)
341{
342  assert(category == fcNormal || category == fcNaN);
343  assert(rhs.partCount() >= partCount());
344
345  APInt::tcAssign(significandParts(), rhs.significandParts(),
346                  partCount());
347}
348
349APFloat &
350APFloat::operator=(const APFloat &rhs)
351{
352  if(this != &rhs) {
353    if(semantics != rhs.semantics) {
354      freeSignificand();
355      initialize(rhs.semantics);
356    }
357    assign(rhs);
358  }
359
360  return *this;
361}
362
363bool
364APFloat::bitwiseIsEqual(const APFloat &rhs) const {
365  if (this == &rhs)
366    return true;
367  if (semantics != rhs.semantics ||
368      category != rhs.category ||
369      sign != rhs.sign)
370    return false;
371  if (semantics==(const llvm::fltSemantics* const)&PPCDoubleDouble &&
372      sign2 != rhs.sign2)
373    return false;
374  if (category==fcZero || category==fcInfinity)
375    return true;
376  else if (category==fcNormal && exponent!=rhs.exponent)
377    return false;
378  else if (semantics==(const llvm::fltSemantics* const)&PPCDoubleDouble &&
379           exponent2!=rhs.exponent2)
380    return false;
381  else {
382    int i= partCount();
383    const integerPart* p=significandParts();
384    const integerPart* q=rhs.significandParts();
385    for (; i>0; i--, p++, q++) {
386      if (*p != *q)
387        return false;
388    }
389    return true;
390  }
391}
392
393APFloat::APFloat(const fltSemantics &ourSemantics, integerPart value)
394{
395  assert(semantics != (const llvm::fltSemantics* const)&PPCDoubleDouble &&
396         "Compile-time arithmetic on PPC long double not supported yet");
397  initialize(&ourSemantics);
398  sign = 0;
399  zeroSignificand();
400  exponent = ourSemantics.precision - 1;
401  significandParts()[0] = value;
402  normalize(rmNearestTiesToEven, lfExactlyZero);
403}
404
405APFloat::APFloat(const fltSemantics &ourSemantics,
406                 fltCategory ourCategory, bool negative)
407{
408  assert(semantics != (const llvm::fltSemantics* const)&PPCDoubleDouble &&
409         "Compile-time arithmetic on PPC long double not supported yet");
410  initialize(&ourSemantics);
411  category = ourCategory;
412  sign = negative;
413  if(category == fcNormal)
414    category = fcZero;
415}
416
417APFloat::APFloat(const fltSemantics &ourSemantics, const char *text)
418{
419  assert(semantics != (const llvm::fltSemantics* const)&PPCDoubleDouble &&
420         "Compile-time arithmetic on PPC long double not supported yet");
421  initialize(&ourSemantics);
422  convertFromString(text, rmNearestTiesToEven);
423}
424
425APFloat::APFloat(const APFloat &rhs)
426{
427  initialize(rhs.semantics);
428  assign(rhs);
429}
430
431APFloat::~APFloat()
432{
433  freeSignificand();
434}
435
436unsigned int
437APFloat::partCount() const
438{
439  return partCountForBits(semantics->precision + 1);
440}
441
442unsigned int
443APFloat::semanticsPrecision(const fltSemantics &semantics)
444{
445  return semantics.precision;
446}
447
448const integerPart *
449APFloat::significandParts() const
450{
451  return const_cast<APFloat *>(this)->significandParts();
452}
453
454integerPart *
455APFloat::significandParts()
456{
457  assert(category == fcNormal || category == fcNaN);
458
459  if(partCount() > 1)
460    return significand.parts;
461  else
462    return &significand.part;
463}
464
465void
466APFloat::zeroSignificand()
467{
468  category = fcNormal;
469  APInt::tcSet(significandParts(), 0, partCount());
470}
471
472/* Increment an fcNormal floating point number's significand.  */
473void
474APFloat::incrementSignificand()
475{
476  integerPart carry;
477
478  carry = APInt::tcIncrement(significandParts(), partCount());
479
480  /* Our callers should never cause us to overflow.  */
481  assert(carry == 0);
482}
483
484/* Add the significand of the RHS.  Returns the carry flag.  */
485integerPart
486APFloat::addSignificand(const APFloat &rhs)
487{
488  integerPart *parts;
489
490  parts = significandParts();
491
492  assert(semantics == rhs.semantics);
493  assert(exponent == rhs.exponent);
494
495  return APInt::tcAdd(parts, rhs.significandParts(), 0, partCount());
496}
497
498/* Subtract the significand of the RHS with a borrow flag.  Returns
499   the borrow flag.  */
500integerPart
501APFloat::subtractSignificand(const APFloat &rhs, integerPart borrow)
502{
503  integerPart *parts;
504
505  parts = significandParts();
506
507  assert(semantics == rhs.semantics);
508  assert(exponent == rhs.exponent);
509
510  return APInt::tcSubtract(parts, rhs.significandParts(), borrow,
511                           partCount());
512}
513
514/* Multiply the significand of the RHS.  If ADDEND is non-NULL, add it
515   on to the full-precision result of the multiplication.  Returns the
516   lost fraction.  */
517lostFraction
518APFloat::multiplySignificand(const APFloat &rhs, const APFloat *addend)
519{
520  unsigned int omsb;        // One, not zero, based MSB.
521  unsigned int partsCount, newPartsCount, precision;
522  integerPart *lhsSignificand;
523  integerPart scratch[4];
524  integerPart *fullSignificand;
525  lostFraction lost_fraction;
526
527  assert(semantics == rhs.semantics);
528
529  precision = semantics->precision;
530  newPartsCount = partCountForBits(precision * 2);
531
532  if(newPartsCount > 4)
533    fullSignificand = new integerPart[newPartsCount];
534  else
535    fullSignificand = scratch;
536
537  lhsSignificand = significandParts();
538  partsCount = partCount();
539
540  APInt::tcFullMultiply(fullSignificand, lhsSignificand,
541                        rhs.significandParts(), partsCount, partsCount);
542
543  lost_fraction = lfExactlyZero;
544  omsb = APInt::tcMSB(fullSignificand, newPartsCount) + 1;
545  exponent += rhs.exponent;
546
547  if(addend) {
548    Significand savedSignificand = significand;
549    const fltSemantics *savedSemantics = semantics;
550    fltSemantics extendedSemantics;
551    opStatus status;
552    unsigned int extendedPrecision;
553
554    /* Normalize our MSB.  */
555    extendedPrecision = precision + precision - 1;
556    if(omsb != extendedPrecision)
557      {
558        APInt::tcShiftLeft(fullSignificand, newPartsCount,
559                           extendedPrecision - omsb);
560        exponent -= extendedPrecision - omsb;
561      }
562
563    /* Create new semantics.  */
564    extendedSemantics = *semantics;
565    extendedSemantics.precision = extendedPrecision;
566
567    if(newPartsCount == 1)
568      significand.part = fullSignificand[0];
569    else
570      significand.parts = fullSignificand;
571    semantics = &extendedSemantics;
572
573    APFloat extendedAddend(*addend);
574    status = extendedAddend.convert(extendedSemantics, rmTowardZero);
575    assert(status == opOK);
576    lost_fraction = addOrSubtractSignificand(extendedAddend, false);
577
578    /* Restore our state.  */
579    if(newPartsCount == 1)
580      fullSignificand[0] = significand.part;
581    significand = savedSignificand;
582    semantics = savedSemantics;
583
584    omsb = APInt::tcMSB(fullSignificand, newPartsCount) + 1;
585  }
586
587  exponent -= (precision - 1);
588
589  if(omsb > precision) {
590    unsigned int bits, significantParts;
591    lostFraction lf;
592
593    bits = omsb - precision;
594    significantParts = partCountForBits(omsb);
595    lf = shiftRight(fullSignificand, significantParts, bits);
596    lost_fraction = combineLostFractions(lf, lost_fraction);
597    exponent += bits;
598  }
599
600  APInt::tcAssign(lhsSignificand, fullSignificand, partsCount);
601
602  if(newPartsCount > 4)
603    delete [] fullSignificand;
604
605  return lost_fraction;
606}
607
608/* Multiply the significands of LHS and RHS to DST.  */
609lostFraction
610APFloat::divideSignificand(const APFloat &rhs)
611{
612  unsigned int bit, i, partsCount;
613  const integerPart *rhsSignificand;
614  integerPart *lhsSignificand, *dividend, *divisor;
615  integerPart scratch[4];
616  lostFraction lost_fraction;
617
618  assert(semantics == rhs.semantics);
619
620  lhsSignificand = significandParts();
621  rhsSignificand = rhs.significandParts();
622  partsCount = partCount();
623
624  if(partsCount > 2)
625    dividend = new integerPart[partsCount * 2];
626  else
627    dividend = scratch;
628
629  divisor = dividend + partsCount;
630
631  /* Copy the dividend and divisor as they will be modified in-place.  */
632  for(i = 0; i < partsCount; i++) {
633    dividend[i] = lhsSignificand[i];
634    divisor[i] = rhsSignificand[i];
635    lhsSignificand[i] = 0;
636  }
637
638  exponent -= rhs.exponent;
639
640  unsigned int precision = semantics->precision;
641
642  /* Normalize the divisor.  */
643  bit = precision - APInt::tcMSB(divisor, partsCount) - 1;
644  if(bit) {
645    exponent += bit;
646    APInt::tcShiftLeft(divisor, partsCount, bit);
647  }
648
649  /* Normalize the dividend.  */
650  bit = precision - APInt::tcMSB(dividend, partsCount) - 1;
651  if(bit) {
652    exponent -= bit;
653    APInt::tcShiftLeft(dividend, partsCount, bit);
654  }
655
656  if(APInt::tcCompare(dividend, divisor, partsCount) < 0) {
657    exponent--;
658    APInt::tcShiftLeft(dividend, partsCount, 1);
659    assert(APInt::tcCompare(dividend, divisor, partsCount) >= 0);
660  }
661
662  /* Long division.  */
663  for(bit = precision; bit; bit -= 1) {
664    if(APInt::tcCompare(dividend, divisor, partsCount) >= 0) {
665      APInt::tcSubtract(dividend, divisor, 0, partsCount);
666      APInt::tcSetBit(lhsSignificand, bit - 1);
667    }
668
669    APInt::tcShiftLeft(dividend, partsCount, 1);
670  }
671
672  /* Figure out the lost fraction.  */
673  int cmp = APInt::tcCompare(dividend, divisor, partsCount);
674
675  if(cmp > 0)
676    lost_fraction = lfMoreThanHalf;
677  else if(cmp == 0)
678    lost_fraction = lfExactlyHalf;
679  else if(APInt::tcIsZero(dividend, partsCount))
680    lost_fraction = lfExactlyZero;
681  else
682    lost_fraction = lfLessThanHalf;
683
684  if(partsCount > 2)
685    delete [] dividend;
686
687  return lost_fraction;
688}
689
690unsigned int
691APFloat::significandMSB() const
692{
693  return APInt::tcMSB(significandParts(), partCount());
694}
695
696unsigned int
697APFloat::significandLSB() const
698{
699  return APInt::tcLSB(significandParts(), partCount());
700}
701
702/* Note that a zero result is NOT normalized to fcZero.  */
703lostFraction
704APFloat::shiftSignificandRight(unsigned int bits)
705{
706  /* Our exponent should not overflow.  */
707  assert((exponent_t) (exponent + bits) >= exponent);
708
709  exponent += bits;
710
711  return shiftRight(significandParts(), partCount(), bits);
712}
713
714/* Shift the significand left BITS bits, subtract BITS from its exponent.  */
715void
716APFloat::shiftSignificandLeft(unsigned int bits)
717{
718  assert(bits < semantics->precision);
719
720  if(bits) {
721    unsigned int partsCount = partCount();
722
723    APInt::tcShiftLeft(significandParts(), partsCount, bits);
724    exponent -= bits;
725
726    assert(!APInt::tcIsZero(significandParts(), partsCount));
727  }
728}
729
730APFloat::cmpResult
731APFloat::compareAbsoluteValue(const APFloat &rhs) const
732{
733  int compare;
734
735  assert(semantics == rhs.semantics);
736  assert(category == fcNormal);
737  assert(rhs.category == fcNormal);
738
739  compare = exponent - rhs.exponent;
740
741  /* If exponents are equal, do an unsigned bignum comparison of the
742     significands.  */
743  if(compare == 0)
744    compare = APInt::tcCompare(significandParts(), rhs.significandParts(),
745                               partCount());
746
747  if(compare > 0)
748    return cmpGreaterThan;
749  else if(compare < 0)
750    return cmpLessThan;
751  else
752    return cmpEqual;
753}
754
755/* Handle overflow.  Sign is preserved.  We either become infinity or
756   the largest finite number.  */
757APFloat::opStatus
758APFloat::handleOverflow(roundingMode rounding_mode)
759{
760  /* Infinity?  */
761  if(rounding_mode == rmNearestTiesToEven
762     || rounding_mode == rmNearestTiesToAway
763     || (rounding_mode == rmTowardPositive && !sign)
764     || (rounding_mode == rmTowardNegative && sign))
765    {
766      category = fcInfinity;
767      return (opStatus) (opOverflow | opInexact);
768    }
769
770  /* Otherwise we become the largest finite number.  */
771  category = fcNormal;
772  exponent = semantics->maxExponent;
773  APInt::tcSetLeastSignificantBits(significandParts(), partCount(),
774                                   semantics->precision);
775
776  return opInexact;
777}
778
779/* Returns TRUE if, when truncating the current number, with BIT the
780   new LSB, with the given lost fraction and rounding mode, the result
781   would need to be rounded away from zero (i.e., by increasing the
782   signficand).  This routine must work for fcZero of both signs, and
783   fcNormal numbers.  */
784bool
785APFloat::roundAwayFromZero(roundingMode rounding_mode,
786                           lostFraction lost_fraction,
787                           unsigned int bit) const
788{
789  /* NaNs and infinities should not have lost fractions.  */
790  assert(category == fcNormal || category == fcZero);
791
792  /* Current callers never pass this so we don't handle it.  */
793  assert(lost_fraction != lfExactlyZero);
794
795  switch(rounding_mode) {
796  default:
797    assert(0);
798
799  case rmNearestTiesToAway:
800    return lost_fraction == lfExactlyHalf || lost_fraction == lfMoreThanHalf;
801
802  case rmNearestTiesToEven:
803    if(lost_fraction == lfMoreThanHalf)
804      return true;
805
806    /* Our zeroes don't have a significand to test.  */
807    if(lost_fraction == lfExactlyHalf && category != fcZero)
808      return APInt::tcExtractBit(significandParts(), bit);
809
810    return false;
811
812  case rmTowardZero:
813    return false;
814
815  case rmTowardPositive:
816    return sign == false;
817
818  case rmTowardNegative:
819    return sign == true;
820  }
821}
822
823APFloat::opStatus
824APFloat::normalize(roundingMode rounding_mode,
825                   lostFraction lost_fraction)
826{
827  unsigned int omsb;                /* One, not zero, based MSB.  */
828  int exponentChange;
829
830  if(category != fcNormal)
831    return opOK;
832
833  /* Before rounding normalize the exponent of fcNormal numbers.  */
834  omsb = significandMSB() + 1;
835
836  if(omsb) {
837    /* OMSB is numbered from 1.  We want to place it in the integer
838       bit numbered PRECISON if possible, with a compensating change in
839       the exponent.  */
840    exponentChange = omsb - semantics->precision;
841
842    /* If the resulting exponent is too high, overflow according to
843       the rounding mode.  */
844    if(exponent + exponentChange > semantics->maxExponent)
845      return handleOverflow(rounding_mode);
846
847    /* Subnormal numbers have exponent minExponent, and their MSB
848       is forced based on that.  */
849    if(exponent + exponentChange < semantics->minExponent)
850      exponentChange = semantics->minExponent - exponent;
851
852    /* Shifting left is easy as we don't lose precision.  */
853    if(exponentChange < 0) {
854      assert(lost_fraction == lfExactlyZero);
855
856      shiftSignificandLeft(-exponentChange);
857
858      return opOK;
859    }
860
861    if(exponentChange > 0) {
862      lostFraction lf;
863
864      /* Shift right and capture any new lost fraction.  */
865      lf = shiftSignificandRight(exponentChange);
866
867      lost_fraction = combineLostFractions(lf, lost_fraction);
868
869      /* Keep OMSB up-to-date.  */
870      if(omsb > (unsigned) exponentChange)
871        omsb -= (unsigned) exponentChange;
872      else
873        omsb = 0;
874    }
875  }
876
877  /* Now round the number according to rounding_mode given the lost
878     fraction.  */
879
880  /* As specified in IEEE 754, since we do not trap we do not report
881     underflow for exact results.  */
882  if(lost_fraction == lfExactlyZero) {
883    /* Canonicalize zeroes.  */
884    if(omsb == 0)
885      category = fcZero;
886
887    return opOK;
888  }
889
890  /* Increment the significand if we're rounding away from zero.  */
891  if(roundAwayFromZero(rounding_mode, lost_fraction, 0)) {
892    if(omsb == 0)
893      exponent = semantics->minExponent;
894
895    incrementSignificand();
896    omsb = significandMSB() + 1;
897
898    /* Did the significand increment overflow?  */
899    if(omsb == (unsigned) semantics->precision + 1) {
900      /* Renormalize by incrementing the exponent and shifting our
901         significand right one.  However if we already have the
902         maximum exponent we overflow to infinity.  */
903      if(exponent == semantics->maxExponent) {
904        category = fcInfinity;
905
906        return (opStatus) (opOverflow | opInexact);
907      }
908
909      shiftSignificandRight(1);
910
911      return opInexact;
912    }
913  }
914
915  /* The normal case - we were and are not denormal, and any
916     significand increment above didn't overflow.  */
917  if(omsb == semantics->precision)
918    return opInexact;
919
920  /* We have a non-zero denormal.  */
921  assert(omsb < semantics->precision);
922  assert(exponent == semantics->minExponent);
923
924  /* Canonicalize zeroes.  */
925  if(omsb == 0)
926    category = fcZero;
927
928  /* The fcZero case is a denormal that underflowed to zero.  */
929  return (opStatus) (opUnderflow | opInexact);
930}
931
932APFloat::opStatus
933APFloat::addOrSubtractSpecials(const APFloat &rhs, bool subtract)
934{
935  switch(convolve(category, rhs.category)) {
936  default:
937    assert(0);
938
939  case convolve(fcNaN, fcZero):
940  case convolve(fcNaN, fcNormal):
941  case convolve(fcNaN, fcInfinity):
942  case convolve(fcNaN, fcNaN):
943  case convolve(fcNormal, fcZero):
944  case convolve(fcInfinity, fcNormal):
945  case convolve(fcInfinity, fcZero):
946    return opOK;
947
948  case convolve(fcZero, fcNaN):
949  case convolve(fcNormal, fcNaN):
950  case convolve(fcInfinity, fcNaN):
951    category = fcNaN;
952    copySignificand(rhs);
953    return opOK;
954
955  case convolve(fcNormal, fcInfinity):
956  case convolve(fcZero, fcInfinity):
957    category = fcInfinity;
958    sign = rhs.sign ^ subtract;
959    return opOK;
960
961  case convolve(fcZero, fcNormal):
962    assign(rhs);
963    sign = rhs.sign ^ subtract;
964    return opOK;
965
966  case convolve(fcZero, fcZero):
967    /* Sign depends on rounding mode; handled by caller.  */
968    return opOK;
969
970  case convolve(fcInfinity, fcInfinity):
971    /* Differently signed infinities can only be validly
972       subtracted.  */
973    if(sign ^ rhs.sign != subtract) {
974      category = fcNaN;
975      // Arbitrary but deterministic value for significand
976      APInt::tcSet(significandParts(), ~0U, partCount());
977      return opInvalidOp;
978    }
979
980    return opOK;
981
982  case convolve(fcNormal, fcNormal):
983    return opDivByZero;
984  }
985}
986
987/* Add or subtract two normal numbers.  */
988lostFraction
989APFloat::addOrSubtractSignificand(const APFloat &rhs, bool subtract)
990{
991  integerPart carry;
992  lostFraction lost_fraction;
993  int bits;
994
995  /* Determine if the operation on the absolute values is effectively
996     an addition or subtraction.  */
997  subtract ^= (sign ^ rhs.sign);
998
999  /* Are we bigger exponent-wise than the RHS?  */
1000  bits = exponent - rhs.exponent;
1001
1002  /* Subtraction is more subtle than one might naively expect.  */
1003  if(subtract) {
1004    APFloat temp_rhs(rhs);
1005    bool reverse;
1006
1007    if (bits == 0) {
1008      reverse = compareAbsoluteValue(temp_rhs) == cmpLessThan;
1009      lost_fraction = lfExactlyZero;
1010    } else if (bits > 0) {
1011      lost_fraction = temp_rhs.shiftSignificandRight(bits - 1);
1012      shiftSignificandLeft(1);
1013      reverse = false;
1014    } else {
1015      lost_fraction = shiftSignificandRight(-bits - 1);
1016      temp_rhs.shiftSignificandLeft(1);
1017      reverse = true;
1018    }
1019
1020    if (reverse) {
1021      carry = temp_rhs.subtractSignificand
1022        (*this, lost_fraction != lfExactlyZero);
1023      copySignificand(temp_rhs);
1024      sign = !sign;
1025    } else {
1026      carry = subtractSignificand
1027        (temp_rhs, lost_fraction != lfExactlyZero);
1028    }
1029
1030    /* Invert the lost fraction - it was on the RHS and
1031       subtracted.  */
1032    if(lost_fraction == lfLessThanHalf)
1033      lost_fraction = lfMoreThanHalf;
1034    else if(lost_fraction == lfMoreThanHalf)
1035      lost_fraction = lfLessThanHalf;
1036
1037    /* The code above is intended to ensure that no borrow is
1038       necessary.  */
1039    assert(!carry);
1040  } else {
1041    if(bits > 0) {
1042      APFloat temp_rhs(rhs);
1043
1044      lost_fraction = temp_rhs.shiftSignificandRight(bits);
1045      carry = addSignificand(temp_rhs);
1046    } else {
1047      lost_fraction = shiftSignificandRight(-bits);
1048      carry = addSignificand(rhs);
1049    }
1050
1051    /* We have a guard bit; generating a carry cannot happen.  */
1052    assert(!carry);
1053  }
1054
1055  return lost_fraction;
1056}
1057
1058APFloat::opStatus
1059APFloat::multiplySpecials(const APFloat &rhs)
1060{
1061  switch(convolve(category, rhs.category)) {
1062  default:
1063    assert(0);
1064
1065  case convolve(fcNaN, fcZero):
1066  case convolve(fcNaN, fcNormal):
1067  case convolve(fcNaN, fcInfinity):
1068  case convolve(fcNaN, fcNaN):
1069    return opOK;
1070
1071  case convolve(fcZero, fcNaN):
1072  case convolve(fcNormal, fcNaN):
1073  case convolve(fcInfinity, fcNaN):
1074    category = fcNaN;
1075    copySignificand(rhs);
1076    return opOK;
1077
1078  case convolve(fcNormal, fcInfinity):
1079  case convolve(fcInfinity, fcNormal):
1080  case convolve(fcInfinity, fcInfinity):
1081    category = fcInfinity;
1082    return opOK;
1083
1084  case convolve(fcZero, fcNormal):
1085  case convolve(fcNormal, fcZero):
1086  case convolve(fcZero, fcZero):
1087    category = fcZero;
1088    return opOK;
1089
1090  case convolve(fcZero, fcInfinity):
1091  case convolve(fcInfinity, fcZero):
1092    category = fcNaN;
1093    // Arbitrary but deterministic value for significand
1094    APInt::tcSet(significandParts(), ~0U, partCount());
1095    return opInvalidOp;
1096
1097  case convolve(fcNormal, fcNormal):
1098    return opOK;
1099  }
1100}
1101
1102APFloat::opStatus
1103APFloat::divideSpecials(const APFloat &rhs)
1104{
1105  switch(convolve(category, rhs.category)) {
1106  default:
1107    assert(0);
1108
1109  case convolve(fcNaN, fcZero):
1110  case convolve(fcNaN, fcNormal):
1111  case convolve(fcNaN, fcInfinity):
1112  case convolve(fcNaN, fcNaN):
1113  case convolve(fcInfinity, fcZero):
1114  case convolve(fcInfinity, fcNormal):
1115  case convolve(fcZero, fcInfinity):
1116  case convolve(fcZero, fcNormal):
1117    return opOK;
1118
1119  case convolve(fcZero, fcNaN):
1120  case convolve(fcNormal, fcNaN):
1121  case convolve(fcInfinity, fcNaN):
1122    category = fcNaN;
1123    copySignificand(rhs);
1124    return opOK;
1125
1126  case convolve(fcNormal, fcInfinity):
1127    category = fcZero;
1128    return opOK;
1129
1130  case convolve(fcNormal, fcZero):
1131    category = fcInfinity;
1132    return opDivByZero;
1133
1134  case convolve(fcInfinity, fcInfinity):
1135  case convolve(fcZero, fcZero):
1136    category = fcNaN;
1137    // Arbitrary but deterministic value for significand
1138    APInt::tcSet(significandParts(), ~0U, partCount());
1139    return opInvalidOp;
1140
1141  case convolve(fcNormal, fcNormal):
1142    return opOK;
1143  }
1144}
1145
1146/* Change sign.  */
1147void
1148APFloat::changeSign()
1149{
1150  /* Look mummy, this one's easy.  */
1151  sign = !sign;
1152}
1153
1154void
1155APFloat::clearSign()
1156{
1157  /* So is this one. */
1158  sign = 0;
1159}
1160
1161void
1162APFloat::copySign(const APFloat &rhs)
1163{
1164  /* And this one. */
1165  sign = rhs.sign;
1166}
1167
1168/* Normalized addition or subtraction.  */
1169APFloat::opStatus
1170APFloat::addOrSubtract(const APFloat &rhs, roundingMode rounding_mode,
1171                       bool subtract)
1172{
1173  opStatus fs;
1174
1175  fs = addOrSubtractSpecials(rhs, subtract);
1176
1177  /* This return code means it was not a simple case.  */
1178  if(fs == opDivByZero) {
1179    lostFraction lost_fraction;
1180
1181    lost_fraction = addOrSubtractSignificand(rhs, subtract);
1182    fs = normalize(rounding_mode, lost_fraction);
1183
1184    /* Can only be zero if we lost no fraction.  */
1185    assert(category != fcZero || lost_fraction == lfExactlyZero);
1186  }
1187
1188  /* If two numbers add (exactly) to zero, IEEE 754 decrees it is a
1189     positive zero unless rounding to minus infinity, except that
1190     adding two like-signed zeroes gives that zero.  */
1191  if(category == fcZero) {
1192    if(rhs.category != fcZero || (sign == rhs.sign) == subtract)
1193      sign = (rounding_mode == rmTowardNegative);
1194  }
1195
1196  return fs;
1197}
1198
1199/* Normalized addition.  */
1200APFloat::opStatus
1201APFloat::add(const APFloat &rhs, roundingMode rounding_mode)
1202{
1203  assert(semantics != (const llvm::fltSemantics* const)&PPCDoubleDouble &&
1204         "Compile-time arithmetic on PPC long double not supported yet");
1205  return addOrSubtract(rhs, rounding_mode, false);
1206}
1207
1208/* Normalized subtraction.  */
1209APFloat::opStatus
1210APFloat::subtract(const APFloat &rhs, roundingMode rounding_mode)
1211{
1212  assert(semantics != (const llvm::fltSemantics* const)&PPCDoubleDouble &&
1213         "Compile-time arithmetic on PPC long double not supported yet");
1214  return addOrSubtract(rhs, rounding_mode, true);
1215}
1216
1217/* Normalized multiply.  */
1218APFloat::opStatus
1219APFloat::multiply(const APFloat &rhs, roundingMode rounding_mode)
1220{
1221  assert(semantics != (const llvm::fltSemantics* const)&PPCDoubleDouble &&
1222         "Compile-time arithmetic on PPC long double not supported yet");
1223  opStatus fs;
1224
1225  sign ^= rhs.sign;
1226  fs = multiplySpecials(rhs);
1227
1228  if(category == fcNormal) {
1229    lostFraction lost_fraction = multiplySignificand(rhs, 0);
1230    fs = normalize(rounding_mode, lost_fraction);
1231    if(lost_fraction != lfExactlyZero)
1232      fs = (opStatus) (fs | opInexact);
1233  }
1234
1235  return fs;
1236}
1237
1238/* Normalized divide.  */
1239APFloat::opStatus
1240APFloat::divide(const APFloat &rhs, roundingMode rounding_mode)
1241{
1242  assert(semantics != (const llvm::fltSemantics* const)&PPCDoubleDouble &&
1243         "Compile-time arithmetic on PPC long double not supported yet");
1244  opStatus fs;
1245
1246  sign ^= rhs.sign;
1247  fs = divideSpecials(rhs);
1248
1249  if(category == fcNormal) {
1250    lostFraction lost_fraction = divideSignificand(rhs);
1251    fs = normalize(rounding_mode, lost_fraction);
1252    if(lost_fraction != lfExactlyZero)
1253      fs = (opStatus) (fs | opInexact);
1254  }
1255
1256  return fs;
1257}
1258
1259/* Normalized remainder.  This is not currently doing TRT.  */
1260APFloat::opStatus
1261APFloat::mod(const APFloat &rhs, roundingMode rounding_mode)
1262{
1263  assert(semantics != (const llvm::fltSemantics* const)&PPCDoubleDouble &&
1264         "Compile-time arithmetic on PPC long double not supported yet");
1265  opStatus fs;
1266  APFloat V = *this;
1267  unsigned int origSign = sign;
1268  fs = V.divide(rhs, rmNearestTiesToEven);
1269  if (fs == opDivByZero)
1270    return fs;
1271
1272  int parts = partCount();
1273  integerPart *x = new integerPart[parts];
1274  fs = V.convertToInteger(x, parts * integerPartWidth, true,
1275                          rmNearestTiesToEven);
1276  if (fs==opInvalidOp)
1277    return fs;
1278
1279  fs = V.convertFromZeroExtendedInteger(x, parts * integerPartWidth, true,
1280                                        rmNearestTiesToEven);
1281  assert(fs==opOK);   // should always work
1282
1283  fs = V.multiply(rhs, rounding_mode);
1284  assert(fs==opOK || fs==opInexact);   // should not overflow or underflow
1285
1286  fs = subtract(V, rounding_mode);
1287  assert(fs==opOK || fs==opInexact);   // likewise
1288
1289  if (isZero())
1290    sign = origSign;    // IEEE754 requires this
1291  delete[] x;
1292  return fs;
1293}
1294
1295/* Normalized fused-multiply-add.  */
1296APFloat::opStatus
1297APFloat::fusedMultiplyAdd(const APFloat &multiplicand,
1298                          const APFloat &addend,
1299                          roundingMode rounding_mode)
1300{
1301  assert(semantics != (const llvm::fltSemantics* const)&PPCDoubleDouble &&
1302         "Compile-time arithmetic on PPC long double not supported yet");
1303  opStatus fs;
1304
1305  /* Post-multiplication sign, before addition.  */
1306  sign ^= multiplicand.sign;
1307
1308  /* If and only if all arguments are normal do we need to do an
1309     extended-precision calculation.  */
1310  if(category == fcNormal
1311     && multiplicand.category == fcNormal
1312     && addend.category == fcNormal) {
1313    lostFraction lost_fraction;
1314
1315    lost_fraction = multiplySignificand(multiplicand, &addend);
1316    fs = normalize(rounding_mode, lost_fraction);
1317    if(lost_fraction != lfExactlyZero)
1318      fs = (opStatus) (fs | opInexact);
1319
1320    /* If two numbers add (exactly) to zero, IEEE 754 decrees it is a
1321       positive zero unless rounding to minus infinity, except that
1322       adding two like-signed zeroes gives that zero.  */
1323    if(category == fcZero && sign != addend.sign)
1324      sign = (rounding_mode == rmTowardNegative);
1325  } else {
1326    fs = multiplySpecials(multiplicand);
1327
1328    /* FS can only be opOK or opInvalidOp.  There is no more work
1329       to do in the latter case.  The IEEE-754R standard says it is
1330       implementation-defined in this case whether, if ADDEND is a
1331       quiet NaN, we raise invalid op; this implementation does so.
1332
1333       If we need to do the addition we can do so with normal
1334       precision.  */
1335    if(fs == opOK)
1336      fs = addOrSubtract(addend, rounding_mode, false);
1337  }
1338
1339  return fs;
1340}
1341
1342/* Comparison requires normalized numbers.  */
1343APFloat::cmpResult
1344APFloat::compare(const APFloat &rhs) const
1345{
1346  assert(semantics != (const llvm::fltSemantics* const)&PPCDoubleDouble &&
1347         "Compile-time arithmetic on PPC long double not supported yet");
1348  cmpResult result;
1349
1350  assert(semantics == rhs.semantics);
1351
1352  switch(convolve(category, rhs.category)) {
1353  default:
1354    assert(0);
1355
1356  case convolve(fcNaN, fcZero):
1357  case convolve(fcNaN, fcNormal):
1358  case convolve(fcNaN, fcInfinity):
1359  case convolve(fcNaN, fcNaN):
1360  case convolve(fcZero, fcNaN):
1361  case convolve(fcNormal, fcNaN):
1362  case convolve(fcInfinity, fcNaN):
1363    return cmpUnordered;
1364
1365  case convolve(fcInfinity, fcNormal):
1366  case convolve(fcInfinity, fcZero):
1367  case convolve(fcNormal, fcZero):
1368    if(sign)
1369      return cmpLessThan;
1370    else
1371      return cmpGreaterThan;
1372
1373  case convolve(fcNormal, fcInfinity):
1374  case convolve(fcZero, fcInfinity):
1375  case convolve(fcZero, fcNormal):
1376    if(rhs.sign)
1377      return cmpGreaterThan;
1378    else
1379      return cmpLessThan;
1380
1381  case convolve(fcInfinity, fcInfinity):
1382    if(sign == rhs.sign)
1383      return cmpEqual;
1384    else if(sign)
1385      return cmpLessThan;
1386    else
1387      return cmpGreaterThan;
1388
1389  case convolve(fcZero, fcZero):
1390    return cmpEqual;
1391
1392  case convolve(fcNormal, fcNormal):
1393    break;
1394  }
1395
1396  /* Two normal numbers.  Do they have the same sign?  */
1397  if(sign != rhs.sign) {
1398    if(sign)
1399      result = cmpLessThan;
1400    else
1401      result = cmpGreaterThan;
1402  } else {
1403    /* Compare absolute values; invert result if negative.  */
1404    result = compareAbsoluteValue(rhs);
1405
1406    if(sign) {
1407      if(result == cmpLessThan)
1408        result = cmpGreaterThan;
1409      else if(result == cmpGreaterThan)
1410        result = cmpLessThan;
1411    }
1412  }
1413
1414  return result;
1415}
1416
1417APFloat::opStatus
1418APFloat::convert(const fltSemantics &toSemantics,
1419                 roundingMode rounding_mode)
1420{
1421  assert(semantics != (const llvm::fltSemantics* const)&PPCDoubleDouble &&
1422         "Compile-time arithmetic on PPC long double not supported yet");
1423  lostFraction lostFraction;
1424  unsigned int newPartCount, oldPartCount;
1425  opStatus fs;
1426
1427  lostFraction = lfExactlyZero;
1428  newPartCount = partCountForBits(toSemantics.precision + 1);
1429  oldPartCount = partCount();
1430
1431  /* Handle storage complications.  If our new form is wider,
1432     re-allocate our bit pattern into wider storage.  If it is
1433     narrower, we ignore the excess parts, but if narrowing to a
1434     single part we need to free the old storage.
1435     Be careful not to reference significandParts for zeroes
1436     and infinities, since it aborts.  */
1437  if (newPartCount > oldPartCount) {
1438    integerPart *newParts;
1439    newParts = new integerPart[newPartCount];
1440    APInt::tcSet(newParts, 0, newPartCount);
1441    if (category==fcNormal || category==fcNaN)
1442      APInt::tcAssign(newParts, significandParts(), oldPartCount);
1443    freeSignificand();
1444    significand.parts = newParts;
1445  } else if (newPartCount < oldPartCount) {
1446    /* Capture any lost fraction through truncation of parts so we get
1447       correct rounding whilst normalizing.  */
1448    if (category==fcNormal)
1449      lostFraction = lostFractionThroughTruncation
1450        (significandParts(), oldPartCount, toSemantics.precision);
1451    if (newPartCount == 1) {
1452        integerPart newPart = 0;
1453        if (category==fcNormal || category==fcNaN)
1454          newPart = significandParts()[0];
1455        freeSignificand();
1456        significand.part = newPart;
1457    }
1458  }
1459
1460  if(category == fcNormal) {
1461    /* Re-interpret our bit-pattern.  */
1462    exponent += toSemantics.precision - semantics->precision;
1463    semantics = &toSemantics;
1464    fs = normalize(rounding_mode, lostFraction);
1465  } else if (category == fcNaN) {
1466    int shift = toSemantics.precision - semantics->precision;
1467    // No normalization here, just truncate
1468    if (shift>0)
1469      APInt::tcShiftLeft(significandParts(), newPartCount, shift);
1470    else if (shift < 0)
1471      APInt::tcShiftRight(significandParts(), newPartCount, -shift);
1472    // gcc forces the Quiet bit on, which means (float)(double)(float_sNan)
1473    // does not give you back the same bits.  This is dubious, and we
1474    // don't currently do it.  You're really supposed to get
1475    // an invalid operation signal at runtime, but nobody does that.
1476    semantics = &toSemantics;
1477    fs = opOK;
1478  } else {
1479    semantics = &toSemantics;
1480    fs = opOK;
1481  }
1482
1483  return fs;
1484}
1485
1486/* Convert a floating point number to an integer according to the
1487   rounding mode.  If the rounded integer value is out of range this
1488   returns an invalid operation exception.  If the rounded value is in
1489   range but the floating point number is not the exact integer, the C
1490   standard doesn't require an inexact exception to be raised.  IEEE
1491   854 does require it so we do that.
1492
1493   Note that for conversions to integer type the C standard requires
1494   round-to-zero to always be used.  */
1495APFloat::opStatus
1496APFloat::convertToInteger(integerPart *parts, unsigned int width,
1497                          bool isSigned,
1498                          roundingMode rounding_mode) const
1499{
1500  assert(semantics != (const llvm::fltSemantics* const)&PPCDoubleDouble &&
1501         "Compile-time arithmetic on PPC long double not supported yet");
1502  lostFraction lost_fraction;
1503  unsigned int msb, partsCount;
1504  int bits;
1505
1506  partsCount = partCountForBits(width);
1507
1508  /* Handle the three special cases first.  We produce
1509     a deterministic result even for the Invalid cases. */
1510  if (category == fcNaN) {
1511    // Neither sign nor isSigned affects this.
1512    APInt::tcSet(parts, 0, partsCount);
1513    return opInvalidOp;
1514  }
1515  if (category == fcInfinity) {
1516    if (!sign && isSigned)
1517      APInt::tcSetLeastSignificantBits(parts, partsCount, width-1);
1518    else if (!sign && !isSigned)
1519      APInt::tcSetLeastSignificantBits(parts, partsCount, width);
1520    else if (sign && isSigned) {
1521      APInt::tcSetLeastSignificantBits(parts, partsCount, 1);
1522      APInt::tcShiftLeft(parts, partsCount, width-1);
1523    } else // sign && !isSigned
1524      APInt::tcSet(parts, 0, partsCount);
1525    return opInvalidOp;
1526  }
1527  if (category == fcZero) {
1528    APInt::tcSet(parts, 0, partsCount);
1529    return opOK;
1530  }
1531
1532  /* Shift the bit pattern so the fraction is lost.  */
1533  APFloat tmp(*this);
1534
1535  bits = (int) semantics->precision - 1 - exponent;
1536
1537  if(bits > 0) {
1538    lost_fraction = tmp.shiftSignificandRight(bits);
1539  } else {
1540    if (-bits >= semantics->precision) {
1541      // Unrepresentably large.
1542      if (!sign && isSigned)
1543        APInt::tcSetLeastSignificantBits(parts, partsCount, width-1);
1544      else if (!sign && !isSigned)
1545        APInt::tcSetLeastSignificantBits(parts, partsCount, width);
1546      else if (sign && isSigned) {
1547        APInt::tcSetLeastSignificantBits(parts, partsCount, 1);
1548        APInt::tcShiftLeft(parts, partsCount, width-1);
1549      } else // sign && !isSigned
1550        APInt::tcSet(parts, 0, partsCount);
1551      return (opStatus)(opOverflow | opInexact);
1552    }
1553    tmp.shiftSignificandLeft(-bits);
1554    lost_fraction = lfExactlyZero;
1555  }
1556
1557  if(lost_fraction != lfExactlyZero
1558     && tmp.roundAwayFromZero(rounding_mode, lost_fraction, 0))
1559    tmp.incrementSignificand();
1560
1561  msb = tmp.significandMSB();
1562
1563  /* Negative numbers cannot be represented as unsigned.  */
1564  if(!isSigned && tmp.sign && msb != -1U)
1565    return opInvalidOp;
1566
1567  /* It takes exponent + 1 bits to represent the truncated floating
1568     point number without its sign.  We lose a bit for the sign, but
1569     the maximally negative integer is a special case.  */
1570  if(msb + 1 > width)                /* !! Not same as msb >= width !! */
1571    return opInvalidOp;
1572
1573  if(isSigned && msb + 1 == width
1574     && (!tmp.sign || tmp.significandLSB() != msb))
1575    return opInvalidOp;
1576
1577  APInt::tcAssign(parts, tmp.significandParts(), partsCount);
1578
1579  if(tmp.sign)
1580    APInt::tcNegate(parts, partsCount);
1581
1582  if(lost_fraction == lfExactlyZero)
1583    return opOK;
1584  else
1585    return opInexact;
1586}
1587
1588/* Convert an unsigned integer SRC to a floating point number,
1589   rounding according to ROUNDING_MODE.  The sign of the floating
1590   point number is not modified.  */
1591APFloat::opStatus
1592APFloat::convertFromUnsignedParts(const integerPart *src,
1593                                  unsigned int srcCount,
1594                                  roundingMode rounding_mode)
1595{
1596  unsigned int omsb, precision, dstCount;
1597  integerPart *dst;
1598  lostFraction lost_fraction;
1599
1600  category = fcNormal;
1601  omsb = APInt::tcMSB(src, srcCount) + 1;
1602  dst = significandParts();
1603  dstCount = partCount();
1604  precision = semantics->precision;
1605
1606  /* We want the most significant PRECISON bits of SRC.  There may not
1607     be that many; extract what we can.  */
1608  if (precision <= omsb) {
1609    exponent = omsb - 1;
1610    lost_fraction = lostFractionThroughTruncation(src, srcCount,
1611                                                  omsb - precision);
1612    APInt::tcExtract(dst, dstCount, src, precision, omsb - precision);
1613  } else {
1614    exponent = precision - 1;
1615    lost_fraction = lfExactlyZero;
1616    APInt::tcExtract(dst, dstCount, src, omsb, 0);
1617  }
1618
1619  return normalize(rounding_mode, lost_fraction);
1620}
1621
1622/* Convert a two's complement integer SRC to a floating point number,
1623   rounding according to ROUNDING_MODE.  ISSIGNED is true if the
1624   integer is signed, in which case it must be sign-extended.  */
1625APFloat::opStatus
1626APFloat::convertFromSignExtendedInteger(const integerPart *src,
1627                                        unsigned int srcCount,
1628                                        bool isSigned,
1629                                        roundingMode rounding_mode)
1630{
1631  assert(semantics != (const llvm::fltSemantics* const)&PPCDoubleDouble &&
1632         "Compile-time arithmetic on PPC long double not supported yet");
1633  opStatus status;
1634
1635  if (isSigned
1636      && APInt::tcExtractBit(src, srcCount * integerPartWidth - 1)) {
1637    integerPart *copy;
1638
1639    /* If we're signed and negative negate a copy.  */
1640    sign = true;
1641    copy = new integerPart[srcCount];
1642    APInt::tcAssign(copy, src, srcCount);
1643    APInt::tcNegate(copy, srcCount);
1644    status = convertFromUnsignedParts(copy, srcCount, rounding_mode);
1645    delete [] copy;
1646  } else {
1647    sign = false;
1648    status = convertFromUnsignedParts(src, srcCount, rounding_mode);
1649  }
1650
1651  return status;
1652}
1653
1654/* FIXME: should this just take a const APInt reference?  */
1655APFloat::opStatus
1656APFloat::convertFromZeroExtendedInteger(const integerPart *parts,
1657                                        unsigned int width, bool isSigned,
1658                                        roundingMode rounding_mode)
1659{
1660  assert(semantics != (const llvm::fltSemantics* const)&PPCDoubleDouble &&
1661         "Compile-time arithmetic on PPC long double not supported yet");
1662  unsigned int partCount = partCountForBits(width);
1663  APInt api = APInt(width, partCount, parts);
1664
1665  sign = false;
1666  if(isSigned && APInt::tcExtractBit(parts, width - 1)) {
1667    sign = true;
1668    api = -api;
1669  }
1670
1671  return convertFromUnsignedParts(api.getRawData(), partCount, rounding_mode);
1672}
1673
1674APFloat::opStatus
1675APFloat::convertFromHexadecimalString(const char *p,
1676                                      roundingMode rounding_mode)
1677{
1678  assert(semantics != (const llvm::fltSemantics* const)&PPCDoubleDouble &&
1679         "Compile-time arithmetic on PPC long double not supported yet");
1680  lostFraction lost_fraction;
1681  integerPart *significand;
1682  unsigned int bitPos, partsCount;
1683  const char *dot, *firstSignificantDigit;
1684
1685  zeroSignificand();
1686  exponent = 0;
1687  category = fcNormal;
1688
1689  significand = significandParts();
1690  partsCount = partCount();
1691  bitPos = partsCount * integerPartWidth;
1692
1693  /* Skip leading zeroes and any (hexa)decimal point.  */
1694  p = skipLeadingZeroesAndAnyDot(p, &dot);
1695  firstSignificantDigit = p;
1696
1697  for(;;) {
1698    integerPart hex_value;
1699
1700    if(*p == '.') {
1701      assert(dot == 0);
1702      dot = p++;
1703    }
1704
1705    hex_value = hexDigitValue(*p);
1706    if(hex_value == -1U) {
1707      lost_fraction = lfExactlyZero;
1708      break;
1709    }
1710
1711    p++;
1712
1713    /* Store the number whilst 4-bit nibbles remain.  */
1714    if(bitPos) {
1715      bitPos -= 4;
1716      hex_value <<= bitPos % integerPartWidth;
1717      significand[bitPos / integerPartWidth] |= hex_value;
1718    } else {
1719      lost_fraction = trailingHexadecimalFraction(p, hex_value);
1720      while(hexDigitValue(*p) != -1U)
1721        p++;
1722      break;
1723    }
1724  }
1725
1726  /* Hex floats require an exponent but not a hexadecimal point.  */
1727  assert(*p == 'p' || *p == 'P');
1728
1729  /* Ignore the exponent if we are zero.  */
1730  if(p != firstSignificantDigit) {
1731    int expAdjustment;
1732
1733    /* Implicit hexadecimal point?  */
1734    if(!dot)
1735      dot = p;
1736
1737    /* Calculate the exponent adjustment implicit in the number of
1738       significant digits.  */
1739    expAdjustment = dot - firstSignificantDigit;
1740    if(expAdjustment < 0)
1741      expAdjustment++;
1742    expAdjustment = expAdjustment * 4 - 1;
1743
1744    /* Adjust for writing the significand starting at the most
1745       significant nibble.  */
1746    expAdjustment += semantics->precision;
1747    expAdjustment -= partsCount * integerPartWidth;
1748
1749    /* Adjust for the given exponent.  */
1750    exponent = totalExponent(p, expAdjustment);
1751  }
1752
1753  return normalize(rounding_mode, lost_fraction);
1754}
1755
1756APFloat::opStatus
1757APFloat::convertFromString(const char *p, roundingMode rounding_mode)
1758{
1759  assert(semantics != (const llvm::fltSemantics* const)&PPCDoubleDouble &&
1760         "Compile-time arithmetic on PPC long double not supported yet");
1761  /* Handle a leading minus sign.  */
1762  if(*p == '-')
1763    sign = 1, p++;
1764  else
1765    sign = 0;
1766
1767  if(p[0] == '0' && (p[1] == 'x' || p[1] == 'X'))
1768    return convertFromHexadecimalString(p + 2, rounding_mode);
1769
1770  assert(0 && "Decimal to binary conversions not yet implemented");
1771  abort();
1772}
1773
1774/* Write out a hexadecimal representation of the floating point value
1775   to DST, which must be of sufficient size, in the C99 form
1776   [-]0xh.hhhhp[+-]d.  Return the number of characters written,
1777   excluding the terminating NUL.
1778
1779   If UPPERCASE, the output is in upper case, otherwise in lower case.
1780
1781   HEXDIGITS digits appear altogether, rounding the value if
1782   necessary.  If HEXDIGITS is 0, the minimal precision to display the
1783   number precisely is used instead.  If nothing would appear after
1784   the decimal point it is suppressed.
1785
1786   The decimal exponent is always printed and has at least one digit.
1787   Zero values display an exponent of zero.  Infinities and NaNs
1788   appear as "infinity" or "nan" respectively.
1789
1790   The above rules are as specified by C99.  There is ambiguity about
1791   what the leading hexadecimal digit should be.  This implementation
1792   uses whatever is necessary so that the exponent is displayed as
1793   stored.  This implies the exponent will fall within the IEEE format
1794   range, and the leading hexadecimal digit will be 0 (for denormals),
1795   1 (normal numbers) or 2 (normal numbers rounded-away-from-zero with
1796   any other digits zero).
1797*/
1798unsigned int
1799APFloat::convertToHexString(char *dst, unsigned int hexDigits,
1800                            bool upperCase, roundingMode rounding_mode) const
1801{
1802  assert(semantics != (const llvm::fltSemantics* const)&PPCDoubleDouble &&
1803         "Compile-time arithmetic on PPC long double not supported yet");
1804  char *p;
1805
1806  p = dst;
1807  if (sign)
1808    *dst++ = '-';
1809
1810  switch (category) {
1811  case fcInfinity:
1812    memcpy (dst, upperCase ? infinityU: infinityL, sizeof infinityU - 1);
1813    dst += sizeof infinityL - 1;
1814    break;
1815
1816  case fcNaN:
1817    memcpy (dst, upperCase ? NaNU: NaNL, sizeof NaNU - 1);
1818    dst += sizeof NaNU - 1;
1819    break;
1820
1821  case fcZero:
1822    *dst++ = '0';
1823    *dst++ = upperCase ? 'X': 'x';
1824    *dst++ = '0';
1825    if (hexDigits > 1) {
1826      *dst++ = '.';
1827      memset (dst, '0', hexDigits - 1);
1828      dst += hexDigits - 1;
1829    }
1830    *dst++ = upperCase ? 'P': 'p';
1831    *dst++ = '0';
1832    break;
1833
1834  case fcNormal:
1835    dst = convertNormalToHexString (dst, hexDigits, upperCase, rounding_mode);
1836    break;
1837  }
1838
1839  *dst = 0;
1840
1841  return dst - p;
1842}
1843
1844/* Does the hard work of outputting the correctly rounded hexadecimal
1845   form of a normal floating point number with the specified number of
1846   hexadecimal digits.  If HEXDIGITS is zero the minimum number of
1847   digits necessary to print the value precisely is output.  */
1848char *
1849APFloat::convertNormalToHexString(char *dst, unsigned int hexDigits,
1850                                  bool upperCase,
1851                                  roundingMode rounding_mode) const
1852{
1853  unsigned int count, valueBits, shift, partsCount, outputDigits;
1854  const char *hexDigitChars;
1855  const integerPart *significand;
1856  char *p;
1857  bool roundUp;
1858
1859  *dst++ = '0';
1860  *dst++ = upperCase ? 'X': 'x';
1861
1862  roundUp = false;
1863  hexDigitChars = upperCase ? hexDigitsUpper: hexDigitsLower;
1864
1865  significand = significandParts();
1866  partsCount = partCount();
1867
1868  /* +3 because the first digit only uses the single integer bit, so
1869     we have 3 virtual zero most-significant-bits.  */
1870  valueBits = semantics->precision + 3;
1871  shift = integerPartWidth - valueBits % integerPartWidth;
1872
1873  /* The natural number of digits required ignoring trailing
1874     insignificant zeroes.  */
1875  outputDigits = (valueBits - significandLSB () + 3) / 4;
1876
1877  /* hexDigits of zero means use the required number for the
1878     precision.  Otherwise, see if we are truncating.  If we are,
1879     find out if we need to round away from zero.  */
1880  if (hexDigits) {
1881    if (hexDigits < outputDigits) {
1882      /* We are dropping non-zero bits, so need to check how to round.
1883         "bits" is the number of dropped bits.  */
1884      unsigned int bits;
1885      lostFraction fraction;
1886
1887      bits = valueBits - hexDigits * 4;
1888      fraction = lostFractionThroughTruncation (significand, partsCount, bits);
1889      roundUp = roundAwayFromZero(rounding_mode, fraction, bits);
1890    }
1891    outputDigits = hexDigits;
1892  }
1893
1894  /* Write the digits consecutively, and start writing in the location
1895     of the hexadecimal point.  We move the most significant digit
1896     left and add the hexadecimal point later.  */
1897  p = ++dst;
1898
1899  count = (valueBits + integerPartWidth - 1) / integerPartWidth;
1900
1901  while (outputDigits && count) {
1902    integerPart part;
1903
1904    /* Put the most significant integerPartWidth bits in "part".  */
1905    if (--count == partsCount)
1906      part = 0;  /* An imaginary higher zero part.  */
1907    else
1908      part = significand[count] << shift;
1909
1910    if (count && shift)
1911      part |= significand[count - 1] >> (integerPartWidth - shift);
1912
1913    /* Convert as much of "part" to hexdigits as we can.  */
1914    unsigned int curDigits = integerPartWidth / 4;
1915
1916    if (curDigits > outputDigits)
1917      curDigits = outputDigits;
1918    dst += partAsHex (dst, part, curDigits, hexDigitChars);
1919    outputDigits -= curDigits;
1920  }
1921
1922  if (roundUp) {
1923    char *q = dst;
1924
1925    /* Note that hexDigitChars has a trailing '0'.  */
1926    do {
1927      q--;
1928      *q = hexDigitChars[hexDigitValue (*q) + 1];
1929    } while (*q == '0');
1930    assert (q >= p);
1931  } else {
1932    /* Add trailing zeroes.  */
1933    memset (dst, '0', outputDigits);
1934    dst += outputDigits;
1935  }
1936
1937  /* Move the most significant digit to before the point, and if there
1938     is something after the decimal point add it.  This must come
1939     after rounding above.  */
1940  p[-1] = p[0];
1941  if (dst -1 == p)
1942    dst--;
1943  else
1944    p[0] = '.';
1945
1946  /* Finally output the exponent.  */
1947  *dst++ = upperCase ? 'P': 'p';
1948
1949  return writeSignedDecimal (dst, exponent);
1950}
1951
1952// For good performance it is desirable for different APFloats
1953// to produce different integers.
1954uint32_t
1955APFloat::getHashValue() const
1956{
1957  if (category==fcZero) return sign<<8 | semantics->precision ;
1958  else if (category==fcInfinity) return sign<<9 | semantics->precision;
1959  else if (category==fcNaN) return 1<<10 | semantics->precision;
1960  else {
1961    uint32_t hash = sign<<11 | semantics->precision | exponent<<12;
1962    const integerPart* p = significandParts();
1963    for (int i=partCount(); i>0; i--, p++)
1964      hash ^= ((uint32_t)*p) ^ (*p)>>32;
1965    return hash;
1966  }
1967}
1968
1969// Conversion from APFloat to/from host float/double.  It may eventually be
1970// possible to eliminate these and have everybody deal with APFloats, but that
1971// will take a while.  This approach will not easily extend to long double.
1972// Current implementation requires integerPartWidth==64, which is correct at
1973// the moment but could be made more general.
1974
1975// Denormals have exponent minExponent in APFloat, but minExponent-1 in
1976// the actual IEEE respresentations.  We compensate for that here.
1977
1978APInt
1979APFloat::convertF80LongDoubleAPFloatToAPInt() const
1980{
1981  assert(semantics == (const llvm::fltSemantics* const)&x87DoubleExtended);
1982  assert (partCount()==2);
1983
1984  uint64_t myexponent, mysignificand;
1985
1986  if (category==fcNormal) {
1987    myexponent = exponent+16383; //bias
1988    mysignificand = significandParts()[0];
1989    if (myexponent==1 && !(mysignificand & 0x8000000000000000ULL))
1990      myexponent = 0;   // denormal
1991  } else if (category==fcZero) {
1992    myexponent = 0;
1993    mysignificand = 0;
1994  } else if (category==fcInfinity) {
1995    myexponent = 0x7fff;
1996    mysignificand = 0x8000000000000000ULL;
1997  } else {
1998    assert(category == fcNaN && "Unknown category");
1999    myexponent = 0x7fff;
2000    mysignificand = significandParts()[0];
2001  }
2002
2003  uint64_t words[2];
2004  words[0] =  (((uint64_t)sign & 1) << 63) |
2005              ((myexponent & 0x7fff) <<  48) |
2006              ((mysignificand >>16) & 0xffffffffffffLL);
2007  words[1] = mysignificand & 0xffff;
2008  return APInt(80, 2, words);
2009}
2010
2011APInt
2012APFloat::convertPPCDoubleDoubleAPFloatToAPInt() const
2013{
2014  assert(semantics == (const llvm::fltSemantics* const)&PPCDoubleDouble);
2015  assert (partCount()==2);
2016
2017  uint64_t myexponent, mysignificand, myexponent2, mysignificand2;
2018
2019  if (category==fcNormal) {
2020    myexponent = exponent + 1023; //bias
2021    myexponent2 = exponent2 + 1023;
2022    mysignificand = significandParts()[0];
2023    mysignificand2 = significandParts()[1];
2024    if (myexponent==1 && !(mysignificand & 0x10000000000000LL))
2025      myexponent = 0;   // denormal
2026    if (myexponent2==1 && !(mysignificand2 & 0x10000000000000LL))
2027      myexponent2 = 0;   // denormal
2028  } else if (category==fcZero) {
2029    myexponent = 0;
2030    mysignificand = 0;
2031    myexponent2 = 0;
2032    mysignificand2 = 0;
2033  } else if (category==fcInfinity) {
2034    myexponent = 0x7ff;
2035    myexponent2 = 0;
2036    mysignificand = 0;
2037    mysignificand2 = 0;
2038  } else {
2039    assert(category == fcNaN && "Unknown category");
2040    myexponent = 0x7ff;
2041    mysignificand = significandParts()[0];
2042    myexponent2 = exponent2;
2043    mysignificand2 = significandParts()[1];
2044  }
2045
2046  uint64_t words[2];
2047  words[0] =  (((uint64_t)sign & 1) << 63) |
2048              ((myexponent & 0x7ff) <<  52) |
2049              (mysignificand & 0xfffffffffffffLL);
2050  words[1] =  (((uint64_t)sign2 & 1) << 63) |
2051              ((myexponent2 & 0x7ff) <<  52) |
2052              (mysignificand2 & 0xfffffffffffffLL);
2053  return APInt(128, 2, words);
2054}
2055
2056APInt
2057APFloat::convertDoubleAPFloatToAPInt() const
2058{
2059  assert(semantics == (const llvm::fltSemantics*)&IEEEdouble);
2060  assert (partCount()==1);
2061
2062  uint64_t myexponent, mysignificand;
2063
2064  if (category==fcNormal) {
2065    myexponent = exponent+1023; //bias
2066    mysignificand = *significandParts();
2067    if (myexponent==1 && !(mysignificand & 0x10000000000000LL))
2068      myexponent = 0;   // denormal
2069  } else if (category==fcZero) {
2070    myexponent = 0;
2071    mysignificand = 0;
2072  } else if (category==fcInfinity) {
2073    myexponent = 0x7ff;
2074    mysignificand = 0;
2075  } else {
2076    assert(category == fcNaN && "Unknown category!");
2077    myexponent = 0x7ff;
2078    mysignificand = *significandParts();
2079  }
2080
2081  return APInt(64, (((((uint64_t)sign & 1) << 63) |
2082                     ((myexponent & 0x7ff) <<  52) |
2083                     (mysignificand & 0xfffffffffffffLL))));
2084}
2085
2086APInt
2087APFloat::convertFloatAPFloatToAPInt() const
2088{
2089  assert(semantics == (const llvm::fltSemantics*)&IEEEsingle);
2090  assert (partCount()==1);
2091
2092  uint32_t myexponent, mysignificand;
2093
2094  if (category==fcNormal) {
2095    myexponent = exponent+127; //bias
2096    mysignificand = *significandParts();
2097    if (myexponent == 1 && !(mysignificand & 0x400000))
2098      myexponent = 0;   // denormal
2099  } else if (category==fcZero) {
2100    myexponent = 0;
2101    mysignificand = 0;
2102  } else if (category==fcInfinity) {
2103    myexponent = 0xff;
2104    mysignificand = 0;
2105  } else {
2106    assert(category == fcNaN && "Unknown category!");
2107    myexponent = 0xff;
2108    mysignificand = *significandParts();
2109  }
2110
2111  return APInt(32, (((sign&1) << 31) | ((myexponent&0xff) << 23) |
2112                    (mysignificand & 0x7fffff)));
2113}
2114
2115// This function creates an APInt that is just a bit map of the floating
2116// point constant as it would appear in memory.  It is not a conversion,
2117// and treating the result as a normal integer is unlikely to be useful.
2118
2119APInt
2120APFloat::convertToAPInt() const
2121{
2122  if (semantics == (const llvm::fltSemantics* const)&IEEEsingle)
2123    return convertFloatAPFloatToAPInt();
2124
2125  if (semantics == (const llvm::fltSemantics* const)&IEEEdouble)
2126    return convertDoubleAPFloatToAPInt();
2127
2128  if (semantics == (const llvm::fltSemantics* const)&PPCDoubleDouble)
2129    return convertPPCDoubleDoubleAPFloatToAPInt();
2130
2131  assert(semantics == (const llvm::fltSemantics* const)&x87DoubleExtended &&
2132         "unknown format!");
2133  return convertF80LongDoubleAPFloatToAPInt();
2134}
2135
2136float
2137APFloat::convertToFloat() const
2138{
2139  assert(semantics == (const llvm::fltSemantics* const)&IEEEsingle);
2140  APInt api = convertToAPInt();
2141  return api.bitsToFloat();
2142}
2143
2144double
2145APFloat::convertToDouble() const
2146{
2147  assert(semantics == (const llvm::fltSemantics* const)&IEEEdouble);
2148  APInt api = convertToAPInt();
2149  return api.bitsToDouble();
2150}
2151
2152/// Integer bit is explicit in this format.  Current Intel book does not
2153/// define meaning of:
2154///  exponent = all 1's, integer bit not set.
2155///  exponent = 0, integer bit set. (formerly "psuedodenormals")
2156///  exponent!=0 nor all 1's, integer bit not set. (formerly "unnormals")
2157void
2158APFloat::initFromF80LongDoubleAPInt(const APInt &api)
2159{
2160  assert(api.getBitWidth()==80);
2161  uint64_t i1 = api.getRawData()[0];
2162  uint64_t i2 = api.getRawData()[1];
2163  uint64_t myexponent = (i1 >> 48) & 0x7fff;
2164  uint64_t mysignificand = ((i1 << 16) &  0xffffffffffff0000ULL) |
2165                          (i2 & 0xffff);
2166
2167  initialize(&APFloat::x87DoubleExtended);
2168  assert(partCount()==2);
2169
2170  sign = i1>>63;
2171  if (myexponent==0 && mysignificand==0) {
2172    // exponent, significand meaningless
2173    category = fcZero;
2174  } else if (myexponent==0x7fff && mysignificand==0x8000000000000000ULL) {
2175    // exponent, significand meaningless
2176    category = fcInfinity;
2177  } else if (myexponent==0x7fff && mysignificand!=0x8000000000000000ULL) {
2178    // exponent meaningless
2179    category = fcNaN;
2180    significandParts()[0] = mysignificand;
2181    significandParts()[1] = 0;
2182  } else {
2183    category = fcNormal;
2184    exponent = myexponent - 16383;
2185    significandParts()[0] = mysignificand;
2186    significandParts()[1] = 0;
2187    if (myexponent==0)          // denormal
2188      exponent = -16382;
2189  }
2190}
2191
2192void
2193APFloat::initFromPPCDoubleDoubleAPInt(const APInt &api)
2194{
2195  assert(api.getBitWidth()==128);
2196  uint64_t i1 = api.getRawData()[0];
2197  uint64_t i2 = api.getRawData()[1];
2198  uint64_t myexponent = (i1 >> 52) & 0x7ff;
2199  uint64_t mysignificand = i1 & 0xfffffffffffffLL;
2200  uint64_t myexponent2 = (i2 >> 52) & 0x7ff;
2201  uint64_t mysignificand2 = i2 & 0xfffffffffffffLL;
2202
2203  initialize(&APFloat::PPCDoubleDouble);
2204  assert(partCount()==2);
2205
2206  sign = i1>>63;
2207  sign2 = i2>>63;
2208  if (myexponent==0 && mysignificand==0) {
2209    // exponent, significand meaningless
2210    // exponent2 and significand2 are required to be 0; we don't check
2211    category = fcZero;
2212  } else if (myexponent==0x7ff && mysignificand==0) {
2213    // exponent, significand meaningless
2214    // exponent2 and significand2 are required to be 0; we don't check
2215    category = fcInfinity;
2216  } else if (myexponent==0x7ff && mysignificand!=0) {
2217    // exponent meaningless.  So is the whole second word, but keep it
2218    // for determinism.
2219    category = fcNaN;
2220    exponent2 = myexponent2;
2221    significandParts()[0] = mysignificand;
2222    significandParts()[1] = mysignificand2;
2223  } else {
2224    category = fcNormal;
2225    // Note there is no category2; the second word is treated as if it is
2226    // fcNormal, although it might be something else considered by itself.
2227    exponent = myexponent - 1023;
2228    exponent2 = myexponent2 - 1023;
2229    significandParts()[0] = mysignificand;
2230    significandParts()[1] = mysignificand2;
2231    if (myexponent==0)          // denormal
2232      exponent = -1022;
2233    else
2234      significandParts()[0] |= 0x10000000000000LL;  // integer bit
2235    if (myexponent2==0)
2236      exponent2 = -1022;
2237    else
2238      significandParts()[1] |= 0x10000000000000LL;  // integer bit
2239  }
2240}
2241
2242void
2243APFloat::initFromDoubleAPInt(const APInt &api)
2244{
2245  assert(api.getBitWidth()==64);
2246  uint64_t i = *api.getRawData();
2247  uint64_t myexponent = (i >> 52) & 0x7ff;
2248  uint64_t mysignificand = i & 0xfffffffffffffLL;
2249
2250  initialize(&APFloat::IEEEdouble);
2251  assert(partCount()==1);
2252
2253  sign = i>>63;
2254  if (myexponent==0 && mysignificand==0) {
2255    // exponent, significand meaningless
2256    category = fcZero;
2257  } else if (myexponent==0x7ff && mysignificand==0) {
2258    // exponent, significand meaningless
2259    category = fcInfinity;
2260  } else if (myexponent==0x7ff && mysignificand!=0) {
2261    // exponent meaningless
2262    category = fcNaN;
2263    *significandParts() = mysignificand;
2264  } else {
2265    category = fcNormal;
2266    exponent = myexponent - 1023;
2267    *significandParts() = mysignificand;
2268    if (myexponent==0)          // denormal
2269      exponent = -1022;
2270    else
2271      *significandParts() |= 0x10000000000000LL;  // integer bit
2272  }
2273}
2274
2275void
2276APFloat::initFromFloatAPInt(const APInt & api)
2277{
2278  assert(api.getBitWidth()==32);
2279  uint32_t i = (uint32_t)*api.getRawData();
2280  uint32_t myexponent = (i >> 23) & 0xff;
2281  uint32_t mysignificand = i & 0x7fffff;
2282
2283  initialize(&APFloat::IEEEsingle);
2284  assert(partCount()==1);
2285
2286  sign = i >> 31;
2287  if (myexponent==0 && mysignificand==0) {
2288    // exponent, significand meaningless
2289    category = fcZero;
2290  } else if (myexponent==0xff && mysignificand==0) {
2291    // exponent, significand meaningless
2292    category = fcInfinity;
2293  } else if (myexponent==0xff && mysignificand!=0) {
2294    // sign, exponent, significand meaningless
2295    category = fcNaN;
2296    *significandParts() = mysignificand;
2297  } else {
2298    category = fcNormal;
2299    exponent = myexponent - 127;  //bias
2300    *significandParts() = mysignificand;
2301    if (myexponent==0)    // denormal
2302      exponent = -126;
2303    else
2304      *significandParts() |= 0x800000; // integer bit
2305  }
2306}
2307
2308/// Treat api as containing the bits of a floating point number.  Currently
2309/// we infer the floating point type from the size of the APInt.  The
2310/// isIEEE argument distinguishes between PPC128 and IEEE128 (not meaningful
2311/// when the size is anything else).
2312void
2313APFloat::initFromAPInt(const APInt& api, bool isIEEE)
2314{
2315  if (api.getBitWidth() == 32)
2316    return initFromFloatAPInt(api);
2317  else if (api.getBitWidth()==64)
2318    return initFromDoubleAPInt(api);
2319  else if (api.getBitWidth()==80)
2320    return initFromF80LongDoubleAPInt(api);
2321  else if (api.getBitWidth()==128 && !isIEEE)
2322    return initFromPPCDoubleDoubleAPInt(api);
2323  else
2324    assert(0);
2325}
2326
2327APFloat::APFloat(const APInt& api, bool isIEEE)
2328{
2329  initFromAPInt(api, isIEEE);
2330}
2331
2332APFloat::APFloat(float f)
2333{
2334  APInt api = APInt(32, 0);
2335  initFromAPInt(api.floatToBits(f));
2336}
2337
2338APFloat::APFloat(double d)
2339{
2340  APInt api = APInt(64, 0);
2341  initFromAPInt(api.doubleToBits(d));
2342}
2343