//===-- lib/Decimal/big-radix-floating-point.h ------------------*- C++ -*-===//
//
// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
// See https://llvm.org/LICENSE.txt for license information.
// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
//
//===----------------------------------------------------------------------===//
#ifndef FORTRAN_DECIMAL_BIG_RADIX_FLOATING_POINT_H_
#define FORTRAN_DECIMAL_BIG_RADIX_FLOATING_POINT_H_
// This is a helper class for use in floating-point conversions between
// binary and decimal representations. It holds a multiple-precision
// integer value using digits of a radix that is a large even power of ten
// (10,000,000,000,000,000 by default, 10**16). These digits are accompanied
// by a signed exponent that denotes multiplication by a power of ten.
// The effective radix point is to the right of the digits (i.e., they do
// not represent a fraction).
//
// The operations supported by this class are limited to those required
// for conversions between binary and decimal representations; it is not
// a general-purpose facility.
#include "flang/Common/bit-population-count.h"
#include "flang/Common/leading-zero-bit-count.h"
#include "flang/Common/uint128.h"
#include "flang/Decimal/binary-floating-point.h"
#include "flang/Decimal/decimal.h"
#include <cinttypes>
#include <limits>
#include <type_traits>
// Some environments, viz. glibc 2.17, allow the macro HUGE
// to leak out of <math.h>.
#undef HUGE
namespace Fortran::decimal {
static constexpr std::uint64_t TenToThe(int power) {
return power <= 0 ? 1 : 10 * TenToThe(power - 1);
}
// 10**(LOG10RADIX + 3) must be < 2**wordbits, and LOG10RADIX must be
// even, so that pairs of decimal digits do not straddle Digits.
// So LOG10RADIX must be 16 or 6.
template <int PREC, int LOG10RADIX = 16> class BigRadixFloatingPointNumber {
public:
using Real = BinaryFloatingPointNumber<PREC>;
static constexpr int log10Radix{LOG10RADIX};
private:
static constexpr std::uint64_t uint64Radix{TenToThe(log10Radix)};
static constexpr int minDigitBits{
64 - common::LeadingZeroBitCount(uint64Radix)};
using Digit = common::HostUnsignedIntType<minDigitBits>;
static constexpr Digit radix{uint64Radix};
static_assert(radix < std::numeric_limits<Digit>::max() / 1000,
"radix is somehow too big");
static_assert(radix > std::numeric_limits<Digit>::max() / 10000,
"radix is somehow too small");
// The base-2 logarithm of the least significant bit that can arise
// in a subnormal IEEE floating-point number.
static constexpr int minLog2AnyBit{
-Real::exponentBias - Real::binaryPrecision};
// The number of Digits needed to represent the smallest subnormal.
static constexpr int maxDigits{3 - minLog2AnyBit / log10Radix};
public:
explicit RT_API_ATTRS BigRadixFloatingPointNumber(
enum FortranRounding rounding = RoundNearest)
: rounding_{rounding} {}
// Converts a binary floating point value.
explicit RT_API_ATTRS BigRadixFloatingPointNumber(
Real, enum FortranRounding = RoundNearest);
RT_API_ATTRS BigRadixFloatingPointNumber &SetToZero() {
isNegative_ = false;
digits_ = 0;
exponent_ = 0;
return *this;
}
RT_API_ATTRS bool IsInteger() const { return exponent_ >= 0; }
// Converts decimal floating-point to binary.
RT_API_ATTRS ConversionToBinaryResult<PREC> ConvertToBinary();
// Parses and converts to binary. Handles leading spaces,
// "NaN", & optionally-signed "Inf". Does not skip internal
// spaces.
// The argument is a reference to a pointer that is left
// pointing to the first character that wasn't parsed.
RT_API_ATTRS ConversionToBinaryResult<PREC> ConvertToBinary(
const char *&, const char *end = nullptr);
// Formats a decimal floating-point number to a user buffer.
// May emit "NaN" or "Inf", or an possibly-signed integer.
// No decimal point is written, but if it were, it would be
// after the last digit; the effective decimal exponent is
// returned as part of the result structure so that it can be
// formatted by the client.
RT_API_ATTRS ConversionToDecimalResult ConvertToDecimal(
char *, std::size_t, enum DecimalConversionFlags, int digits) const;
// Discard decimal digits not needed to distinguish this value
// from the decimal encodings of two others (viz., the nearest binary
// floating-point numbers immediately below and above this one).
// The last decimal digit may not be uniquely determined in all
// cases, and will be the mean value when that is so (e.g., if
// last decimal digit values 6-8 would all work, it'll be a 7).
// This minimization necessarily assumes that the value will be
// emitted and read back into the same (or less precise) format
// with default rounding to the nearest value.
RT_API_ATTRS void Minimize(
BigRadixFloatingPointNumber &&less, BigRadixFloatingPointNumber &&more);
template <typename STREAM> STREAM &Dump(STREAM &) const;
private:
RT_API_ATTRS BigRadixFloatingPointNumber(
const BigRadixFloatingPointNumber &that)
: digits_{that.digits_}, exponent_{that.exponent_},
isNegative_{that.isNegative_}, rounding_{that.rounding_} {
for (int j{0}; j < digits_; ++j) {
digit_[j] = that.digit_[j];
}
}
RT_API_ATTRS bool IsZero() const {
// Don't assume normalization.
for (int j{0}; j < digits_; ++j) {
if (digit_[j] != 0) {
return false;
}
}
return true;
}
// Predicate: true when 10*value would cause a carry.
// (When this happens during decimal-to-binary conversion,
// there are more digits in the input string than can be
// represented precisely.)
RT_API_ATTRS bool IsFull() const {
return digits_ == digitLimit_ && digit_[digits_ - 1] >= radix / 10;
}
// Sets *this to an unsigned integer value.
// Returns any remainder.
template <typename UINT> RT_API_ATTRS UINT SetTo(UINT n) {
static_assert(
std::is_same_v<UINT, common::uint128_t> || std::is_unsigned_v<UINT>);
SetToZero();
while (n != 0) {
auto q{n / 10u};
if (n != q * 10) {
break;
}
++exponent_;
n = q;
}
if constexpr (sizeof n < sizeof(Digit)) {
if (n != 0) {
digit_[digits_++] = n;
}
return 0;
} else {
while (n != 0 && digits_ < digitLimit_) {
auto q{n / radix};
digit_[digits_++] = static_cast<Digit>(n - q * radix);
n = q;
}
return n;
}
}
RT_API_ATTRS int RemoveLeastOrderZeroDigits() {
int remove{0};
if (digits_ > 0 && digit_[0] == 0) {
while (remove < digits_ && digit_[remove] == 0) {
++remove;
}
if (remove >= digits_) {
digits_ = 0;
} else if (remove > 0) {
#if defined __GNUC__ && __GNUC__ < 8
// (&& j + remove < maxDigits) was added to avoid GCC < 8 build failure
// on -Werror=array-bounds. This can be removed if -Werror is disable.
for (int j{0}; j + remove < digits_ && (j + remove < maxDigits); ++j) {
#else
for (int j{0}; j + remove < digits_; ++j) {
#endif
digit_[j] = digit_[j + remove];
}
digits_ -= remove;
}
}
return remove;
}
RT_API_ATTRS void RemoveLeadingZeroDigits() {
while (digits_ > 0 && digit_[digits_ - 1] == 0) {
--digits_;
}
}
RT_API_ATTRS void Normalize() {
RemoveLeadingZeroDigits();
exponent_ += RemoveLeastOrderZeroDigits() * log10Radix;
}
// This limited divisibility test only works for even divisors of the radix,
// which is fine since it's only ever used with 2 and 5.
template <int N> RT_API_ATTRS bool IsDivisibleBy() const {
static_assert(N > 1 && radix % N == 0, "bad modulus");
return digits_ == 0 || (digit_[0] % N) == 0;
}
template <unsigned DIVISOR> RT_API_ATTRS int DivideBy() {
Digit remainder{0};
for (int j{digits_ - 1}; j >= 0; --j) {
Digit q{digit_[j] / DIVISOR};
Digit nrem{digit_[j] - DIVISOR * q};
digit_[j] = q + (radix / DIVISOR) * remainder;
remainder = nrem;
}
return remainder;
}
RT_API_ATTRS void DivideByPowerOfTwo(int twoPow) { // twoPow <= log10Radix
Digit remainder{0};
auto mask{(Digit{1} << twoPow) - 1};
auto coeff{radix >> twoPow};
for (int j{digits_ - 1}; j >= 0; --j) {
auto nrem{digit_[j] & mask};
digit_[j] = (digit_[j] >> twoPow) + coeff * remainder;
remainder = nrem;
}
}
// Returns true on overflow
RT_API_ATTRS bool DivideByPowerOfTwoInPlace(int twoPow) {
if (digits_ > 0) {
while (twoPow > 0) {
int chunk{twoPow > log10Radix ? log10Radix : twoPow};
if ((digit_[0] & ((Digit{1} << chunk) - 1)) == 0) {
DivideByPowerOfTwo(chunk);
twoPow -= chunk;
continue;
}
twoPow -= chunk;
if (digit_[digits_ - 1] >> chunk != 0) {
if (digits_ == digitLimit_) {
return true; // overflow
}
digit_[digits_++] = 0;
}
auto remainder{digit_[digits_ - 1]};
exponent_ -= log10Radix;
auto coeff{radix >> chunk}; // precise; radix is (5*2)**log10Radix
auto mask{(Digit{1} << chunk) - 1};
for (int j{digits_ - 1}; j >= 1; --j) {
digit_[j] = (digit_[j - 1] >> chunk) + coeff * remainder;
remainder = digit_[j - 1] & mask;
}
digit_[0] = coeff * remainder;
}
}
return false; // no overflow
}
RT_API_ATTRS int AddCarry(int position = 0, int carry = 1) {
for (; position < digits_; ++position) {
Digit v{digit_[position] + carry};
if (v < radix) {
digit_[position] = v;
return 0;
}
digit_[position] = v - radix;
carry = 1;
}
if (digits_ < digitLimit_) {
digit_[digits_++] = carry;
return 0;
}
Normalize();
if (digits_ < digitLimit_) {
digit_[digits_++] = carry;
return 0;
}
return carry;
}
RT_API_ATTRS void Decrement() {
for (int j{0}; digit_[j]-- == 0; ++j) {
digit_[j] = radix - 1;
}
}
template <int N> RT_API_ATTRS int MultiplyByHelper(int carry = 0) {
for (int j{0}; j < digits_; ++j) {
auto v{N * digit_[j] + carry};
carry = v / radix;
digit_[j] = v - carry * radix; // i.e., v % radix
}
return carry;
}
template <int N> RT_API_ATTRS int MultiplyBy(int carry = 0) {
if (int newCarry{MultiplyByHelper<N>(carry)}) {
return AddCarry(digits_, newCarry);
} else {
return 0;
}
}
template <int N> RT_API_ATTRS int MultiplyWithoutNormalization() {
if (int carry{MultiplyByHelper<N>(0)}) {
if (digits_ < digitLimit_) {
digit_[digits_++] = carry;
return 0;
} else {
return carry;
}
} else {
return 0;
}
}
RT_API_ATTRS void LoseLeastSignificantDigit(); // with rounding
RT_API_ATTRS void PushCarry(int carry) {
if (digits_ == maxDigits && RemoveLeastOrderZeroDigits() == 0) {
LoseLeastSignificantDigit();
digit_[digits_ - 1] += carry;
} else {
digit_[digits_++] = carry;
}
}
// Adds another number and then divides by two.
// Assumes same exponent and sign.
// Returns true when the result has effectively been rounded down.
RT_API_ATTRS bool Mean(const BigRadixFloatingPointNumber &);
// Parses a floating-point number; leaves the pointer reference
// argument pointing at the next character after what was recognized.
// The "end" argument can be left null if the caller is sure that the
// string is properly terminated with an addressable character that
// can't be in a valid floating-point character.
RT_API_ATTRS bool ParseNumber(const char *&, bool &inexact, const char *end);
using Raw = typename Real::RawType;
constexpr RT_API_ATTRS Raw SignBit() const {
return Raw{isNegative_} << (Real::bits - 1);
}
constexpr RT_API_ATTRS Raw Infinity() const {
Raw result{static_cast<Raw>(Real::maxExponent)};
result <<= Real::significandBits;
result |= SignBit();
if constexpr (Real::bits == 80) { // x87
result |= Raw{1} << 63;
}
return result;
}
constexpr RT_API_ATTRS Raw NaN(bool isQuiet = true) {
Raw result{Real::maxExponent};
result <<= Real::significandBits;
result |= SignBit();
if constexpr (Real::bits == 80) { // x87
result |= Raw{isQuiet ? 3u : 2u} << 62;
} else {
Raw quiet{isQuiet ? Raw{2} : Raw{1}};
quiet <<= Real::significandBits - 2;
result |= quiet;
}
return result;
}
constexpr RT_API_ATTRS Raw HUGE() const {
Raw result{static_cast<Raw>(Real::maxExponent)};
result <<= Real::significandBits;
result |= SignBit();
return result - 1; // decrement exponent, set all significand bits
}
Digit digit_[maxDigits]; // in little-endian order: digit_[0] is LSD
int digits_{0}; // # of elements in digit_[] array; zero when zero
int digitLimit_{maxDigits}; // precision clamp
int exponent_{0}; // signed power of ten
bool isNegative_{false};
enum FortranRounding rounding_ { RoundNearest };
};
} // namespace Fortran::decimal
#endif
Codebase Browser
llvm