/usr/include/boost/multiprecision
Edit: /usr/include/boost/multiprecision/cpp_dec_float.hpp (141498B)
///////////////////////////////////////////////////////////////////////////////
// Copyright Christopher Kormanyos 2002 - 2013.
// Copyright 2011 -2013 John Maddock. Distributed under the Boost
// Software License, Version 1.0. (See accompanying file
// LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
//
// This work is based on an earlier work:
// "Algorithm 910: A Portable C++ Multiple-Precision System for Special-Function Calculations",
// in ACM TOMS, {VOL 37, ISSUE 4, (February 2011)} (C) ACM, 2011. http://doi.acm.org/10.1145/1916461.1916469
//
// Note that there are no "noexcept" specifications on the functions in this file: there are too many
// calls to lexical_cast (and similar) to easily analyse the code for correctness. So until compilers
// can detect noexcept misuse at compile time, the only realistic option is to simply not use it here.
//
#ifndef BOOST_MP_CPP_DEC_FLOAT_BACKEND_HPP
#define BOOST_MP_CPP_DEC_FLOAT_BACKEND_HPP
#include
#include
#include
#ifndef BOOST_NO_CXX11_HDR_ARRAY
#include
#else
#include
#endif
#include
#include
#include
#include
#include
#include
//
// Headers required for Boost.Math integration:
//
#include
//
// Some includes we need from Boost.Math, since we rely on that library to provide these functions:
//
#include
#include
#include
#include
#include
#include
#ifdef BOOST_MSVC
#pragma warning(push)
#pragma warning(disable : 6326) // comparison of two constants
#endif
namespace boost {
namespace multiprecision {
namespace backends {
template
class cpp_dec_float;
} // namespace backends
template
struct number_category > : public mpl::int_
{};
namespace backends {
template
class cpp_dec_float
{
private:
static const boost::int32_t cpp_dec_float_digits10_setting = Digits10;
// We need at least 16-bits in the exponent type to do anything sensible:
BOOST_STATIC_ASSERT_MSG(boost::is_signed::value, "ExponentType must be a signed built in integer type.");
BOOST_STATIC_ASSERT_MSG(sizeof(ExponentType) > 1, "ExponentType is too small.");
public:
typedef mpl::list signed_types;
typedef mpl::list unsigned_types;
typedef mpl::list float_types;
typedef ExponentType exponent_type;
static const boost::int32_t cpp_dec_float_radix = 10L;
static const boost::int32_t cpp_dec_float_digits10_limit_lo = 9L;
static const boost::int32_t cpp_dec_float_digits10_limit_hi = boost::integer_traits::const_max - 100;
static const boost::int32_t cpp_dec_float_digits10 = ((cpp_dec_float_digits10_setting < cpp_dec_float_digits10_limit_lo) ? cpp_dec_float_digits10_limit_lo : ((cpp_dec_float_digits10_setting > cpp_dec_float_digits10_limit_hi) ? cpp_dec_float_digits10_limit_hi : cpp_dec_float_digits10_setting));
static const ExponentType cpp_dec_float_max_exp10 = (static_cast(1) << (std::numeric_limits::digits - 5));
static const ExponentType cpp_dec_float_min_exp10 = -cpp_dec_float_max_exp10;
static const ExponentType cpp_dec_float_max_exp = cpp_dec_float_max_exp10;
static const ExponentType cpp_dec_float_min_exp = cpp_dec_float_min_exp10;
BOOST_STATIC_ASSERT((cpp_dec_float::cpp_dec_float_max_exp10 == -cpp_dec_float::cpp_dec_float_min_exp10));
private:
static const boost::int32_t cpp_dec_float_elem_digits10 = 8L;
static const boost::int32_t cpp_dec_float_elem_mask = 100000000L;
BOOST_STATIC_ASSERT(0 == cpp_dec_float_max_exp10 % cpp_dec_float_elem_digits10);
// There are three guard limbs.
// 1) The first limb has 'play' from 1...8 decimal digits.
// 2) The last limb also has 'play' from 1...8 decimal digits.
// 3) One limb can get lost when justifying after multiply,
// as only half of the triangle is multiplied and a carry
// from below is missing.
static const boost::int32_t cpp_dec_float_elem_number_request = static_cast((cpp_dec_float_digits10 / cpp_dec_float_elem_digits10) + (((cpp_dec_float_digits10 % cpp_dec_float_elem_digits10) != 0) ? 1 : 0));
// The number of elements needed (with a minimum of two) plus three added guard limbs.
static const boost::int32_t cpp_dec_float_elem_number = static_cast(((cpp_dec_float_elem_number_request < 2L) ? 2L : cpp_dec_float_elem_number_request) + 3L);
public:
static const boost::int32_t cpp_dec_float_total_digits10 = static_cast(cpp_dec_float_elem_number * cpp_dec_float_elem_digits10);
private:
typedef enum enum_fpclass_type
{
cpp_dec_float_finite,
cpp_dec_float_inf,
cpp_dec_float_NaN
} fpclass_type;
#ifndef BOOST_NO_CXX11_HDR_ARRAY
typedef typename mpl::if_,
std::array,
detail::dynamic_array >::type array_type;
#else
typedef typename mpl::if_,
boost::array,
detail::dynamic_array >::type array_type;
#endif
array_type data;
ExponentType exp;
bool neg;
fpclass_type fpclass;
boost::int32_t prec_elem;
//
// Special values constructor:
//
cpp_dec_float(fpclass_type c) : data(),
exp(static_cast(0)),
neg(false),
fpclass(c),
prec_elem(cpp_dec_float_elem_number) {}
//
// Static data initializer:
//
struct initializer
{
initializer()
{
cpp_dec_float::nan();
cpp_dec_float::inf();
(cpp_dec_float::min)();
(cpp_dec_float::max)();
cpp_dec_float::zero();
cpp_dec_float::one();
cpp_dec_float::two();
cpp_dec_float::half();
cpp_dec_float::double_min();
cpp_dec_float::double_max();
//cpp_dec_float::long_double_max();
//cpp_dec_float::long_double_min();
cpp_dec_float::long_long_max();
cpp_dec_float::long_long_min();
cpp_dec_float::ulong_long_max();
cpp_dec_float::eps();
cpp_dec_float::pow2(0);
}
void do_nothing() {}
};
static initializer init;
struct long_double_initializer
{
long_double_initializer()
{
cpp_dec_float::long_double_max();
cpp_dec_float::long_double_min();
}
void do_nothing() {}
};
static long_double_initializer linit;
public:
// Constructors
cpp_dec_float() BOOST_MP_NOEXCEPT_IF(noexcept(array_type())) : data(),
exp(static_cast(0)),
neg(false),
fpclass(cpp_dec_float_finite),
prec_elem(cpp_dec_float_elem_number) {}
cpp_dec_float(const char* s) : data(),
exp(static_cast(0)),
neg(false),
fpclass(cpp_dec_float_finite),
prec_elem(cpp_dec_float_elem_number)
{
*this = s;
}
template
cpp_dec_float(I i, typename enable_if >::type* = 0) : data(),
exp(static_cast(0)),
neg(false),
fpclass(cpp_dec_float_finite),
prec_elem(cpp_dec_float_elem_number)
{
from_unsigned_long_long(i);
}
template
cpp_dec_float(I i, typename enable_if >::type* = 0) : data(),
exp(static_cast(0)),
neg(false),
fpclass(cpp_dec_float_finite),
prec_elem(cpp_dec_float_elem_number)
{
if (i < 0)
{
from_unsigned_long_long(boost::multiprecision::detail::unsigned_abs(i));
negate();
}
else
from_unsigned_long_long(i);
}
cpp_dec_float(const cpp_dec_float& f) BOOST_MP_NOEXCEPT_IF(noexcept(array_type(std::declval()))) : data(f.data),
exp(f.exp),
neg(f.neg),
fpclass(f.fpclass),
prec_elem(f.prec_elem) {}
template
cpp_dec_float(const cpp_dec_float& f, typename enable_if_c::type* = 0) : data(),
exp(f.exp),
neg(f.neg),
fpclass(static_cast(static_cast(f.fpclass))),
prec_elem(cpp_dec_float_elem_number)
{
std::copy(f.data.begin(), f.data.begin() + f.prec_elem, data.begin());
}
template
explicit cpp_dec_float(const cpp_dec_float& f, typename disable_if_c::type* = 0) : data(),
exp(f.exp),
neg(f.neg),
fpclass(static_cast(static_cast(f.fpclass))),
prec_elem(cpp_dec_float_elem_number)
{
// TODO: this doesn't round!
std::copy(f.data.begin(), f.data.begin() + prec_elem, data.begin());
}
template
cpp_dec_float(const F val, typename enable_if_c::value
#ifdef BOOST_HAS_FLOAT128
&& !boost::is_same::value
#endif
>::type* = 0) : data(),
exp(static_cast(0)),
neg(false),
fpclass(cpp_dec_float_finite),
prec_elem(cpp_dec_float_elem_number)
{
*this = val;
}
cpp_dec_float(const double mantissa, const ExponentType exponent);
std::size_t hash() const
{
std::size_t result = 0;
for (int i = 0; i < prec_elem; ++i)
boost::hash_combine(result, data[i]);
boost::hash_combine(result, exp);
boost::hash_combine(result, neg);
boost::hash_combine(result, fpclass);
return result;
}
// Specific special values.
static const cpp_dec_float& nan()
{
static const cpp_dec_float val(cpp_dec_float_NaN);
init.do_nothing();
return val;
}
static const cpp_dec_float& inf()
{
static const cpp_dec_float val(cpp_dec_float_inf);
init.do_nothing();
return val;
}
static const cpp_dec_float&(max)()
{
init.do_nothing();
static cpp_dec_float val_max = std::string("1.0e" + boost::multiprecision::detail::itos(cpp_dec_float_max_exp10)).c_str();
return val_max;
}
static const cpp_dec_float&(min)()
{
init.do_nothing();
static cpp_dec_float val_min = std::string("1.0e" + boost::multiprecision::detail::itos(cpp_dec_float_min_exp10)).c_str();
return val_min;
}
static const cpp_dec_float& zero()
{
init.do_nothing();
static cpp_dec_float val(static_cast(0u));
return val;
}
static const cpp_dec_float& one()
{
init.do_nothing();
static cpp_dec_float val(static_cast(1u));
return val;
}
static const cpp_dec_float& two()
{
init.do_nothing();
static cpp_dec_float val(static_cast(2u));
return val;
}
static const cpp_dec_float& half()
{
init.do_nothing();
static cpp_dec_float val(0.5L);
return val;
}
static const cpp_dec_float& double_min()
{
init.do_nothing();
static cpp_dec_float val((std::numeric_limits::min)());
return val;
}
static const cpp_dec_float& double_max()
{
init.do_nothing();
static cpp_dec_float val((std::numeric_limits::max)());
return val;
}
static const cpp_dec_float& long_double_min()
{
linit.do_nothing();
#ifdef BOOST_MATH_NO_LONG_DOUBLE_MATH_FUNCTIONS
static cpp_dec_float val(static_cast((std::numeric_limits::min)()));
#else
static cpp_dec_float val((std::numeric_limits::min)());
#endif
return val;
}
static const cpp_dec_float& long_double_max()
{
linit.do_nothing();
#ifdef BOOST_MATH_NO_LONG_DOUBLE_MATH_FUNCTIONS
static cpp_dec_float val(static_cast((std::numeric_limits::max)()));
#else
static cpp_dec_float val((std::numeric_limits::max)());
#endif
return val;
}
static const cpp_dec_float& long_long_max()
{
init.do_nothing();
static cpp_dec_float val((std::numeric_limits::max)());
return val;
}
static const cpp_dec_float& long_long_min()
{
init.do_nothing();
static cpp_dec_float val((std::numeric_limits::min)());
return val;
}
static const cpp_dec_float& ulong_long_max()
{
init.do_nothing();
static cpp_dec_float val((std::numeric_limits::max)());
return val;
}
static const cpp_dec_float& eps()
{
init.do_nothing();
static cpp_dec_float val(1.0, 1 - static_cast(cpp_dec_float_digits10));
return val;
}
// Basic operations.
cpp_dec_float& operator=(const cpp_dec_float& v) BOOST_MP_NOEXCEPT_IF(noexcept(std::declval() = std::declval()))
{
data = v.data;
exp = v.exp;
neg = v.neg;
fpclass = v.fpclass;
prec_elem = v.prec_elem;
return *this;
}
template
cpp_dec_float& operator=(const cpp_dec_float& f)
{
exp = f.exp;
neg = f.neg;
fpclass = static_cast(static_cast(f.fpclass));
unsigned elems = (std::min)(f.prec_elem, cpp_dec_float_elem_number);
std::copy(f.data.begin(), f.data.begin() + elems, data.begin());
std::fill(data.begin() + elems, data.end(), 0);
prec_elem = cpp_dec_float_elem_number;
return *this;
}
cpp_dec_float& operator=(boost::long_long_type v)
{
if (v < 0)
{
from_unsigned_long_long(1u - boost::ulong_long_type(v + 1)); // Avoid undefined behaviour in negation of minimum value for long long
negate();
}
else
from_unsigned_long_long(v);
return *this;
}
cpp_dec_float& operator=(boost::ulong_long_type v)
{
from_unsigned_long_long(v);
return *this;
}
template
typename boost::enable_if_c::value, cpp_dec_float&>::type operator=(Float v);
cpp_dec_float& operator=(const char* v)
{
rd_string(v);
return *this;
}
cpp_dec_float& operator+=(const cpp_dec_float& v);
cpp_dec_float& operator-=(const cpp_dec_float& v);
cpp_dec_float& operator*=(const cpp_dec_float& v);
cpp_dec_float& operator/=(const cpp_dec_float& v);
cpp_dec_float& add_unsigned_long_long(const boost::ulong_long_type n)
{
cpp_dec_float t;
t.from_unsigned_long_long(n);
return *this += t;
}
cpp_dec_float& sub_unsigned_long_long(const boost::ulong_long_type n)
{
cpp_dec_float t;
t.from_unsigned_long_long(n);
return *this -= t;
}
cpp_dec_float& mul_unsigned_long_long(const boost::ulong_long_type n);
cpp_dec_float& div_unsigned_long_long(const boost::ulong_long_type n);
// Elementary primitives.
cpp_dec_float& calculate_inv();
cpp_dec_float& calculate_sqrt();
void negate()
{
if (!iszero())
neg = !neg;
}
// Comparison functions
bool isnan BOOST_PREVENT_MACRO_SUBSTITUTION() const { return (fpclass == cpp_dec_float_NaN); }
bool isinf BOOST_PREVENT_MACRO_SUBSTITUTION() const { return (fpclass == cpp_dec_float_inf); }
bool isfinite BOOST_PREVENT_MACRO_SUBSTITUTION() const { return (fpclass == cpp_dec_float_finite); }
bool iszero() const
{
return ((fpclass == cpp_dec_float_finite) && (data[0u] == 0u));
}
bool isone() const;
bool isint() const;
bool isneg() const { return neg; }
// Operators pre-increment and pre-decrement
cpp_dec_float& operator++()
{
return *this += one();
}
cpp_dec_float& operator--()
{
return *this -= one();
}
std::string str(boost::intmax_t digits, std::ios_base::fmtflags f) const;
int compare(const cpp_dec_float& v) const;
template
int compare(const V& v) const
{
cpp_dec_float t;
t = v;
return compare(t);
}
void swap(cpp_dec_float& v)
{
data.swap(v.data);
std::swap(exp, v.exp);
std::swap(neg, v.neg);
std::swap(fpclass, v.fpclass);
std::swap(prec_elem, v.prec_elem);
}
double extract_double() const;
long double extract_long_double() const;
boost::long_long_type extract_signed_long_long() const;
boost::ulong_long_type extract_unsigned_long_long() const;
void extract_parts(double& mantissa, ExponentType& exponent) const;
cpp_dec_float extract_integer_part() const;
void precision(const boost::int32_t prec_digits)
{
if (prec_digits >= cpp_dec_float_total_digits10)
{
prec_elem = cpp_dec_float_elem_number;
}
else
{
const boost::int32_t elems = static_cast(static_cast((prec_digits + (cpp_dec_float_elem_digits10 / 2)) / cpp_dec_float_elem_digits10) + static_cast(((prec_digits % cpp_dec_float_elem_digits10) != 0) ? 1 : 0));
prec_elem = (std::min)(cpp_dec_float_elem_number, (std::max)(elems, static_cast(2)));
}
}
static cpp_dec_float pow2(boost::long_long_type i);
ExponentType order() const
{
const bool bo_order_is_zero = ((!(isfinite)()) || (data[0] == static_cast(0u)));
//
// Binary search to find the order of the leading term:
//
ExponentType prefix = 0;
if (data[0] >= 100000UL)
{
if (data[0] >= 10000000UL)
{
if (data[0] >= 100000000UL)
{
if (data[0] >= 1000000000UL)
prefix = 9;
else
prefix = 8;
}
else
prefix = 7;
}
else
{
if (data[0] >= 1000000UL)
prefix = 6;
else
prefix = 5;
}
}
else
{
if (data[0] >= 1000UL)
{
if (data[0] >= 10000UL)
prefix = 4;
else
prefix = 3;
}
else
{
if (data[0] >= 100)
prefix = 2;
else if (data[0] >= 10)
prefix = 1;
}
}
return (bo_order_is_zero ? static_cast(0) : static_cast(exp + prefix));
}
template
void serialize(Archive& ar, const unsigned int /*version*/)
{
for (unsigned i = 0; i < data.size(); ++i)
ar& boost::make_nvp("digit", data[i]);
ar& boost::make_nvp("exponent", exp);
ar& boost::make_nvp("sign", neg);
ar& boost::make_nvp("class-type", fpclass);
ar& boost::make_nvp("precision", prec_elem);
}
private:
static bool data_elem_is_non_zero_predicate(const boost::uint32_t& d) { return (d != static_cast(0u)); }
static bool data_elem_is_non_nine_predicate(const boost::uint32_t& d) { return (d != static_cast(cpp_dec_float::cpp_dec_float_elem_mask - 1)); }
static bool char_is_nonzero_predicate(const char& c) { return (c != static_cast('0')); }
void from_unsigned_long_long(const boost::ulong_long_type u);
int cmp_data(const array_type& vd) const;
static boost::uint32_t mul_loop_uv(boost::uint32_t* const u, const boost::uint32_t* const v, const boost::int32_t p);
static boost::uint32_t mul_loop_n(boost::uint32_t* const u, boost::uint32_t n, const boost::int32_t p);
static boost::uint32_t div_loop_n(boost::uint32_t* const u, boost::uint32_t n, const boost::int32_t p);
bool rd_string(const char* const s);
template
friend class cpp_dec_float;
};
template
typename cpp_dec_float::initializer cpp_dec_float::init;
template
typename cpp_dec_float::long_double_initializer cpp_dec_float::linit;
template
const boost::int32_t cpp_dec_float::cpp_dec_float_radix;
template
const boost::int32_t cpp_dec_float::cpp_dec_float_digits10_setting;
template
const boost::int32_t cpp_dec_float::cpp_dec_float_digits10_limit_lo;
template
const boost::int32_t cpp_dec_float::cpp_dec_float_digits10_limit_hi;
template
const boost::int32_t cpp_dec_float::cpp_dec_float_digits10;
template
const ExponentType cpp_dec_float::cpp_dec_float_max_exp;
template
const ExponentType cpp_dec_float::cpp_dec_float_min_exp;
template
const ExponentType cpp_dec_float::cpp_dec_float_max_exp10;
template
const ExponentType cpp_dec_float::cpp_dec_float_min_exp10;
template
const boost::int32_t cpp_dec_float::cpp_dec_float_elem_digits10;
template
const boost::int32_t cpp_dec_float::cpp_dec_float_elem_number_request;
template
const boost::int32_t cpp_dec_float::cpp_dec_float_elem_number;
template
const boost::int32_t cpp_dec_float::cpp_dec_float_elem_mask;
template
cpp_dec_float& cpp_dec_float::operator+=(const cpp_dec_float& v)
{
if ((isnan)())
{
return *this;
}
if ((isinf)())
{
if ((v.isinf)() && (isneg() != v.isneg()))
{
*this = nan();
}
return *this;
}
if (iszero())
{
return operator=(v);
}
if ((v.isnan)() || (v.isinf)())
{
*this = v;
return *this;
}
// Get the offset for the add/sub operation.
static const ExponentType max_delta_exp = static_cast((cpp_dec_float_elem_number - 1) * cpp_dec_float_elem_digits10);
const ExponentType ofs_exp = static_cast(exp - v.exp);
// Check if the operation is out of range, requiring special handling.
if (v.iszero() || (ofs_exp > max_delta_exp))
{
// Result is *this unchanged since v is negligible compared to *this.
return *this;
}
else if (ofs_exp < -max_delta_exp)
{
// Result is *this = v since *this is negligible compared to v.
return operator=(v);
}
// Do the add/sub operation.
typename array_type::iterator p_u = data.begin();
typename array_type::const_iterator p_v = v.data.begin();
bool b_copy = false;
const boost::int32_t ofs = static_cast(static_cast(ofs_exp) / cpp_dec_float_elem_digits10);
array_type n_data;
if (neg == v.neg)
{
// Add v to *this, where the data array of either *this or v
// might have to be treated with a positive, negative or zero offset.
// The result is stored in *this. The data are added one element
// at a time, each element with carry.
if (ofs >= static_cast(0))
{
std::copy(v.data.begin(), v.data.end() - static_cast(ofs), n_data.begin() + static_cast(ofs));
std::fill(n_data.begin(), n_data.begin() + static_cast(ofs), static_cast(0u));
p_v = n_data.begin();
}
else
{
std::copy(data.begin(), data.end() - static_cast(-ofs), n_data.begin() + static_cast(-ofs));
std::fill(n_data.begin(), n_data.begin() + static_cast(-ofs), static_cast(0u));
p_u = n_data.begin();
b_copy = true;
}
// Addition algorithm
boost::uint32_t carry = static_cast(0u);
for (boost::int32_t j = static_cast(cpp_dec_float_elem_number - static_cast(1)); j >= static_cast(0); j--)
{
boost::uint32_t t = static_cast(static_cast(p_u[j] + p_v[j]) + carry);
carry = t / static_cast(cpp_dec_float_elem_mask);
p_u[j] = static_cast(t - static_cast(carry * static_cast(cpp_dec_float_elem_mask)));
}
if (b_copy)
{
data = n_data;
exp = v.exp;
}
// There needs to be a carry into the element -1 of the array data
if (carry != static_cast(0u))
{
std::copy_backward(data.begin(), data.end() - static_cast(1u), data.end());
data[0] = carry;
exp += static_cast(cpp_dec_float_elem_digits10);
}
}
else
{
// Subtract v from *this, where the data array of either *this or v
// might have to be treated with a positive, negative or zero offset.
if ((ofs > static_cast(0)) || ((ofs == static_cast(0)) && (cmp_data(v.data) > static_cast(0))))
{
// In this case, |u| > |v| and ofs is positive.
// Copy the data of v, shifted down to a lower value
// into the data array m_n. Set the operand pointer p_v
// to point to the copied, shifted data m_n.
std::copy(v.data.begin(), v.data.end() - static_cast(ofs), n_data.begin() + static_cast(ofs));
std::fill(n_data.begin(), n_data.begin() + static_cast(ofs), static_cast(0u));
p_v = n_data.begin();
}
else
{
if (ofs != static_cast(0))
{
// In this case, |u| < |v| and ofs is negative.
// Shift the data of u down to a lower value.
std::copy_backward(data.begin(), data.end() - static_cast(-ofs), data.end());
std::fill(data.begin(), data.begin() + static_cast(-ofs), static_cast(0u));
}
// Copy the data of v into the data array n_data.
// Set the u-pointer p_u to point to m_n and the
// operand pointer p_v to point to the shifted
// data m_data.
n_data = v.data;
p_u = n_data.begin();
p_v = data.begin();
b_copy = true;
}
boost::int32_t j;
// Subtraction algorithm
boost::int32_t borrow = static_cast(0);
for (j = static_cast(cpp_dec_float_elem_number - static_cast(1)); j >= static_cast(0); j--)
{
boost::int32_t t = static_cast(static_cast(static_cast(p_u[j]) - static_cast(p_v[j])) - borrow);
// Underflow? Borrow?
if (t < static_cast(0))
{
// Yes, underflow and borrow
t += static_cast(cpp_dec_float_elem_mask);
borrow = static_cast(1);
}
else
{
borrow = static_cast(0);
}
p_u[j] = static_cast(static_cast(t) % static_cast(cpp_dec_float_elem_mask));
}
if (b_copy)
{
data = n_data;
exp = v.exp;
neg = v.neg;
}
// Is it necessary to justify the data?
const typename array_type::const_iterator first_nonzero_elem = std::find_if(data.begin(), data.end(), data_elem_is_non_zero_predicate);
if (first_nonzero_elem != data.begin())
{
if (first_nonzero_elem == data.end())
{
// This result of the subtraction is exactly zero.
// Reset the sign and the exponent.
neg = false;
exp = static_cast(0);
}
else
{
// Justify the data
const std::size_t sj = static_cast(std::distance(data.begin(), first_nonzero_elem));
std::copy(data.begin() + static_cast(sj), data.end(), data.begin());
std::fill(data.end() - sj, data.end(), static_cast(0u));
exp -= static_cast(sj * static_cast(cpp_dec_float_elem_digits10));
}
}
}
// Handle underflow.
if (iszero())
return (*this = zero());
// Check for potential overflow.
const bool b_result_might_overflow = (exp >= static_cast(cpp_dec_float_max_exp10));
// Handle overflow.
if (b_result_might_overflow)
{
const bool b_result_is_neg = neg;
neg = false;
if (compare((cpp_dec_float::max)()) > 0)
*this = inf();
neg = b_result_is_neg;
}
return *this;
}
template
cpp_dec_float& cpp_dec_float::operator-=(const cpp_dec_float& v)
{
// Use *this - v = -(-*this + v).
negate();
*this += v;
negate();
return *this;
}
template
cpp_dec_float& cpp_dec_float::operator*=(const cpp_dec_float& v)
{
// Evaluate the sign of the result.
const bool b_result_is_neg = (neg != v.neg);
// Artificially set the sign of the result to be positive.
neg = false;
// Handle special cases like zero, inf and NaN.
const bool b_u_is_inf = (isinf)();
const bool b_v_is_inf = (v.isinf)();
const bool b_u_is_zero = iszero();
const bool b_v_is_zero = v.iszero();
if (((isnan)() || (v.isnan)()) || (b_u_is_inf && b_v_is_zero) || (b_v_is_inf && b_u_is_zero))
{
*this = nan();
return *this;
}
if (b_u_is_inf || b_v_is_inf)
{
*this = inf();
if (b_result_is_neg)
negate();
return *this;
}
if (b_u_is_zero || b_v_is_zero)
{
return *this = zero();
}
// Check for potential overflow or underflow.
const bool b_result_might_overflow = ((exp + v.exp) >= static_cast(cpp_dec_float_max_exp10));
const bool b_result_might_underflow = ((exp + v.exp) <= static_cast(cpp_dec_float_min_exp10));
// Set the exponent of the result.
exp += v.exp;
const boost::int32_t prec_mul = (std::min)(prec_elem, v.prec_elem);
const boost::uint32_t carry = mul_loop_uv(data.data(), v.data.data(), prec_mul);
// Handle a potential carry.
if (carry != static_cast(0u))
{
exp += cpp_dec_float_elem_digits10;
// Shift the result of the multiplication one element to the right...
std::copy_backward(data.begin(),
data.begin() + static_cast(prec_elem - static_cast(1)),
data.begin() + static_cast(prec_elem));
// ... And insert the carry.
data.front() = carry;
}
// Handle overflow.
if (b_result_might_overflow && (compare((cpp_dec_float::max)()) > 0))
{
*this = inf();
}
// Handle underflow.
if (b_result_might_underflow && (compare((cpp_dec_float::min)()) < 0))
{
*this = zero();
return *this;
}
// Set the sign of the result.
neg = b_result_is_neg;
return *this;
}
template
cpp_dec_float& cpp_dec_float::operator/=(const cpp_dec_float& v)
{
if (iszero())
{
if ((v.isnan)())
{
return *this = v;
}
else if (v.iszero())
{
return *this = nan();
}
}
const bool u_and_v_are_finite_and_identical = ((isfinite)() && (fpclass == v.fpclass) && (exp == v.exp) && (cmp_data(v.data) == static_cast(0)));
if (u_and_v_are_finite_and_identical)
{
if (neg != v.neg)
{
*this = one();
negate();
}
else
*this = one();
return *this;
}
else
{
cpp_dec_float t(v);
t.calculate_inv();
return operator*=(t);
}
}
template
cpp_dec_float& cpp_dec_float::mul_unsigned_long_long(const boost::ulong_long_type n)
{
// Multiply *this with a constant boost::ulong_long_type.
// Evaluate the sign of the result.
const bool b_neg = neg;
// Artificially set the sign of the result to be positive.
neg = false;
// Handle special cases like zero, inf and NaN.
const bool b_u_is_inf = (isinf)();
const bool b_n_is_zero = (n == static_cast(0));
if ((isnan)() || (b_u_is_inf && b_n_is_zero))
{
return (*this = nan());
}
if (b_u_is_inf)
{
*this = inf();
if (b_neg)
negate();
return *this;
}
if (iszero() || b_n_is_zero)
{
// Multiplication by zero.
return *this = zero();
}
if (n >= static_cast(cpp_dec_float_elem_mask))
{
neg = b_neg;
cpp_dec_float t;
t = n;
return operator*=(t);
}
if (n == static_cast(1u))
{
neg = b_neg;
return *this;
}
// Set up the multiplication loop.
const boost::uint32_t nn = static_cast(n);
const boost::uint32_t carry = mul_loop_n(data.data(), nn, prec_elem);
// Handle the carry and adjust the exponent.
if (carry != static_cast(0u))
{
exp += static_cast(cpp_dec_float_elem_digits10);
// Shift the result of the multiplication one element to the right.
std::copy_backward(data.begin(),
data.begin() + static_cast(prec_elem - static_cast(1)),
data.begin() + static_cast(prec_elem));
data.front() = static_cast(carry);
}
// Check for potential overflow.
const bool b_result_might_overflow = (exp >= cpp_dec_float_max_exp10);
// Handle overflow.
if (b_result_might_overflow && (compare((cpp_dec_float::max)()) > 0))
{
*this = inf();
}
// Set the sign.
neg = b_neg;
return *this;
}
template
cpp_dec_float& cpp_dec_float::div_unsigned_long_long(const boost::ulong_long_type n)
{
// Divide *this by a constant boost::ulong_long_type.
// Evaluate the sign of the result.
const bool b_neg = neg;
// Artificially set the sign of the result to be positive.
neg = false;
// Handle special cases like zero, inf and NaN.
if ((isnan)())
{
return *this;
}
if ((isinf)())
{
*this = inf();
if (b_neg)
negate();
return *this;
}
if (n == static_cast