#include "decimal.h"
#include "integer.h"
+#include "leading-zero-bit-count.h"
+#include "../common/idioms.h"
+#include <iostream> // TODO pmk rm
+bool pmk{false};
namespace Fortran::evaluate::value {
-static std::ostream *debug{nullptr};
+static std::ostream *debug{nullptr}; // pmk constexpr
template<typename REAL>
std::ostream &Decimal<REAL>::Dump(std::ostream &o) const {
return o << " e" << exponent_ << '\n';
}
-template<typename REAL> void Decimal<REAL>::FromReal(const Real &x) {
+template<typename REAL> void Decimal<REAL>::FromReal(const REAL &x) {
if (x.IsNegative()) {
FromReal(x.Negate());
isNegative_ = true;
return;
}
if (debug) {
- *debug << "FromReal(" << x.DumpHexadecimal() << ")\n";
+ *debug << "FromReal(" << x.DumpHexadecimal() << ") bits " << Real::bits
+ << '\n';
}
if (x.IsZero()) {
return;
if (debug) {
*debug << "second twoPow " << twoPow << ", lshift " << lshift << '\n';
}
- SetTo(x.GetFraction().SHIFTL(lshift));
+ using Word = typename Real::Word;
+ Word word{Word::ConvertUnsigned(x.GetFraction()).value};
+ SetTo(word.SHIFTL(lshift));
+ if (debug) {
+ Dump(*debug);
+ }
for (; twoPow > 0 && IsDivisibleBy<5>(); --twoPow) {
DivideBy<5>();
using Real = REAL;
using Word = typename Real::Word;
- void SetTo(Word n) {
- word_ = n;
+ void SetTo(std::uint64_t n) {
exponent_ = 0;
+ if constexpr (Word::bits >= 8 * sizeof n) {
+ word_ = n;
+ } else {
+ int shift{64 - LeadingZeroBitCount(n) - Word::bits};
+ if (shift <= 0) {
+ word_ = n;
+ } else {
+ word_ = n >> shift;
+ exponent_ = shift;
+ bool sticky{n << (64 - shift) != 0};
+ if (sticky) {
+ word_ = word_.IOR(Word{1});
+ }
+ }
+ }
}
void MultiplyAndAdd(std::uint32_t n, std::uint32_t plus = 0) {
sticky |= product.lower.BTEST(0);
product.lower = product.lower.DSHIFTR(product.upper, 1);
product.upper = product.upper.SHIFTR(1);
+ ++exponent_;
}
word_ = product.lower;
if (sticky) {
- word_ = word_.IOR(word_.MASKR(1));
+ word_ = word_.IOR(Word{1});
}
}
void AdjustExponent(int by = -1) { exponent_ += by; }
- Real ToReal(bool isNegative = false) const;
+ ValueWithRealFlags<Real> ToReal(
+ bool isNegative = false, Rounding rounding = Rounding::TiesToEven) const;
private:
Word word_{0};
}
template<typename REAL>
-REAL IntermediateFloat<REAL>::ToReal(bool isNegative) const {
+ValueWithRealFlags<REAL> IntermediateFloat<REAL>::ToReal(
+ bool isNegative, Rounding rounding) const {
if (word_.IsNegative()) {
+ // word_ represents an unsigned quantity, so shift it down if the MSB is set
IntermediateFloat shifted;
- shifted.word_ =
- word_.SHIFTR(1).IOR(word_.IAND(word_.MASKR(1))); // sticky bit
+ Word sticky{word_.IAND(Word{1})};
+ shifted.word_ = word_.SHIFTR(1).IOR(sticky);
shifted.exponent_ = exponent_ + 1;
- return shifted.ToReal(false);
+ if (debug) {
+ shifted.Dump(*debug << "IntermediateFloat::ToReal: shifted: ") << '\n';
+ }
+ return shifted.ToReal(isNegative, rounding);
+ }
+ ValueWithRealFlags<Real> result;
+ if (isNegative) {
+ result = Real::FromInteger(word_.Negate().value, rounding);
+ } else {
+ result = Real::FromInteger(word_, rounding);
+ }
+ if (debug) {
+ *debug << "IntermediateFloat::ToReal: after FromInteger: "
+ << result.value.DumpHexadecimal() << " * 2**" << exponent_ << '\n';
}
- Real result = Real::FromInteger(word_).value;
int expo{exponent_};
while (expo + Real::exponentBias < 1) {
- Real twoPow{Word{1}.SHIFTL(Real::significandBits)};
- result = result.Multiply(twoPow).value;
+ Real twoPow{Word{1}.SHIFTL(Real::significandBits)}; // min normal value
+ result.value = result.value.Multiply(twoPow).AccumulateFlags(result.flags);
expo += Real::exponentBias - 1;
+ if (debug) {
+ *debug << "IntermediateFloat::ToReal: reduced: "
+ << result.value.DumpHexadecimal() << " * 2**" << expo << '\n';
+ }
}
while (expo + Real::exponentBias >= Real::maxExponent) {
Real twoPow{Word{Real::maxExponent - 1}.SHIFTL(Real::significandBits)};
- result = result.Multiply(twoPow).value;
+ result.value = result.value.Multiply(twoPow).AccumulateFlags(result.flags);
expo += Real::maxExponent - 1 - Real::exponentBias;
+ if (debug) {
+ *debug << "IntermediateFloat::ToReal: magnified: "
+ << result.value.DumpHexadecimal() << " * 2**" << expo << '\n';
+ }
}
Real twoPow{Word{expo + Real::exponentBias}.SHIFTL(Real::significandBits)};
- if (isNegative) {
- twoPow = twoPow.Negate();
+ if (debug) {
+ *debug << "IntermediateFloat::ToReal: twoPow: " << twoPow.DumpHexadecimal()
+ << '\n';
}
- return result.Multiply(twoPow).value;
+ result.value = result.value.Multiply(twoPow).AccumulateFlags(result.flags);
+ return result;
}
-template<typename REAL> REAL Decimal<REAL>::ToReal(const char *&p) {
+template<typename REAL>
+ValueWithRealFlags<REAL> Decimal<REAL>::ToReal(
+ const char *&p, Rounding rounding) {
+ if (std::string{p} == "3.4028234663852885981170418348451692544e38_4") {
+ debug = &std::cerr;
+ pmk = true;
+ } else {
+ debug = nullptr;
+ pmk = false;
+ } // pmk
+ if (debug) {
+ *debug << "ToReal('" << p << "')\n";
+ }
while (*p == ' ') {
++p;
}
++exponent_;
}
} else {
- MultiplyBy<10>(c - '0');
+ int carry{MultiplyBy<10>(c - '0')};
+ CHECK(carry == 0);
if (decimalPoint) {
--exponent_;
}
}
+ if (debug) Dump(*debug << "ToReal in loop, p at '" << p << "'\n'");
}
switch (*p) {
case 'D':
case 'q':
case 'Q':
- ++p;
- bool negExpo{*p == '-'};
- if (*p == '+' || *p == '-') {
+ bool negExpo{*++p == '-'};
+ if (negExpo || *p == '+') {
++p;
}
char *q;
}
if (debug) {
- Dump(*debug << "ToReal start: ");
+ Dump(*debug << "ToReal start, p at '" << p << "'\n");
}
if (IsZero()) {
- return Real{}.Negate(); // -0.0
+ ValueWithRealFlags<Real> result;
+ if (isNegative_) {
+ result.value = Real{}.Negate(); // -0.0
+ }
+ return result;
}
if (exponent_ < 0) {
// result. The most significant digit can be moved directly;
// lesser-order digits require transfer of carries.
if (exponent_ >= -(digits_ - 1) * log10Quintillion) {
+ if (debug) {
+ Dump(*debug << "converting integer part:");
+ }
f.SetTo(digit_[--digits_]);
+ if (debug) {
+ f.Dump(*debug << "after converting top digit: ") << '\n';
+ Dump(*debug);
+ }
while (exponent_ > -digits_ * log10Quintillion) {
digitLimit_ = digits_;
int carry{MultiplyBy<10>()};
f.MultiplyAndAdd(10, carry);
--exponent_;
+ if (debug) {
+ f.Dump(*debug << "foot of loop after carry " << carry << ": ") << '\n';
+ Dump(*debug);
+ }
}
}
Dump(*debug);
}
- return f.ToReal(isNegative_);
+ auto result{f.ToReal(isNegative_, rounding)};
+ return result;
}
template<typename REAL>
#ifndef FORTRAN_EVALUATE_DECIMAL_H_
#define FORTRAN_EVALUATE_DECIMAL_H_
+#include "common.h"
+#include "integer.h"
#include "real.h"
#include <cinttypes>
#include <limits>
first_ = 0;
exponent_ = 0;
}
+
void FromReal(const Real &);
- Real ToReal(const char *&); // arg left pointing to first unparsed char
+
+ // Convert a character representation of a floating-point value to
+ // the underlying Real type. The reference argument is a pointer that
+ // is left pointing to the first character that wasn't included.
+ ValueWithRealFlags<Real> ToReal(
+ const char *&, Rounding rounding = Rounding::TiesToEven);
+
std::string ToString(int maxDigits = 1000000) const;
private:
++exponent_;
n = qr.quotient;
}
- while (!n.IsZero() && digits_ < digitLimit_) {
- auto qr{n.DivideUnsigned(quintillion)};
- digit_[digits_++] = qr.remainder.ToUInt64();
- if (digits_ == first_ + 1 && digit_[first_] == 0) {
- ++first_;
+ if constexpr (INT::bits < 60) {
+ // n is necessarily less than a quintillion
+ if (!n.IsZero()) {
+ digit_[digits_++] = n.ToUInt64();
}
- n = qr.quotient;
+ return 0;
+ } else {
+ while (!n.IsZero() && digits_ < digitLimit_) {
+ auto qr{n.DivideUnsigned(quintillion)};
+ digit_[digits_++] = qr.remainder.ToUInt64();
+ if (digits_ == first_ + 1 && digit_[first_] == 0) {
+ ++first_;
+ }
+ n = qr.quotient;
+ }
+ return n;
}
- return n;
}
int RemoveLeastOrderZeroDigits() {
// limitations under the License.
#include "real.h"
+#include "decimal.h"
#include "int-power.h"
#include "../common/idioms.h"
#include "../parser/characters.h"
result.value = y; // x + +/-Inf -> +/-Inf
return result;
}
- std::uint64_t exponent{Exponent()};
- std::uint64_t yExponent{y.Exponent()};
+ int exponent{Exponent()};
+ int yExponent{y.Exponent()};
if (exponent < yExponent) {
// y is larger in magnitude; simplify by reversing operands
return y.Add(*this, rounding);
}
template<typename W, int P, bool IM>
-RealFlags Real<W, P, IM>::Normalize(bool negative, std::uint64_t exponent,
+RealFlags Real<W, P, IM>::Normalize(bool negative, int exponent,
const Fraction &fraction, Rounding rounding, RoundingBits *roundingBits) {
- std::uint64_t lshift = fraction.LEADZ();
+ int lshift{fraction.LEADZ()};
if (lshift == fraction.bits /* fraction is zero */ &&
(roundingBits == nullptr || roundingBits->empty())) {
// No fraction, no rounding bits -> +/-0.0
template<typename W, int P, bool IM>
RealFlags Real<W, P, IM>::Round(
Rounding rounding, const RoundingBits &bits, bool multiply) {
- std::uint64_t origExponent{Exponent()};
+ int origExponent{Exponent()};
RealFlags flags;
bool inexact{!bits.empty()};
if (inexact) {
bits.MustRound(rounding, IsNegative(), word_.BTEST(0) /* is odd */)) {
typename Fraction::ValueWithCarry sum{
GetFraction().AddUnsigned(Fraction{}, true)};
- std::uint64_t newExponent{origExponent};
+ int newExponent{origExponent};
if (sum.carry) {
// The fraction was all ones before rounding; sum.value is now zero
sum.value = sum.value.IBSET(precision - 1);
template<typename W, int P, bool IM>
void Real<W, P, IM>::NormalizeAndRound(ValueWithRealFlags<Real> &result,
- bool isNegative, std::uint64_t exponent, const Fraction &fraction,
- Rounding rounding, RoundingBits roundingBits, bool multiply) {
+ bool isNegative, int exponent, const Fraction &fraction, Rounding rounding,
+ RoundingBits roundingBits, bool multiply) {
result.flags |= result.value.Normalize(
isNegative, exponent, fraction, rounding, &roundingBits);
result.flags |= result.value.Round(rounding, roundingBits, multiply);
template<typename W, int P, bool IM>
ValueWithRealFlags<Real<W, P, IM>> Real<W, P, IM>::Read(
const char *&p, Rounding rounding) {
- ValueWithRealFlags<Real> result;
- while (*p == ' ') {
- ++p;
- }
- bool isNegative{*p == '-'};
- if (*p == '-' || *p == '+') {
- ++p;
- }
- Word integer{0};
- int decimalExponent{0};
- bool full{false};
- bool inFraction{false};
- for (;; ++p) {
- if (*p == '.') {
- if (inFraction) {
- break;
- }
- inFraction = true;
- } else if (!parser::IsDecimalDigit(*p)) {
- break;
- } else if (full) {
- if (!inFraction) {
- ++decimalExponent;
- }
- } else {
- auto times10{integer.MultiplyUnsigned(Word{10})};
- if (!times10.upper.IsZero()) {
- full = true;
- if (!inFraction) {
- ++decimalExponent;
- }
- } else {
- auto augmented{times10.lower.AddUnsigned(Word{*p - '0'})};
- if (augmented.carry) {
- full = true;
- if (!inFraction) {
- ++decimalExponent;
- }
- } else {
- integer = augmented.value;
- if (inFraction) {
- --decimalExponent;
- }
- }
- }
- }
- }
-
- if (parser::IsLetter(*p)) {
- ++p;
- bool negExpo{*p == '-'};
- if (*p == '+' || *p == '-') {
- ++p;
- }
- auto expo{Integer<32>::ReadUnsigned(p)};
- std::int64_t expoVal{expo.value.ToInt64()};
- if (expo.overflow) {
- expoVal = std::numeric_limits<std::int32_t>::max();
- } else if (negExpo) {
- expoVal *= -1;
- }
- decimalExponent += expoVal;
- }
-
- int binaryExponent{exponentBias + bits - 1};
- if (integer.IsZero()) {
- decimalExponent = 0;
- } else {
- int leadz{integer.LEADZ()};
- binaryExponent -= leadz;
- integer = integer.SHIFTL(leadz);
- }
- for (; decimalExponent > 0; --decimalExponent) {
- auto times5{integer.MultiplyUnsigned(Word{5})};
- ++binaryExponent;
- integer = times5.lower;
- for (; !times5.upper.IsZero(); times5.upper = times5.upper.SHIFTR(1)) {
- ++binaryExponent;
- integer = integer.SHIFTR(1);
- if (times5.upper.BTEST(0)) {
- integer = integer.IBSET(bits - 1);
- }
- }
- }
- for (; decimalExponent < 0; ++decimalExponent) {
- auto div5{integer.DivideUnsigned(Word{5})};
- --binaryExponent;
- integer = div5.quotient;
- std::uint8_t lost = div5.remainder.ToUInt64() * 0x33;
- while (!integer.BTEST(bits - 1)) {
- integer = integer.SHIFTL(1);
- if (lost & 0x80) {
- integer = integer.IBSET(0);
- }
- lost <<= 1;
- --binaryExponent;
- }
- }
-
- RoundingBits roundingBits;
- for (int j{0}; bits - j > precision; ++j) {
- roundingBits.ShiftRight(integer.BTEST(0));
- integer = integer.SHIFTR(1);
- }
-
- Fraction fraction{Fraction::ConvertUnsigned(integer).value};
- while (binaryExponent < 1) {
- if (fraction.IsZero()) {
- binaryExponent = 0;
- break;
- } else {
- ++binaryExponent;
- roundingBits.ShiftRight(fraction.BTEST(0));
- fraction = fraction.SHIFTR(1);
- }
- }
-
- NormalizeAndRound(
- result, isNegative, binaryExponent, fraction, rounding, roundingBits);
- return result;
+ Decimal<Real> decimal;
+ return decimal.ToReal(p, rounding);
}
template<typename W, int P, bool IM>
template<typename W, int P, bool IM>
std::ostream &Real<W, P, IM>::AsFortran(std::ostream &o, int kind) const {
- ValueWithRealFlags<ScaledDecimal> scaled{AsScaledDecimal()};
- if (scaled.flags.test(RealFlag::InvalidArgument)) {
+ if (IsNotANumber()) {
o << "(0._" << kind << "/0.)";
- } else if (scaled.flags.test(RealFlag::Overflow)) {
+ } else if (IsInfinite()) {
if (IsNegative()) {
o << "(-1._" << kind << "/0.)";
} else {
o << "(1._" << kind << "/0.)";
}
} else {
- if (scaled.value.negative) {
- o << "(-";
- }
- std::string digits{scaled.value.integer.UnsignedDecimal()};
- int exponent = scaled.value.decimalExponent + digits.size() - 1;
- o << digits[0] << '.' << digits.substr(1);
- if (exponent != 0) {
- o << 'e' << exponent;
+ bool parenthesize{IsNegative()};
+ if (parenthesize) {
+ o << "(";
}
+ Decimal<Real> decimal;
+ decimal.FromReal(*this);
+ o << decimal.ToString();
o << '_' << kind;
- if (scaled.value.negative) {
+ if (parenthesize) {
o << ')';
}
}
return o;
}
-template<typename W, int P, bool IM>
-auto Real<W, P, IM>::AsScaledDecimal(Rounding rounding) const
- -> ValueWithRealFlags<ScaledDecimal> {
-
- ValueWithRealFlags<ScaledDecimal> result;
- if (IsNotANumber()) {
- result.flags.set(RealFlag::InvalidArgument);
- return result;
- }
- if (IsInfinite()) {
- result.flags.set(RealFlag::Overflow);
- result.value.integer = Word::HUGE();
- return result;
- }
- if (IsNegative()) {
- if (rounding == Rounding::Up) {
- rounding = Rounding::Down;
- } else if (rounding == Rounding::Down) {
- rounding = Rounding::Up;
- }
- result = Negate().AsScaledDecimal(rounding);
- result.value.negative = true;
- return result;
- }
- if (IsZero()) {
- return result;
- }
-
- Word fraction{Word::ConvertUnsigned(GetFraction()).value};
- Word asInt{fraction.SHIFTL(bits - precision)};
- int twoPower{UnbiasedExponent() - bits + 1};
-
- // The original Real, a finite positive number, is now represented as a
- // left-justified integer and a binary exponent. In other words, asInt
- // has its sign bit set and the value of "asInt * 2.0**twoPower" equals
- // the original Real exactly; both asInt and twoPower are integers.
- //
- // To generate the scaled decimal result, we multiply or divide asInt by
- // 2**twoPower iteratively. Whenever the result of a multiplication or
- // division would overflow (or lose precision, respectively), we scale
- // asInt by dividing it by ten (or multiplying, respectively). Since
- // 10==2*5, the scaling is actually by a factor of five since the factor
- // of two can be accommodated by shifting. These scalings are recorded
- // in a decimal exponent. Once all of the "2.0**twoPower" scaling is
- // done, we have a value that is represented as
- // "asInt * 10.0**decimalExponent", and that is the penultimate result.
-
- static constexpr Word five{5};
- if (twoPower > 0) {
- // Multiply asInt by 2**twoPower.
- static constexpr Word two{2};
- for (; twoPower > 0; --twoPower) {
- if (asInt.BTEST(bits - 1)) {
- // asInt * 2 would overflow: divide by five and increment the
- // decimal exponent (effectively dividing by ten and then multiplying
- // by two).
- auto qr{asInt.DivideUnsigned(five)};
- asInt = qr.quotient;
- ++result.value.decimalExponent;
- if (twoPower > 1 &&
- qr.remainder.CompareUnsigned(two) == Ordering::Greater) {
- --twoPower;
- asInt = asInt.SHIFTL(1).IBSET(0);
- }
- } else {
- asInt = asInt.SHIFTL(1);
- }
- }
- } else {
- // Divide asInt by 2**(-twoPower).
- std::uint32_t lower{0};
- for (; twoPower < 0; ++twoPower) {
- auto times5{asInt.MultiplyUnsigned(five)};
- if (!times5.upper.IsZero()) {
- // asInt is too big to need scaling, just shift it down.
- lower >>= 1;
- if (asInt.BTEST(0)) {
- lower |= 1 << 31;
- }
- asInt = asInt.SHIFTR(1);
- } else {
- std::uint64_t lowerTimes5{lower * static_cast<std::uint64_t>(5)};
- std::uint32_t round = lowerTimes5 >> 32;
- auto rounded{times5.lower.AddUnsigned(Word{round})};
- if (rounded.carry) {
- // asInt is still too big to need scaling (rounding would overflow)
- lower >>= 1;
- if (asInt.BTEST(0)) {
- lower |= 1 << 31;
- }
- asInt = asInt.SHIFTR(1);
- } else {
- // asInt is small enough to be scaled; do so.
- --result.value.decimalExponent;
- lower = lowerTimes5;
- asInt = rounded.value;
- }
- }
- }
- }
-
- // Canonicalize to the minimum integer value: the final result will not
- // be a multiple of ten.
- result.value.integer = asInt;
- if (!result.value.integer.IsZero()) {
- while (true) {
- auto qr{result.value.integer.DivideUnsigned(10)};
- if (!qr.remainder.IsZero()) {
- break;
- }
- ++result.value.decimalExponent;
- result.value.integer = qr.quotient;
- }
- }
-
- return result;
-}
-
template class Real<Integer<16>, 11>;
template class Real<Integer<32>, 24>;
template class Real<Integer<64>, 53>;
#include <ostream>
#include <string>
+#include <iostream> // TODO pmk rm
+extern bool pmk;
+
// Some environments, viz. clang on Darwin, allow the macro HUGE
// to leak out of <math.h> even when it is never directly included.
#undef HUGE
template<typename WORD, int PRECISION, bool IMPLICIT_MSB = true> class Real {
public:
using Word = WORD;
+ static constexpr int bits{Word::bits};
static constexpr int precision{PRECISION};
using Fraction = Integer<precision>; // all bits made explicit
-
- static constexpr int bits{Word::bits};
static constexpr bool implicitMSB{IMPLICIT_MSB};
static constexpr int significandBits{precision - implicitMSB};
static constexpr int exponentBits{bits - significandBits - 1 /*sign*/};
static_assert(precision > 0);
static_assert(exponentBits > 1);
- static constexpr std::uint64_t maxExponent{(1 << exponentBits) - 1};
- static constexpr std::uint64_t exponentBias{maxExponent / 2};
+ static_assert(exponentBits <= 16);
+ static constexpr int maxExponent{(1 << exponentBits) - 1};
+ static constexpr int exponentBias{maxExponent / 2};
// Decimal precision of a binary floating-point representation is
// actually the same as the base-5 precision, as factors of two
// Calculate "precision*0.43" with integer arithmetic so as to be constexpr.
static constexpr int decimalDigits{(precision * 43) / 100};
- // Associates a decimal exponent with an integral value for formatting.
- // TODO pmk remove
- struct ScaledDecimal {
- bool negative{false};
- Word integer; // unsigned
- int decimalExponent{0}; // Exxx
- };
-
template<typename W, int P, bool I> friend class Real;
constexpr Real() {} // +0.0
const Real &, Rounding rounding = defaultRounding) const;
template<typename INT> constexpr INT EXPONENT() const {
- std::uint64_t exponent{Exponent()};
+ int exponent{Exponent()};
if (exponent == maxExponent) {
return INT::HUGE();
} else {
- int result = exponent - exponentBias;
+ int result{exponent - exponentBias};
if (IsDenormal()) {
++result;
}
return {}; // all bits zero -> +0.0
}
ValueWithRealFlags<Real> result;
- std::uint64_t exponent{exponentBias + absN.bits - leadz - 1};
+ int exponent{exponentBias + absN.bits - leadz - 1};
int bitsNeeded{absN.bits - (leadz + implicitMSB)};
int bitsLost{bitsNeeded - significandBits};
+ if (pmk)
+ std::cerr << "pmk real.h exponent " << exponent << " bitsLost "
+ << bitsLost << " rounding " << static_cast<int>(rounding)
+ << '\n';
if (bitsLost <= 0) {
Fraction fraction{Fraction::ConvertUnsigned(absN).value};
result.flags |= result.value.Normalize(
} else {
Fraction fraction{Fraction::ConvertUnsigned(absN.SHIFTR(bitsLost)).value};
result.flags |= result.value.Normalize(isNegative, exponent, fraction);
+ if (pmk)
+ std::cerr << "pmk Normalized " << result.value.DumpHexadecimal()
+ << '\n';
RoundingBits roundingBits{absN, bitsLost};
result.flags |= result.value.Round(rounding, roundingBits);
+ if (pmk)
+ std::cerr << "pmk Rounded " << result.value.DumpHexadecimal() << '\n';
}
return result;
}
} else if (IsInfinite()) {
result.flags.set(RealFlag::Overflow);
} else {
- std::uint64_t exponent{Exponent()};
+ int exponent{Exponent()};
if (exponent < exponentBias) { // |x| < 1.0
result.value.Normalize(IsNegative(), 0, Fraction{}); // +/-0.0
} else {
- constexpr std::uint64_t noClipExponent{exponentBias + precision - 1};
+ constexpr int noClipExponent{exponentBias + precision - 1};
if (int clip = noClipExponent - exponent; clip > 0) {
result.value.word_ = result.value.word_.IAND(Word::MASKR(clip).NOT());
}
return result;
}
bool isNegative{IsNegative()};
- std::uint64_t exponent{Exponent()};
+ int exponent{Exponent()};
Fraction fraction{GetFraction()};
if (exponent >= maxExponent || // +/-Inf
exponent >= exponentBias + result.value.bits) { // too big
}
} else {
// finite number |x| >= 1.0
- constexpr std::uint64_t noShiftExponent{exponentBias + precision - 1};
+ constexpr int noShiftExponent{exponentBias + precision - 1};
if (exponent < noShiftExponent) {
int rshift = noShiftExponent - exponent;
if (!fraction.IBITS(0, rshift).IsZero()) {
absX = x.Negate();
}
ValueWithRealFlags<Real> result;
- std::uint64_t exponent{exponentBias + x.Exponent() - A::exponentBias};
+ int exponent{exponentBias + x.Exponent() - A::exponentBias};
int bitsLost{A::precision - precision};
typename A::Fraction xFraction{x.GetFraction()};
if (bitsLost <= 0) {
return result;
}
- // Represents the number as "J*(10**K)" where J and K are integers.
- ValueWithRealFlags<ScaledDecimal> AsScaledDecimal(
- Rounding rounding = defaultRounding) const;
-
constexpr Word RawBits() const { return word_; }
// Extracts "raw" biased exponent field.
- constexpr std::uint64_t Exponent() const {
+ constexpr int Exponent() const {
return word_.IBITS(significandBits, exponentBits).ToUInt64();
}
- // Extracts unbiased exponent value.
- constexpr int UnbiasedExponent() const {
- int exponent =
- word_.IBITS(significandBits, exponentBits).ToUInt64() - exponentBias;
- if (IsDenormal()) {
- ++exponent;
- }
- return exponent;
- }
-
// Extracts the fraction; any implied bit is made explicit.
constexpr Fraction GetFraction() const {
Fraction result{Fraction::ConvertUnsigned(word_).value};
if constexpr (!implicitMSB) {
return result;
} else {
- std::uint64_t exponent{Exponent()};
+ int exponent{Exponent()};
if (exponent > 0 && exponent < maxExponent) {
return result.IBSET(significandBits);
} else {
}
}
+ // Extracts unbiased exponent value.
+ constexpr int UnbiasedExponent() const {
+ int exponent =
+ word_.IBITS(significandBits, exponentBits).ToUInt64() - exponentBias;
+ if (IsDenormal()) {
+ ++exponent;
+ }
+ return exponent;
+ }
+
static ValueWithRealFlags<Real> Read(
const char *&, Rounding rounding = defaultRounding);
std::string DumpHexadecimal() const;
return Significand::ConvertUnsigned(word_).value;
}
- constexpr std::int64_t CombineExponents(const Real &y, bool forDivide) const {
- std::int64_t exponent = Exponent(), yExponent = y.Exponent();
+ constexpr int CombineExponents(const Real &y, bool forDivide) const {
+ int exponent = Exponent(), yExponent = y.Exponent();
// A zero exponent field value has the same weight as 1.
exponent += !exponent;
yExponent += !yExponent;
// The value is a number, and a zero fraction means a zero value (i.e.,
// a maximal exponent and zero fraction doesn't signify infinity, although
// this member function will detect overflow and encode infinities).
- RealFlags Normalize(bool negative, std::uint64_t exponent,
- const Fraction &fraction, Rounding rounding = defaultRounding,
+ RealFlags Normalize(bool negative, int exponent, const Fraction &fraction,
+ Rounding rounding = defaultRounding,
RoundingBits *roundingBits = nullptr);
// Rounds a result, if necessary, in place.
RealFlags Round(Rounding, const RoundingBits &, bool multiply = false);
static void NormalizeAndRound(ValueWithRealFlags<Real> &result,
- bool isNegative, std::uint64_t exponent, const Fraction &, Rounding,
- RoundingBits, bool multiply = false);
+ bool isNegative, int exponent, const Fraction &, Rounding, RoundingBits,
+ bool multiply = false);
Word word_{}; // an Integer<>
};
#include <cmath>
#include <cstdio>
#include <cstdlib>
+#include <sstream>
#include <type_traits>
using namespace Fortran::evaluate;
return (x & 0x7fffffffffffffff) == 0x7ff0000000000000;
}
+inline bool IsNegative(std::uint32_t x) { return (x & 0x80000000) != 0; }
+
+inline bool IsNegative(std::uint64_t x) {
+ return (x & 0x8000000000000000) != 0;
+}
+
inline std::uint32_t NormalizeNaN(std::uint32_t x) {
if (IsNaN(x)) {
x = 0x7fe00000;
MATCH(IsInfinite(rj), x.IsInfinite())
("%d IsInfinite(0x%llx)", pass, static_cast<long long>(rj));
- if (rounding == Rounding::TiesToEven) {
- auto scaled{x.AsScaledDecimal()};
+ if (rounding == Rounding::TiesToEven ||
+ rounding == Rounding::ToZero) { // pmk
+ int kind{REAL::bits / 8};
+ std::stringstream ss, css;
+ x.AsFortran(ss, kind);
+ std::string s{ss.str()};
if (IsNaN(rj)) {
- TEST(scaled.flags.test(RealFlag::InvalidArgument))
+ css << "(0._" << kind << "/0.)";
+ MATCH(css.str(), s)
("%d invalid(0x%llx)", pass, static_cast<long long>(rj));
} else if (IsInfinite(rj)) {
- TEST(scaled.flags.test(RealFlag::Overflow))
+ css << '(';
+ if (IsNegative(rj)) {
+ css << '-';
+ }
+ css << "1._" << kind << "/0.)";
+ MATCH(css.str(), s)
("%d overflow(0x%llx)", pass, static_cast<long long>(rj));
} else {
- auto integer{scaled.value.integer.ToUInt64()};
- MATCH(x.IsNegative(), scaled.value.negative)
- ("%d IsNegative(0x%llx)", pass, static_cast<long long>(rj));
- char buffer[128];
- const char *p = buffer;
- snprintf(buffer, sizeof buffer, "%c%llu.0E%d",
- "+-"[scaled.value.negative],
- static_cast<unsigned long long>(integer),
- scaled.value.decimalExponent);
+ const char *p = s.data();
+ if (*p == '(') {
+ ++p;
+ }
auto readBack{REAL::Read(p, rounding)};
MATCH(rj, readBack.value.RawBits().ToUInt64())
- ("%d scaled decimal 0x%llx %s", pass, static_cast<long long>(rj),
- buffer);
+ ("%d Read(AsFortran()) 0x%llx %s %g", pass,
+ static_cast<long long>(rj), s.data(), static_cast<double>(fj));
+ MATCH('_', *p)
+ ("%d Read(AsFortran()) 0x%llx %s %d", pass,
+ static_cast<long long>(rj), s.data(),
+ static_cast<int>(p - s.data()));
}
}
}
}
}
+FailureDetailPrinter Match(const char *file, int line, const std::string &want,
+ const char *gots, const std::string &got) {
+ return Match(file, line, want.data(), gots, got);
+}
+
FailureDetailPrinter Compare(const char *file, int line, const char *xs,
const char *rel, const char *ys, unsigned long long x,
unsigned long long y) {
const char *gots, unsigned long long got);
FailureDetailPrinter Match(const char *file, int line, const char *want,
const char *gots, const std::string &got);
+FailureDetailPrinter Match(const char *file, int line, const std::string &want,
+ const char *gots, const std::string &got);
FailureDetailPrinter Compare(const char *file, int line, const char *xs,
const char *rel, const char *ys, unsigned long long x,
unsigned long long y);