More tests and syntax

This commit is contained in:
ckormanyos
2025-07-19 22:30:57 +02:00
parent 618132d3bd
commit e4aada5489
6 changed files with 201 additions and 67 deletions
+47 -28
View File
@@ -168,6 +168,9 @@ namespace backends {
// The type of the floating-point constituents should adhere to IEEE754.
// This class has been tested with floats having single-precision (4 byte),
// double-precision (8 byte) and quad precision (16 byte, such as GCC's __float128).
// Although the constituent parts (a0 and a1) satisfy |a1| <= (1 / 2) * ulp(a0),
// the composite type does not adhere to these strict error bounds. Its error
// bounds are larger.
template <typename FloatingPointType>
class cpp_double_fp_backend
@@ -184,7 +187,7 @@ class cpp_double_fp_backend
|| ((cpp_df_qf_detail::ccmath::numeric_limits<float_type>::digits == 53) && std::numeric_limits<float_type>::is_specialized && std::numeric_limits<float_type>::is_iec559)
|| ((cpp_df_qf_detail::ccmath::numeric_limits<float_type>::digits == 64) && std::numeric_limits<float_type>::is_specialized && std::numeric_limits<float_type>::is_iec559)
|| (cpp_df_qf_detail::ccmath::numeric_limits<float_type>::digits == 113)
}, "Error: float_type does not fulfil the backend requirements of cpp_double_fp_backend"
}, "Error: float_type does not fulfill the backend requirements of cpp_double_fp_backend"
);
using rep_type = cpp_df_qf_detail::pair<float_type, float_type>;
@@ -317,7 +320,7 @@ class cpp_double_fp_backend
&& ((static_cast<int>(sizeof(SignedIntegralType) * 8) - 1) > cpp_df_qf_detail::ccmath::numeric_limits<float_type>::digits))>::type const* = nullptr>
constexpr cpp_double_fp_backend(SignedIntegralType n)
{
const bool is_neg { n < SignedIntegralType { INT8_C(0) } };
const bool is_neg { (n < SignedIntegralType { INT8_C(0) }) };
using local_unsigned_integral_type = typename boost::multiprecision::detail::make_unsigned<SignedIntegralType>::type;
@@ -436,7 +439,14 @@ class cpp_double_fp_backend
constexpr auto iszero_unchecked() const noexcept -> bool { return (data.first == float_type { 0.0F }); }
constexpr auto is_one() const noexcept -> bool { return ((data.second == float_type { 0.0F }) && (data.first == float_type { 1.0F })); }
constexpr auto is_one() const noexcept -> bool
{
return
(
(data.second == float_type { 0.0F })
&& (data.first == float_type { 1.0F })
);
}
constexpr auto negate() -> void
{
@@ -714,7 +724,8 @@ class cpp_double_fp_backend
// The division algorithm has been taken from Victor Shoup,
// package WinNTL-5_3_2. It might originally be related to the
// K. Briggs work. The algorithm has been significantly simplified
// while still atempting to retain proper rounding corrections.
// while still attempting to retain proper rounding corrections.
// Checks for overflow and underflow have been added.
const float_type C { data.first / v.data.first };
@@ -762,11 +773,11 @@ class cpp_double_fp_backend
u = cpp_df_qf_detail::ccmath::unsafe::fma(hc, hv, -U);
const float_type tv { v.data.first - hv };
u = cpp_df_qf_detail::ccmath::unsafe::fma(hc, tv, u);
{
const float_type tv { v.data.first - hv };
u = cpp_df_qf_detail::ccmath::unsafe::fma(hc, tv, u);
const float_type tc { C - hc };
u = cpp_df_qf_detail::ccmath::unsafe::fma(tc, hv, u);
@@ -1028,7 +1039,8 @@ class cpp_double_fp_backend
// The multiplication algorithm has been taken from Victor Shoup,
// package WinNTL-5_3_2. It might originally be related to the
// K. Briggs work. The algorithm has been significantly simplified
// while still atempting to retain proper rounding corrections.
// while still attempting to retain proper rounding corrections.
// Checks for overflow and underflow have been added.
float_type C { cpp_df_qf_detail::split_maker<float_type>::value * data.first };
@@ -1081,19 +1093,26 @@ class cpp_double_fp_backend
hv = c - float_type { c - v.data.first };
}
const float_type tv { v.data.first - hv };
{
const float_type tv { v.data.first - hv };
float_type t1 { cpp_df_qf_detail::ccmath::unsafe::fma(hu, hv, -C) };
const float_type
t1
{
cpp_df_qf_detail::ccmath::unsafe::fma
(
hu,
tv,
cpp_df_qf_detail::ccmath::unsafe::fma(hu, hv, -C)
)
};
t1 = cpp_df_qf_detail::ccmath::unsafe::fma(hu, tv, t1);
const float_type tu { data.first - hu };
const float_type tu { data.first - hu };
t1 = cpp_df_qf_detail::ccmath::unsafe::fma(tu, hv, t1);
c = cpp_df_qf_detail::ccmath::unsafe::fma(tu, tv, t1)
+ (data.first * v.data.second)
+ (data.second * v.data.first);
c = cpp_df_qf_detail::ccmath::unsafe::fma(tu, tv, cpp_df_qf_detail::ccmath::unsafe::fma(tu, hv, t1))
+ (data.first * v.data.second)
+ (data.second * v.data.first);
}
// Perform even more simplifications compared to Victor Shoup.
data.first = C + c;
@@ -1116,13 +1135,13 @@ class cpp_double_fp_backend
template <typename FloatingPointType>
auto cpp_double_fp_backend<FloatingPointType>::rd_string(const char* pstr) -> bool
{
using local_double_fp_type = cpp_double_fp_backend<FloatingPointType>;
cpp_bin_float_read_write_type f_bin { pstr };
const auto fpc = fpclassify(f_bin);
const int fpc { fpclassify(f_bin) };
const auto is_definitely_nan = (fpc == FP_NAN);
const bool is_definitely_nan { (fpc == FP_NAN) };
using local_double_fp_type = cpp_double_fp_backend<FloatingPointType>;
if (is_definitely_nan)
{
@@ -1562,7 +1581,7 @@ constexpr auto eval_pow(cpp_double_fp_backend<FloatingPointType>& result, const
template <typename FloatingPointType,
typename IntegralType>
constexpr auto eval_pow(cpp_double_fp_backend<FloatingPointType>& result, const cpp_double_fp_backend<FloatingPointType>& x, IntegralType p) -> typename ::std::enable_if<boost::multiprecision::detail::is_integral<IntegralType>::value, void>::type
constexpr auto eval_pow(cpp_double_fp_backend<FloatingPointType>& result, const cpp_double_fp_backend<FloatingPointType>& x, IntegralType p) -> typename ::std::enable_if<::boost::multiprecision::detail::is_integral<IntegralType>::value, void>::type
{
const int fpc { eval_fpclassify(x) };
@@ -2503,11 +2522,11 @@ auto hash_value(const cpp_double_fp_backend<FloatingPointType>& a) -> ::std::siz
using backends::cpp_double_fp_backend;
using cpp_double_float = number<cpp_double_fp_backend<float>, boost::multiprecision::et_off>;
using cpp_double_double = number<cpp_double_fp_backend<double>, boost::multiprecision::et_off>;
using cpp_double_long_double = number<cpp_double_fp_backend<long double>, boost::multiprecision::et_off>;
using cpp_double_float = number<cpp_double_fp_backend<float>, ::boost::multiprecision::et_off>;
using cpp_double_double = number<cpp_double_fp_backend<double>, ::boost::multiprecision::et_off>;
using cpp_double_long_double = number<cpp_double_fp_backend<long double>, ::boost::multiprecision::et_off>;
#ifdef BOOST_MP_CPP_DOUBLE_FP_HAS_FLOAT128
using cpp_double_float128 = number<cpp_double_fp_backend<::boost::float128_type>, boost::multiprecision::et_off>;
using cpp_double_float128 = number<cpp_double_fp_backend<::boost::float128_type>, ::boost::multiprecision::et_off>;
#endif
} } // namespace boost::multiprecision
@@ -672,7 +672,7 @@ inline void eval_pow(T& result, const T& x, const T& a)
}
else
{
result = std::numeric_limits<number<T, et_on> >::infinity().backend();
result = std::numeric_limits<number<T, et_on> >::infinity().backend(); // LCOV_EXCL_LINE There is no NaN for this number type.
}
}
else BOOST_IF_CONSTEXPR (std::numeric_limits<number<T, et_on> >::has_quiet_NaN)
@@ -170,7 +170,7 @@ void eval_sin(T& result, const T& x)
using si_type = typename boost::multiprecision::detail::canonical<std::int32_t, T>::type ;
using ui_type = typename boost::multiprecision::detail::canonical<std::uint32_t, T>::type;
using fp_type = typename std::tuple_element<0, typename T::float_types>::type ;
using fp_type = typename std::tuple_element<0, typename T::float_types>::type;
switch (eval_fpclassify(x))
{
+95 -22
View File
@@ -17,14 +17,11 @@
#include <test.hpp>
#include <boost/detail/lightweight_test.hpp>
#include <array>
#include <ctime>
#include <random>
#include <boost/math/constants/constants.hpp>
#if !defined(TEST_MPF_50) && !defined(TEST_MPF) && !defined(TEST_BACKEND) && !defined(TEST_CPP_DEC_FLOAT) && !defined(TEST_MPFR) && !defined(TEST_MPFR_50) && !defined(TEST_MPFI_50) && !defined(TEST_FLOAT128) && !defined(TEST_CPP_BIN_FLOAT) && !defined(TEST_CPP_DOUBLE_FLOAT)
#define TEST_MPF_50
//# define TEST_MPF
//#define TEST_MPF
#define TEST_BACKEND
#define TEST_CPP_DEC_FLOAT
#define TEST_MPFI_50
@@ -68,6 +65,13 @@
#include <boost/multiprecision/cpp_double_fp.hpp>
#endif
#include <array>
#include <ctime>
#include <random>
template<typename FloatType> auto my_zero() -> FloatType&;
template<typename FloatType> auto my_one() -> FloatType&;
template <class T>
void test()
{
@@ -77,7 +81,9 @@ void test()
BOOST_IF_CONSTEXPR (std::numeric_limits<T>::is_specialized && (std::numeric_limits<T>::max_exponent10 > 4000))
{
static const std::array<const char*, 51u> data =
using table_data_array_type = std::array<const char*, 51u>;
static const table_data_array_type table_data =
{{
"1.00000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000",
"9.47747587596218770242116751705184563668845029215054154915126374673142219159548534317576897266130328412495991561490384353e76",
@@ -132,26 +138,31 @@ void test()
"6.83336127500041943234365059231968669406267422759442985746460610830503287734479988530512309065240678799786759250323660701e3848",
}};
T pi = static_cast<T>("3.141592653589793238462643383279502884197169399375105820974944592307816406286208998628034825342117067982148086513282306647093844609550582231725359408128481117450284102701938521105559644622948954930381964428810975665933446128475648233786783165271201909145648566923460348610454326648213393607260249141273724587006606315588174881520920962829254091715364367892590360011330530548820466521384146951941511609");
const T pi = static_cast<T>("3.141592653589793238462643383279502884197169399375105820974944592307816406286208998628034825342117067982148086513282306647093844609550582231725359408128481117450284102701938521105559644622948954930381964428810975665933446128475648233786783165271201909145648566923460348610454326648213393607260249141273724587006606315588174881520920962829254091715364367892590360011330530548820466521384146951941511609");
for (unsigned k = 0; k < data.size(); k++)
for (unsigned k = 0u; k < static_cast<unsigned>(std::tuple_size<table_data_array_type>::value); ++k)
{
T val = exp(sqrt((pi * (100 * k)) * (100 * k)));
T e = relative_error(val, T(data[k]));
T val = exp(sqrt((pi * (100u * k)) * (100u * k)));
T e = relative_error(val, T(table_data[k]));
unsigned err = e.template convert_to<unsigned>();
if (err > max_err)
{
max_err = err;
}
val = exp(-sqrt((pi * (100 * k)) * (100 * k)));
e = relative_error(val, T(1 / T(data[k])));
val = exp(-sqrt((pi * (100u * k)) * (100u * k)));
e = relative_error(val, T(1 / T(table_data[k])));
err = e.template convert_to<unsigned>();
if (err > max_err)
{
max_err = err;
}
}
std::cout << "Max error was: " << max_err << std::endl;
std::cout << "Max error table_data was: " << max_err << std::endl;
#if defined(BOOST_INTEL) && defined(TEST_FLOAT128)
BOOST_TEST(max_err < 40000);
#elif defined(TEST_CPP_BIN_FLOAT)
@@ -163,7 +174,10 @@ void test()
using std::ldexp;
static const std::array<std::array<T, 2>, 12> exact_data =
using exact_data_pair_type = std::array<T, 2>;
using exact_data_array_type = std::array<exact_data_pair_type, 12>;
static const exact_data_array_type exact_data =
{{
{{ldexp(1.0, -50), static_cast<T>("1.000000000000000888178419700125626769357944978675730736309709508287711059579809241499236575743374705946980126761002249532899810120773330275600867188192232364653350140592709328919189811425580736404494785")}},
{{ldexp(1.0, -20), static_cast<T>("1.00000095367477115374544678824955687428365188553281789775169686343569285229334215539516690752571791280462887427635269562079697496032436580742164524046357050365736415701568566320292733574692386949329504")}},
@@ -175,42 +189,83 @@ void test()
{{10.5, static_cast<T>("36315.5026742466377389120269013166179689315579671275857607480190550842856628099187749764427758174866310742771977376827511779240563292392413797025035252944088473059613508585937585174897482367249804575829")}},
{{25, static_cast<T>("7.20048993373858725241613514661261579152235338133952787362213864472320593107782569745000325654258093194727871848859163683530731259683905547347259346394203424381962487428148616548084393082148632398178103e10")}},
{{31.25, static_cast<T>("3.72994612957188849046766396046821396700589012875701157893019118883826370993674081486706667149871508642909416337810227575115729423667289322729813729072376711271966343652718317681492250090327804073781915e13")}},
// N[Log[39614081257132168796771975168], 201]
{{static_cast<T>("65.8489821531948043946370515385267739671725127642242491414646009018723940871209979825570160646597753164901406969542151447001244223970224029181037214053322163834191192286952863430791686692697089797567183"), static_cast<T>("39614081257132168796771975168.0")}},
// N[Log[19807040628566084398385987584], 201]
{{static_cast<T>("65.1558349726348590852198194170685973990970126298639938873439208923790004651513032669511527376633566289481392159336444589664389021612642723610710506536971404214883916578669149078888616306458173062855949"), static_cast<T>("19807040628566084398385987584.0")}},
}};
max_err = 0;
for (unsigned k = 0; k < exact_data.size(); k++)
for (unsigned k = 0u; k < static_cast<unsigned>(std::tuple_size<exact_data_array_type>::value); ++k)
{
T val = exp(exact_data[k][0]);
T e = relative_error(val, exact_data[k][1]);
T val = exp(exact_data[k][0]);
T e = relative_error(val, exact_data[k][1]);
unsigned err = e.template convert_to<unsigned>();
if (err > max_err)
{
max_err = err;
}
val = exp(-exact_data[k][0]);
e = relative_error(val, T(1 / exact_data[k][1]));
err = e.template convert_to<unsigned>();
if (err > max_err)
{
max_err = err;
}
}
std::cout << "Max error was: " << max_err << std::endl;
std::cout << "Max error exact_data was: " << max_err << std::endl;
BOOST_TEST(max_err < 60);
BOOST_TEST(exp(T(0)) == 1);
if (!boost::multiprecision::is_interval_number<T>::value)
BOOST_IF_CONSTEXPR (!boost::multiprecision::is_interval_number<T>::value)
{
std::mt19937_64 gen { };
gen.seed(static_cast<typename std::mt19937_64::result_type>(std::clock()));
std::uniform_real_distribution<float>
dist
(
static_cast<float>(1.01L),
static_cast<float>(1.04L)
);
for (int index = 0; index < 8; ++index)
{
static_cast<void>(index);
const T val_zero { ::my_zero<T>() * dist(gen) };
const T exp_zero = exp(val_zero);
BOOST_CHECK(exp_zero == ::my_one<T>());
}
BOOST_IF_CONSTEXPR (std::numeric_limits<T>::is_specialized)
{
for (int index = 0; index < 8; ++index)
{
static_cast<void>(index);
const T val_one(static_cast<int>(::my_one<T>() * dist(gen)));
const T exp_one = exp(val_one);
BOOST_CHECK_CLOSE_FRACTION(exp_one, boost::math::constants::e<T>(), std::numeric_limits<T>::epsilon() * 4);
}
}
T bug_case = -1.05 * log((std::numeric_limits<T>::max)());
for (unsigned i = 0U; bug_case > -20 / std::numeric_limits<T>::epsilon(); ++i, bug_case *= 1.05)
{
static_cast<void>(i);
if (std::numeric_limits<T>::has_infinity)
BOOST_IF_CONSTEXPR (std::numeric_limits<T>::has_infinity)
{
BOOST_CHECK_EQUAL(exp(bug_case), 0);
}
@@ -290,3 +345,21 @@ int main()
#endif
return boost::report_errors();
}
template<typename FloatType> auto my_zero() -> FloatType&
{
using float_type = FloatType;
static float_type my_val_zero(0);
return my_val_zero;
}
template<typename FloatType> auto my_one() -> FloatType&
{
using float_type = FloatType;
static float_type my_val_one(1);
return my_val_one;
}
+1
View File
@@ -25,6 +25,7 @@
//# define TEST_MPF
#define TEST_BACKEND
#define TEST_CPP_DEC_FLOAT
#define TEST_FLOAT128
#define TEST_MPFR_50
#define TEST_MPFI_50
#define TEST_CPP_BIN_FLOAT
+56 -15
View File
@@ -625,6 +625,47 @@ namespace local
return result_is_ok;
}
template<class FloatType>
bool test_edges_trig()
{
using float_type = FloatType;
auto dis =
std::uniform_real_distribution<float>
{
static_cast<float>(1.01F),
static_cast<float>(1.04F)
};
bool result_is_ok { true };
using std::isinf;
using std::isnan;
using std::signbit;
BOOST_IF_CONSTEXPR(std::numeric_limits<float_type>::has_quiet_NaN)
{
std::mt19937_64 gen { UINT64_C(0xF00DCAFEDEADBEEF) };
for(auto index = static_cast<unsigned>(UINT8_C(0)); index < static_cast<unsigned>(UINT8_C(8)); ++index)
{
static_cast<void>(index);
float_type flt_nan = std::numeric_limits<float_type>::quiet_NaN() * dis(gen);
const float_type atan_nan = atan(flt_nan);
const bool result_atan_nan_is_ok { (boost::multiprecision::isnan)(atan_nan) && ((boost::multiprecision::fpclassify)(atan_nan) == FP_NAN) };
BOOST_TEST(result_atan_nan_is_ok);
result_is_ok = (result_atan_nan_is_ok && result_is_ok);
}
}
return result_is_ok;
}
template<class OtherFloatType>
auto test_convert_and_back_caller(const float epsilon_factor = 0.0F) -> bool
{
@@ -790,27 +831,27 @@ namespace local
auto main() -> int
{
using bin_float_backend_type = boost::multiprecision::cpp_bin_float<50>;
using dec_float_backend_type = boost::multiprecision::cpp_dec_float<50>;
using bin_float_type = boost::multiprecision::number<bin_float_backend_type, boost::multiprecision::et_off>;
using dec_float_type = boost::multiprecision::number<dec_float_backend_type, boost::multiprecision::et_off>;
{
using float_backend_type = boost::multiprecision::cpp_bin_float<50>;
std::cout << "Testing type: " << typeid(bin_float_type).name() << std::endl;
using float_type = boost::multiprecision::number<float_backend_type, boost::multiprecision::et_off>;
std::cout << "Testing type: " << typeid(float_type).name() << std::endl;
static_cast<void>(local::test_edges<float_type>());
static_cast<void>(local::test_edges<bin_float_type>());
static_cast<void>(local::test_edges_trig<bin_float_type>());
}
{
using float_backend_type = boost::multiprecision::cpp_dec_float<50>;
std::cout << "Testing type: " << typeid(dec_float_type).name() << std::endl;
using float_type = boost::multiprecision::number<float_backend_type, boost::multiprecision::et_off>;
std::cout << "Testing type: " << typeid(float_type).name() << std::endl;
static_cast<void>(local::test_edges<float_type>());
local::test_cpp_dec_float_rd_ovf_unf<float_type>();
local::test_convert_and_back<double, float_type>(0.0F);
local::test_cpp_dec_float_frexp_edge<float_type>();
static_cast<void>(local::test_edges<dec_float_type>());
static_cast<void>(local::test_edges_trig<dec_float_type>());
local::test_cpp_dec_float_rd_ovf_unf<dec_float_type>();
local::test_convert_and_back<double, dec_float_type>(0.0F);
local::test_cpp_dec_float_frexp_edge<dec_float_type>();
}
{