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