boost/math
2
0
Fork 0
mirror of https://github.com/boostorg/math.git synced 2026-01-19 04:22:09 +00:00
Code Issues Projects Releases Wiki Activity

Merge branch 'boostorg:develop' into develop

Browse Source
This commit is contained in:
Maksym Zhelyeznyakov
2025-08-16 15:00:26 +02:00
committed by GitHub
parent b5429d4dd5 7535ee105b
commit 0b319a3cfe
21 changed files with 554 additions and 53 deletions
Download Patch File Download Diff File Expand all files Collapse all files

1
doc/constants/constants.qbk
View File

@@ -340,6 +340,7 @@ This section lists the mathematical constants, their use(s) (and sometimes ratio
[[pi_minus_three] [[pi]-3] [0.141593] [] ]
[[four_minus_pi] [4 -[pi]] [0.858407] [] ]
[[pi_pow_e] [[pi][super e]] [22.4591] [] ]
[[log_pi] [ln([pi])] [1.14473] [] ]
[[pi_sqr] [[pi][super 2]] [9.86960] [] ]
[[pi_sqr_div_six] [[pi][super 2]/6] [1.64493] [] ]

9
doc/math.qbk
View File

@@ -339,6 +339,10 @@ and use the function's name as the link text.]
[def __hermite [link math_toolkit.sf_poly.hermite hermite]]
[def __cardinal_b_splines [link math_toolkit.sf_poly.cardinal_b_splines cardinal_b_splines]]
[/logistic functions]
[def __logit [link math_toolkit.logistic.logit logit]]
[def __logistic_sigmoid [link math_toolkit.logistic.logistic_sigmoid logistic_sigmoid]]
[/Misc]
[def __expint [link math_toolkit.expint.expint_i expint]]
[def __spherical_harmonic [link math_toolkit.sf_poly.sph_harm spherical_harmonic]]
@@ -666,6 +670,11 @@ and as a CD ISBN 0-9504833-2-X 978-0-9504833-2-0, Classification 519.2-dc22.
[include sf/jacobi.qbk]
[endsect] [/section:sf_poly Polynomials]
[section:logistic Logistic Functions]
[include sf/logit.qbk]
[include sf/logistic_sigmoid.qbk]
[endsect] [/section:logistic Logistic Functions]
[section:bessel Bessel Functions]
[include sf/bessel_introduction.qbk]
[include sf/bessel_jy.qbk]

30
doc/sf/logistic_sigmoid.qbk Normal file
View File

@@ -0,0 +1,30 @@
[/
Copyright Matt Borland 2025
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).
]
[section:logistic_sigmoid logistic_sigmoid]
[h4 Synopsis]
#include <boost/math/special_functions/logistic_sigmoid.hpp>
namespace boost { namespace math {
template <typename RealType, typename Policy>
RealType logistic_sigmoid(RealType x, const Policy&);
template <typename RealType>
RealType logistic_sigmoid(RealType x)
}} // namespaces
[h4 Description]
Returns the [@https://en.wikipedia.org/wiki/Logistic_function logistic sigmoid function]
This function is broadly useful, and is used to compute the CDF of the logistic distribution.
It is also sometimes referred to as expit as it is the inverse of the logit function.
[endsect]

30
doc/sf/logit.qbk Normal file
View File

@@ -0,0 +1,30 @@
[/
Copyright Matt Borland 2025
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).
]
[section:logit logit]
[h4 Synopsis]
#include <boost/math/special_functions/logit.hpp>
namespace boost { namespace math {
template <typename RealType, typename Policy>
RealType logit(RealType x, const Policy&);
template <typename RealType>
RealType logit(RealType x)
}} // namespaces
[h4 Description]
Returns the [@https://en.wikipedia.org/wiki/Logit logit function]
This function is broadly useful, and is used to compute the Quantile of the logistic distribution.
The inverse of this function is the logistic_sigmoid.
[endsect]

8
include/boost/math/constants/calculate_constants.hpp
View File

@@ -95,6 +95,14 @@ inline T constant_one_div_two_pi<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_S
return 1 / two_pi<T, policies::policy<policies::digits2<N> > >();
}
template <class T>
template<int N>
inline T constant_log_pi<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC((std::integral_constant<int, N>)))
{
BOOST_MATH_STD_USING
return log(pi<T, policies::policy<policies::digits2<N> > >());
}
template <class T>
template<int N>
inline T constant_root_pi<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC((std::integral_constant<int, N>)))

3
include/boost/math/constants/constants.hpp
View File

@@ -280,6 +280,7 @@ namespace constants {
BOOST_DEFINE_MATH_CONSTANT(four_thirds_pi, 4.188790204786390984616857844372670512e+00, "4.18879020478639098461685784437267051226289253250014109463325945641042187504827866483737976712282275730953078202e+00")
BOOST_DEFINE_MATH_CONSTANT(one_div_two_pi, 1.591549430918953357688837633725143620e-01, "1.59154943091895335768883763372514362034459645740456448747667344058896797634226535090113802766253085956072842727e-01")
BOOST_DEFINE_MATH_CONSTANT(one_div_root_two_pi, 3.989422804014326779399460599343818684e-01, "3.98942280401432677939946059934381868475858631164934657665925829670657925899301838501252333907306936430302558863e-01")
BOOST_DEFINE_MATH_CONSTANT(log_pi, 1.144729885849400174143427351353058711e+00, "1.14472988584940017414342735135305871164729481291531157151362307147213776988482607978362327027548970770200981223e+00")
BOOST_DEFINE_MATH_CONSTANT(root_pi, 1.772453850905516027298167483341145182e+00, "1.77245385090551602729816748334114518279754945612238712821380778985291128459103218137495065673854466541622682362e+00")
BOOST_DEFINE_MATH_CONSTANT(root_half_pi, 1.253314137315500251207882642405522626e+00, "1.25331413731550025120788264240552262650349337030496915831496178817114682730392098747329791918902863305800498633e+00")
BOOST_DEFINE_MATH_CONSTANT(root_two_pi, 2.506628274631000502415765284811045253e+00, "2.50662827463100050241576528481104525300698674060993831662992357634229365460784197494659583837805726611600997267e+00")
@@ -357,5 +358,3 @@ BOOST_MATH_GPU_ENABLED inline constexpr T tau() { return two_pi<T>(); }
#endif
#endif // BOOST_MATH_CONSTANTS_CONSTANTS_INCLUDED

48
include/boost/math/distributions/logistic.hpp
View File

@@ -19,6 +19,8 @@
#include <boost/math/constants/constants.hpp>
#include <boost/math/policies/policy.hpp>
#include <boost/math/policies/error_handling.hpp>
#include <boost/math/special_functions/logit.hpp>
#include <boost/math/special_functions/logistic_sigmoid.hpp>
namespace boost { namespace math {
@@ -144,13 +146,10 @@ namespace boost { namespace math {
{
return result;
}
BOOST_MATH_STD_USING
RealType power = (location - x) / scale;
if(power > tools::log_max_value<RealType>())
return 0;
if(power < -tools::log_max_value<RealType>())
return 1;
return 1 / (1 + exp(power));
using promoted_real_type = typename policies::evaluation<RealType, Policy>::type;
promoted_real_type power = (static_cast<promoted_real_type>(x) - static_cast<promoted_real_type>(location)) / static_cast<promoted_real_type>(scale);
return logistic_sigmoid(power, policies::make_forwarding_policy_t<Policy>());
}
template <class RealType, class Policy>
@@ -159,7 +158,7 @@ namespace boost { namespace math {
RealType scale = dist.scale();
RealType location = dist.location();
RealType result = 0; // of checks.
constexpr auto function = "boost::math::cdf(const logistic_distribution<%1%>&, %1%)";
constexpr auto function = "boost::math::logcdf(const logistic_distribution<%1%>&, %1%)";
if(false == detail::check_scale(function, scale, &result, Policy()))
{
return result;
@@ -199,7 +198,6 @@ namespace boost { namespace math {
template <class RealType, class Policy>
BOOST_MATH_GPU_ENABLED inline RealType quantile(const logistic_distribution<RealType, Policy>& dist, const RealType& p)
{
BOOST_MATH_STD_USING
RealType location = dist.location();
RealType scale = dist.scale();
@@ -221,21 +219,13 @@ namespace boost { namespace math {
{
return policies::raise_overflow_error<RealType>(function,"probability argument is 1, must be >0 and <1",Policy());
}
//Expressions to try
//return location+scale*log(p/(1-p));
//return location+scale*log1p((2*p-1)/(1-p));
//return location - scale*log( (1-p)/p);
//return location - scale*log1p((1-2*p)/p);
//return -scale*log(1/p-1) + location;
return location - scale * log((1 - p) / p);
return location + scale * logit(p, Policy());
} // RealType quantile(const logistic_distribution<RealType, Policy>& dist, const RealType& p)
template <class RealType, class Policy>
BOOST_MATH_GPU_ENABLED inline RealType cdf(const complemented2_type<logistic_distribution<RealType, Policy>, RealType>& c)
{
BOOST_MATH_STD_USING
RealType location = c.dist.location();
RealType scale = c.dist.scale();
RealType x = c.param;
@@ -259,12 +249,10 @@ namespace boost { namespace math {
{
return result;
}
RealType power = (x - location) / scale;
if(power > tools::log_max_value<RealType>())
return 0;
if(power < -tools::log_max_value<RealType>())
return 1;
return 1 / (1 + exp(power));
using promoted_real_type = typename policies::evaluation<RealType, Policy>::type;
promoted_real_type power = (static_cast<promoted_real_type>(location) - static_cast<promoted_real_type>(x)) / static_cast<promoted_real_type>(scale);
return logistic_sigmoid(power, policies::make_forwarding_policy_t<Policy>());
}
template <class RealType, class Policy>
@@ -274,7 +262,7 @@ namespace boost { namespace math {
RealType location = c.dist.location();
RealType scale = c.dist.scale();
RealType x = c.param;
constexpr auto function = "boost::math::cdf(const complement(logistic_distribution<%1%>&), %1%)";
constexpr auto function = "boost::math::logcdf(const complement(logistic_distribution<%1%>&), %1%)";
RealType result = 0;
if(false == detail::check_scale(function, scale, &result, Policy()))
@@ -306,7 +294,6 @@ namespace boost { namespace math {
template <class RealType, class Policy>
BOOST_MATH_GPU_ENABLED inline RealType quantile(const complemented2_type<logistic_distribution<RealType, Policy>, RealType>& c)
{
BOOST_MATH_STD_USING
RealType scale = c.dist.scale();
RealType location = c.dist.location();
constexpr auto function = "boost::math::quantile(const complement(logistic_distribution<%1%>&), %1%)";
@@ -328,15 +315,8 @@ namespace boost { namespace math {
{
return policies::raise_overflow_error<RealType>(function,"probability argument is 0, but must be >0 and <1",Policy());
}
//Expressions to try
//return location+scale*log((1-q)/q);
return location + scale * log((1 - q) / q);
//return location-scale*log(q/(1-q));
//return location-scale*log1p((2*q-1)/(1-q));
//return location+scale*log(1/q-1);
//return location+scale*log1p(1/q-2);
return location - scale * logit(q, Policy());
}
template <class RealType, class Policy>

19
include/boost/math/policies/policy.hpp
View File

@@ -1003,6 +1003,25 @@ struct is_noexcept_error_policy
&& (t8::value != throw_on_error) && (t8::value != user_error));
};
// Generate a forwarding policy to stop further promotion from occurring
// For example if a special function for float promotes to double, we don't want the next
// function in the call chain to then promote to long double
template <typename Policy>
struct make_forwarding_policy
{
using type =
typename policies::normalise<
Policy,
policies::promote_float<false>,
policies::promote_double<false>,
policies::discrete_quantile<>,
policies::assert_undefined<>
>::type;
};
template <typename Policy>
using make_forwarding_policy_t = typename make_forwarding_policy<Policy>::type;
}}} // namespaces
#endif // BOOST_MATH_POLICY_HPP

2
include/boost/math/special_functions/detail/bessel_jn.hpp
View File

@@ -84,7 +84,7 @@ BOOST_MATH_GPU_ENABLED T bessel_jn(int n, T x, const Policy& pol)
current = value;
}
}
else if((x < 1) || (n > x * x / 4) || (x < 5))
else if((x < 5) || (n > x * x / 4))
{
return factor * bessel_j_small_z_series(T(n), x, pol);
}

5
include/boost/math/special_functions/detail/bessel_jy.hpp
View File

@@ -327,7 +327,10 @@ namespace boost { namespace math {
// x is positive until reflection
W = T(2) / (x * pi<T>()); // Wronskian
T Yv_scale = 1;
if(((kind & need_y) == 0) && ((x < 1) || (v > x * x / 4) || (x < 5)))
const bool kind_does_not_need_y { ((kind & need_y) == 0) };
if(kind_does_not_need_y && ((x < 5) || (v > x * x / 4)))
{
//
// This series will actually converge rapidly for all small

46
include/boost/math/special_functions/logistic_sigmoid.hpp Normal file
View File

@@ -0,0 +1,46 @@
// Copyright Matt Borland 2025.
// Use, modification and distribution are subject to 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)
#ifndef BOOST_MATH_SF_EXPIT_HPP
#define BOOST_MATH_SF_EXPIT_HPP
#include <boost/math/policies/policy.hpp>
#include <boost/math/tools/precision.hpp>
#include <cmath>
namespace boost {
namespace math {
template <typename RealType, typename Policy>
RealType logistic_sigmoid(RealType x, const Policy&)
{
BOOST_MATH_STD_USING
using promoted_real_type = typename policies::evaluation<RealType, Policy>::type;
if(-x >= tools::log_max_value<RealType>())
{
return static_cast<RealType>(0);
}
if(-x <= -tools::log_max_value<RealType>())
{
return static_cast<RealType>(1);
}
const auto res {static_cast<RealType>(1 / (1 + exp(static_cast<promoted_real_type>(-x))))};
return res;
}
template <typename RealType>
RealType logistic_sigmoid(RealType x)
{
return logistic_sigmoid(x, policies::policy<>());
}
} // namespace math
} // namespace boost
#endif // BOOST_MATH_SF_EXPIT_HPP

56
include/boost/math/special_functions/logit.hpp Normal file
View File

@@ -0,0 +1,56 @@
// Copyright Matt Borland 2025.
// Use, modification and distribution are subject to 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)
#ifndef BOOST_MATH_SF_LOGIT_HPP
#define BOOST_MATH_SF_LOGIT_HPP
#include <boost/math/tools/config.hpp>
#include <boost/math/policies/policy.hpp>
#include <boost/math/policies/error_handling.hpp>
#include <cmath>
#include <cfenv>
namespace boost {
namespace math {
template <typename RealType, typename Policy>
RealType logit(RealType p, const Policy&)
{
BOOST_MATH_STD_USING
using std::atanh;
using promoted_real_type = typename policies::evaluation<RealType, Policy>::type;
if (p < tools::min_value<RealType>())
{
return -policies::raise_overflow_error<RealType>("logit", "sub-normals will overflow ln(x/(1-x))", Policy());
}
static const RealType crossover {RealType{1}/4};
const auto promoted_p {static_cast<promoted_real_type>(p)};
RealType result {};
if (p > crossover)
{
result = static_cast<RealType>(2 * atanh(2 * promoted_p - 1));
}
else
{
result = static_cast<RealType>(log(promoted_p / (1 - promoted_p)));
}
return result;
}
template <typename RealType>
RealType logit(RealType p)
{
return logit(p, policies::policy<>());
}
} // namespace math
} // namespace boost
#endif // BOOST_MATH_SF_LOGIT_HPP

2
test/Jamfile.v2
View File

@@ -590,6 +590,8 @@ test-suite special_fun :
[ run test_sinc.cpp /boost/test//boost_unit_test_framework pch_light ]
[ run test_fibonacci.cpp /boost/test//boost_unit_test_framework ]
[ run test_prime.cpp /boost/test//boost_unit_test_framework ]
[ run test_logistic_sigmoid.cpp ]
[ run test_logit.cpp ]
;
test-suite distribution_tests :

33
test/git_issue_1294.cpp Normal file
View File

@@ -0,0 +1,33 @@
// Copyright 2025 Matt Borland
// Use, modification and distribution are subject to 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)
//
// See: https://github.com/boostorg/math/issues/1294
#include <boost/math/distributions/logistic.hpp>
#include "math_unit_test.hpp"
int main()
{
using namespace boost::math::policies;
using boost::math::logistic_distribution;
typedef policy<
promote_float<true>,
promote_double<true>
> with_promotion;
constexpr double p = 2049.0/4096;
constexpr double ref = 9.76562577610225755e-04;
logistic_distribution<double, with_promotion> dist_promote;
const double x = quantile(dist_promote, p);
// Previously we had: 9.76562577610170027e-04
// Which is an ULP distance of 256
CHECK_ULP_CLOSE(x, ref, 1);
return boost::math::test::report_errors();
}

2
test/math_unit_test.hpp
View File

@@ -58,7 +58,7 @@ bool check_mollified_close(Real expected, Real computed, Real tol, std::string c
}
using std::max;
using std::abs;
Real denom = (max)(abs(expected), Real(1));
Real denom = (max)(Real(abs(expected)), Real(1));
Real mollified_relative_error = abs(expected - computed)/denom;
if (mollified_relative_error > tol)
{

4
test/test_bessel_j.cpp
View File

@@ -345,7 +345,3 @@ BOOST_AUTO_TEST_CASE( test_main )
"to pass.</note>" << std::endl;
#endif
}

95
test/test_bessel_j.hpp
View File

@@ -1,4 +1,5 @@
// (C) Copyright John Maddock 2007.
// (C) Copyright John Maddock 2007 - 2025.
// (C) Copyright Christopher Kormanyos 2025.
// Use, modification and distribution are subject to 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)
@@ -131,6 +132,8 @@ void do_test_sph_bessel_j(const T& data, const char* type_name, const char* test
template <class T>
void test_bessel(T, const char* name)
{
using table_entry_type = typename table_type<T>::type;
//
// The actual test data is rather verbose, so it's in a separate file
//
@@ -138,7 +141,7 @@ void test_bessel(T, const char* name)
// three items, input value a, input value b and erf(a, b):
//
// function values calculated on http://functions.wolfram.com/
static const std::array<std::array<typename table_type<T>::type, 3>, 8> j0_data = {{
static const std::array<std::array<table_entry_type, 3>, 8> j0_data = {{
{ { SC_(0.0), SC_(0.0), SC_(1.0) } },
{ { SC_(0.0), SC_(1.0), SC_(0.7651976865579665514497175261026632209093) } },
{ { SC_(0.0), SC_(-2.0), SC_(0.2238907791412356680518274546499486258252) } },
@@ -148,7 +151,7 @@ void test_bessel(T, const char* name)
{{ SC_(0.0), SC_(1e-10), SC_(0.999999999999999999997500000000000000000) }},
{{ SC_(0.0), SC_(-1e+01), SC_(-0.2459357644513483351977608624853287538296) }},
}};
static const std::array<std::array<typename table_type<T>::type, 3>, 6> j0_tricky = {{
static const std::array<std::array<table_entry_type, 3>, 6> j0_tricky = {{
// Big numbers make the accuracy of std::sin the limiting factor:
{ { SC_(0.0), SC_(1e+03), SC_(0.02478668615242017456133073111569370878617) } },
{ { SC_(0.0), SC_(1e+05), SC_(-0.001719201116235972192570601477073201747532) } },
@@ -159,7 +162,7 @@ void test_bessel(T, const char* name)
{ { SC_(0.0), SC_(11.791534423828125) /*T(12364320.0) / (1024 * 1024)*/, SC_(-3.53017140778223781420794006033810387155048392363051866610931e-9) } }
}};
static const std::array<std::array<typename table_type<T>::type, 3>, 8> j1_data = {{
static const std::array<std::array<table_entry_type, 3>, 8> j1_data = {{
{ { SC_(1.0), SC_(0.0), SC_(0.0) } },
{ { SC_(1.0), SC_(1.0), SC_(0.4400505857449335159596822037189149131274) } },
{ { SC_(1.0), SC_(-2.0), SC_(-0.5767248077568733872024482422691370869203) } },
@@ -169,7 +172,7 @@ void test_bessel(T, const char* name)
{ { SC_(1.0), SC_(1e-10), SC_(4.999999999999999999993750000000000000000e-11) } },
{ { SC_(1.0), SC_(-1e+01), SC_(-4.347274616886143666974876802585928830627e-02) } },
}};
static const std::array<std::array<typename table_type<T>::type, 3>, 5> j1_tricky = {{
static const std::array<std::array<table_entry_type, 3>, 5> j1_tricky = {{
// Big numbers make the accuracy of std::sin the limiting factor:
{ { SC_(1.0), SC_(1e+03), SC_(4.728311907089523917576071901216916285418e-03) } },
{ { SC_(1.0), SC_(1e+05), SC_(1.846757562882567716362123967114215743694e-03) } },
@@ -179,7 +182,7 @@ void test_bessel(T, const char* name)
{ { SC_(1.0), SC_(10.1734676361083984375) /*T(10667654) / (1024 * 1024)*/, SC_(1.24591331097191900488116495350277530373473085499043086981229e-7) } },
}};
static const std::array<std::array<typename table_type<T>::type, 3>, 17> jn_data = {{
static const std::array<std::array<table_entry_type, 3>, 17> jn_data = {{
// This first one is a modified test case from https://svn.boost.org/trac/boost/ticket/2733
{ { SC_(-1.0), SC_(1.25), SC_(-0.510623260319880467069474837274910375352924050139633057168856) } },
{ { SC_(2.0), SC_(0.0), SC_(0.0) } },
@@ -211,7 +214,7 @@ void test_bessel(T, const char* name)
do_test_cyl_bessel_j_int<T>(j1_tricky, name, "Bessel J1: Mathworld Data (tricky cases) (Integer Version)");
do_test_cyl_bessel_j_int<T>(jn_data, name, "Bessel JN: Mathworld Data (Integer Version)");
static const std::array<std::array<typename table_type<T>::type, 3>, 20> jv_data = {{
static const std::array<std::array<table_entry_type, 3>, 20> jv_data = {{
//SC_(-2.4), {{ SC_(0.0), std::numeric_limits<T>::infinity() }},
{ { SC_(22.5), SC_(0.0), SC_(0.0) } },
{ { SC_(2.3994140625) /*2457.0 / 1024*/, SC_(0.0009765625) /* 1 / 1024*/, SC_(3.80739920118603335646474073457326714709615200130620574875292e-9) } },
@@ -236,7 +239,7 @@ void test_bessel(T, const char* name)
{{ SC_(-8.5), SC_(12.566370614359172953850573533118011536788677597500423283899778369231265625144835994512139301368468271928592346053) /*Pi * 4*/, SC_(-0.257086543428224355151772807588810984369026142375675714560864) }},
}};
do_test_cyl_bessel_j<T>(jv_data, name, "Bessel J: Mathworld Data");
static const std::array<std::array<typename table_type<T>::type, 3>, 4> jv_large_data = {{
static const std::array<std::array<table_entry_type, 3>, 4> jv_large_data = {{
// Bug report https://svn.boost.org/trac/boost/ticket/5560:
{{ SC_(-0.5), SC_(1.2458993688871959419388378518880931736878259938089494331010226962863582408064841833232475731084062642684629e-206) /*static_cast<T>(std::ldexp(0.5, -683))*/, SC_(7.14823099969225685526188875418476476336424046896822867989728e102) }},
{ { SC_(256.0), SC_(512.0), SC_(0.00671672065717513246956991122723250578101154313313749938944675) } },
@@ -261,6 +264,81 @@ void test_bessel(T, const char* name)
//
// Some special cases:
//
// Specific tests for Git-issue1292.
{
using local_ctrl_array_type = std::array<table_entry_type, std::size_t { UINT8_C(43) }>;
// Table[N[BesselJ[3, n/10], 40], {n, 9, 51, 1}]
static const local_ctrl_array_type ctrl_data =
{{
SC_(0.01443402847586617545767791623904539755731), SC_(0.01956335398266840591890532162175150825451), SC_(0.02569452861246328174726417617756888741432), SC_(0.03287433692499494270867882730165246683837),
SC_(0.04113582571991693187673486447516908751463), SC_(0.05049771328895129623567992727476043273558), SC_(0.06096395114113963064394955997646387979571), SC_(0.07252344333261900300034928368068248675877),
SC_(0.08514992694801526415321095754253909148633), SC_(0.09880201565861918291536618746528733463749), SC_(0.1134234066389601112649841240858617923591), SC_(0.1289432494744020510987933329692398352700),
SC_(0.1452766740542063665759023355570418120750), SC_(0.1623254728332874543121706910035271736854), SC_(0.1799789312775334540800304157279732327839), SC_(0.1981147987975668248498434552081155790183),
SC_(0.2166003910391135247666890035159637217168), SC_(0.2352938130489638091015220916013129483423), SC_(0.2540452915872273499615464996563039918262), SC_(0.2726986037216204380267188592437356599939),
SC_(0.2910925878291867784836313080855848616815), SC_(0.3090627222552516436182601949468331494291), SC_(0.3264427561473409695937042738575781129080), SC_(0.3430663764006682009386373318558777864023),
SC_(0.3587688942275418259451574456258027163924), SC_(0.3733889346000900583527754127339797472980), SC_(0.3867701117168813668578718121131100327218), SC_(0.3987626737105880326848194417650226836608),
SC_(0.4092251000454309977422936498249743734653), SC_(0.4180256354477855744864458808409348352597), SC_(0.4250437447674560017637404058105525727991), SC_(0.4301714738756219403581834788533355563393),
SC_(0.4333147025616927046073022200802734463060), SC_(0.4343942763872007823091130214493427347554), SC_(0.4333470055809823422144251313032973397899), SC_(0.4301265203055088083605755042771532591535),
SC_(0.4247039729774556002468140098011553543390), SC_(0.4170685797734672711167804755454067582755), SC_(0.4072279949807128989552790124633945783765), SC_(0.3952085134465309348696666123753181072022),
SC_(0.3810550980268886849843356923521907577982), SC_(0.3648312306136669944635769493587219791343), SC_(0.3466185870197064968846647990300282094299)
}};
const T tolerance { 128 * boost::math::tools::epsilon<T>() };
int n_val { 9 };
for(const T& ctrl: ctrl_data)
{
const T x_val { static_cast<T>(static_cast<T>(n_val) / 10) };
const T jn_val { boost::math::cyl_bessel_j(3, x_val) };
++n_val;
BOOST_CHECK_CLOSE_FRACTION(jn_val, ctrl, tolerance);
}
}
// More specific tests for Git-issue1292.
{
using local_ctrl_array_type = std::array<table_entry_type, std::size_t { UINT8_C(43) }>;
// Table[N[BesselJ[31 / 10, n/10], 40], {n, 9, 51, 1}]
static const local_ctrl_array_type ctrl_data =
{{
SC_(0.01175139795214295170487105485346781171863), SC_(0.01610092560641321584451701371378836908343), SC_(0.02135659148701787280713314897691402425560), SC_(0.02757316602094671775387912184375438913864),
SC_(0.03479372470942323424236519547108889796879), SC_(0.04304884612770317349181083608393216357649), SC_(0.05235595486839477302545989170267853515412), SC_(0.06271881444009116864560457340917173135492),
SC_(0.07412717428914531462554459978389560408816), SC_(0.08655657401374878955779092959288219537482), SC_(0.09996830657138141192006478769855556316995), SC_(0.1143095409011066041623799431143340837007),
SC_(0.1295136029333754366226206053365805922136), SC_(0.1455004124749064242445094865982308151833), SC_(0.1621770719619318324763655518038365467780), SC_(0.1794386015949089605595848208470049166853),
SC_(0.1971688139222575484595939469080319406355), SC_(0.2152413195483352761119610246805962611319), SC_(0.2335206543185642795903443565627832597552), SC_(0.2518635170977563283446184976264924077045),
SC_(0.2701201061202457234361463802484836927464), SC_(0.2881355408650536940851830262098172944660), SC_(0.3057513555072237123913927136839853900197), SC_(0.3228070492275195774383850177821175151772),
SC_(0.3391416780352496831991017115498371039332), SC_(0.3545954722799571667706920253581728149226), SC_(0.3690114637024459111815176585531943143085), SC_(0.3822371057078773326211699424651752096338),
SC_(0.3941258705356621024731438508712990073773), SC_(0.4045388071531719615179931506197074798366), SC_(0.4133460440118967760814988295083688541739), SC_(0.4204282212729730452147271372029342292548),
SC_(0.4256778377298582292448895379334671298297), SC_(0.4290004984236535775073755577697874033289), SC_(0.4303160498540635572666996674883665298558), SC_(0.4295595907277138677145484170304736319185),
SC_(0.4266823473457215250850626209433240464185), SC_(0.4216524040030023831246842156638018726147), SC_(0.4144552801406973126588418824207056743782), SC_(0.4050943474472003316564108983454900189419),
SC_(0.3935910816286283019261303507307813773886), SC_(0.3799851451515116989978414766085547417497), SC_(0.3643342988837623802358278078893847672838)
}};
const T tolerance { 128 * boost::math::tools::epsilon<T>() };
const T vu_val { static_cast<T>(static_cast<T>(31) / 10) };
int n_val { 9 };
for(const T& ctrl : ctrl_data)
{
const T x_val { static_cast<T>(static_cast<T>(n_val) / 10) };
const T jn_val { boost::math::cyl_bessel_j(vu_val, x_val) };
++n_val;
BOOST_CHECK_CLOSE_FRACTION(jn_val, ctrl, tolerance);
}
}
BOOST_CHECK_EQUAL(boost::math::sph_bessel(0, T(0)), T(1));
BOOST_CHECK_EQUAL(boost::math::sph_bessel(1, T(0)), T(0));
BOOST_CHECK_EQUAL(boost::math::sph_bessel(100000, T(0)), T(0));
@@ -301,4 +379,3 @@ void test_bessel(T, const char* name)
BOOST_CHECK_EQUAL(boost::math::cyl_bessel_j(T(2), -std::numeric_limits<T>::infinity()), T(0));
}
}

2
test/test_constant_generate.cpp
View File

@@ -48,6 +48,7 @@ int main()
BOOST_CONSTANTS_GENERATE(four_thirds_pi);
BOOST_CONSTANTS_GENERATE(one_div_two_pi);
BOOST_CONSTANTS_GENERATE(one_div_root_two_pi);
BOOST_CONSTANTS_GENERATE(log_pi);
BOOST_CONSTANTS_GENERATE(root_pi);
BOOST_CONSTANTS_GENERATE(root_half_pi);
BOOST_CONSTANTS_GENERATE(root_two_pi);
@@ -189,4 +190,3 @@ Output
*/

14
test/test_constants.cpp
View File

@@ -103,6 +103,7 @@ void test_spots(RealType)
BOOST_MATH_STD_USING
CHECK_ULP_CLOSE(3.14159265358979323846264338327950288419716939937510L, pi<RealType>(), 2);
CHECK_ULP_CLOSE(log(3.14159265358979323846264338327950288419716939937510L), log_pi<RealType>(), 2);
CHECK_ULP_CLOSE(sqrt(3.14159265358979323846264338327950288419716939937510L), root_pi<RealType>(), 2);
CHECK_ULP_CLOSE(sqrt(3.14159265358979323846264338327950288419716939937510L/2), root_half_pi<RealType>(), 2);
CHECK_ULP_CLOSE(sqrt(3.14159265358979323846264338327950288419716939937510L * 2), root_two_pi<RealType>(), 2);
@@ -160,6 +161,7 @@ void test_spots(RealType)
CHECK_ULP_CLOSE(1 / (3.14159265358979323846264338327950288419716939937510L), one_div_pi<RealType>(), 2);
CHECK_ULP_CLOSE(2 / (3.14159265358979323846264338327950288419716939937510L), two_div_pi<RealType>(), 2);
CHECK_ULP_CLOSE(1 / (2 * 3.14159265358979323846264338327950288419716939937510L), one_div_two_pi<RealType>(), 2);
CHECK_ULP_CLOSE(log(3.14159265358979323846264338327950288419716939937510L), log_pi<RealType>(), 2);
CHECK_ULP_CLOSE(sqrt(3.14159265358979323846264338327950288419716939937510L), root_pi<RealType>(), 2);
CHECK_ULP_CLOSE(sqrt(3.14159265358979323846264338327950288419716939937510L / 2), root_half_pi<RealType>(), 2);
CHECK_ULP_CLOSE(sqrt(2 * 3.14159265358979323846264338327950288419716939937510L), root_two_pi<RealType>(), 2);
@@ -242,6 +244,7 @@ void test_float_spots()
BOOST_MATH_STD_USING
CHECK_ULP_CLOSE(static_cast<float>(3.14159265358979323846264338327950288419716939937510F), pi, 2);
CHECK_ULP_CLOSE(static_cast<float>(log(3.14159265358979323846264338327950288419716939937510F)), log_pi, 2);
CHECK_ULP_CLOSE(static_cast<float>(sqrt(3.14159265358979323846264338327950288419716939937510F)), root_pi, 2);
CHECK_ULP_CLOSE(static_cast<float>(sqrt(3.14159265358979323846264338327950288419716939937510F/2)), root_half_pi, 2);
CHECK_ULP_CLOSE(static_cast<float>(sqrt(3.14159265358979323846264338327950288419716939937510F * 2)), root_two_pi, 2);
@@ -293,6 +296,7 @@ void test_float_spots()
CHECK_ULP_CLOSE(static_cast<float>(1 / (3.14159265358979323846264338327950288419716939937510F)), one_div_pi, 2);
CHECK_ULP_CLOSE(static_cast<float>(2 / (3.14159265358979323846264338327950288419716939937510F)), two_div_pi, 2);
CHECK_ULP_CLOSE(static_cast<float>(1 / (2 * 3.14159265358979323846264338327950288419716939937510F)), one_div_two_pi, 2);
CHECK_ULP_CLOSE(static_cast<float>(log(3.14159265358979323846264338327950288419716939937510F)), log_pi, 2);
CHECK_ULP_CLOSE(static_cast<float>(sqrt(3.14159265358979323846264338327950288419716939937510F)), root_pi, 2);
CHECK_ULP_CLOSE(static_cast<float>(sqrt(3.14159265358979323846264338327950288419716939937510F / 2)), root_half_pi, 2);
CHECK_ULP_CLOSE(static_cast<float>(sqrt(2 * 3.14159265358979323846264338327950288419716939937510F)), root_two_pi, 2);
@@ -362,6 +366,7 @@ void test_f32_spots()
BOOST_MATH_STD_USING
CHECK_ULP_CLOSE(static_cast<std::float32_t>(3.14159265358979323846264338327950288419716939937510F32), pi, 2);
CHECK_ULP_CLOSE(static_cast<std::float32_t>(log(3.14159265358979323846264338327950288419716939937510F32)), log_pi, 2);
CHECK_ULP_CLOSE(static_cast<std::float32_t>(sqrt(3.14159265358979323846264338327950288419716939937510F32)), root_pi, 2);
CHECK_ULP_CLOSE(static_cast<std::float32_t>(sqrt(3.14159265358979323846264338327950288419716939937510F32/2)), root_half_pi, 2);
CHECK_ULP_CLOSE(static_cast<std::float32_t>(sqrt(3.14159265358979323846264338327950288419716939937510F32 * 2)), root_two_pi, 2);
@@ -413,6 +418,7 @@ void test_f32_spots()
CHECK_ULP_CLOSE(static_cast<std::float32_t>(1 / (3.14159265358979323846264338327950288419716939937510F32)), one_div_pi, 2);
CHECK_ULP_CLOSE(static_cast<std::float32_t>(2 / (3.14159265358979323846264338327950288419716939937510F32)), two_div_pi, 2);
CHECK_ULP_CLOSE(static_cast<std::float32_t>(1 / (2 * 3.14159265358979323846264338327950288419716939937510F32)), one_div_two_pi, 2);
CHECK_ULP_CLOSE(static_cast<std::float32_t>(log(3.14159265358979323846264338327950288419716939937510F32)), log_pi, 2);
CHECK_ULP_CLOSE(static_cast<std::float32_t>(sqrt(3.14159265358979323846264338327950288419716939937510F32)), root_pi, 2);
CHECK_ULP_CLOSE(static_cast<std::float32_t>(sqrt(3.14159265358979323846264338327950288419716939937510F32 / 2)), root_half_pi, 2);
CHECK_ULP_CLOSE(static_cast<std::float32_t>(sqrt(2 * 3.14159265358979323846264338327950288419716939937510F32)), root_two_pi, 2);
@@ -482,6 +488,7 @@ void test_double_spots()
BOOST_MATH_STD_USING
CHECK_ULP_CLOSE(static_cast<double>(3.14159265358979323846264338327950288419716939937510), pi, 2);
CHECK_ULP_CLOSE(static_cast<double>(log(3.14159265358979323846264338327950288419716939937510)), log_pi, 2);
CHECK_ULP_CLOSE(static_cast<double>(sqrt(3.14159265358979323846264338327950288419716939937510)), root_pi, 2);
CHECK_ULP_CLOSE(static_cast<double>(sqrt(3.14159265358979323846264338327950288419716939937510/2)), root_half_pi, 2);
CHECK_ULP_CLOSE(static_cast<double>(sqrt(3.14159265358979323846264338327950288419716939937510 * 2)), root_two_pi, 2);
@@ -533,6 +540,7 @@ void test_double_spots()
CHECK_ULP_CLOSE(static_cast<double>(1 / (3.14159265358979323846264338327950288419716939937510)), one_div_pi, 2);
CHECK_ULP_CLOSE(static_cast<double>(2 / (3.14159265358979323846264338327950288419716939937510)), two_div_pi, 2);
CHECK_ULP_CLOSE(static_cast<double>(1 / (2 * 3.14159265358979323846264338327950288419716939937510)), one_div_two_pi, 2);
CHECK_ULP_CLOSE(static_cast<double>(log(3.14159265358979323846264338327950288419716939937510)), log_pi, 2);
CHECK_ULP_CLOSE(static_cast<double>(sqrt(3.14159265358979323846264338327950288419716939937510)), root_pi, 2);
CHECK_ULP_CLOSE(static_cast<double>(sqrt(3.14159265358979323846264338327950288419716939937510 / 2)), root_half_pi, 2);
CHECK_ULP_CLOSE(static_cast<double>(sqrt(2 * 3.14159265358979323846264338327950288419716939937510)), root_two_pi, 2);
@@ -601,6 +609,7 @@ void test_f64_spots()
BOOST_MATH_STD_USING
CHECK_ULP_CLOSE(static_cast<std::float64_t>(3.14159265358979323846264338327950288419716939937510F64), pi, 2);
CHECK_ULP_CLOSE(static_cast<std::float64_t>(log(3.14159265358979323846264338327950288419716939937510F64)), log_pi, 2);
CHECK_ULP_CLOSE(static_cast<std::float64_t>(sqrt(3.14159265358979323846264338327950288419716939937510F64)), root_pi, 2);
CHECK_ULP_CLOSE(static_cast<std::float64_t>(sqrt(3.14159265358979323846264338327950288419716939937510F64/2)), root_half_pi, 2);
CHECK_ULP_CLOSE(static_cast<std::float64_t>(sqrt(3.14159265358979323846264338327950288419716939937510F64 * 2)), root_two_pi, 2);
@@ -652,6 +661,7 @@ void test_f64_spots()
CHECK_ULP_CLOSE(static_cast<std::float64_t>(1 / (3.14159265358979323846264338327950288419716939937510F64)), one_div_pi, 2);
CHECK_ULP_CLOSE(static_cast<std::float64_t>(2 / (3.14159265358979323846264338327950288419716939937510F64)), two_div_pi, 2);
CHECK_ULP_CLOSE(static_cast<std::float64_t>(1 / (2 * 3.14159265358979323846264338327950288419716939937510F64)), one_div_two_pi, 2);
CHECK_ULP_CLOSE(static_cast<std::float64_t>(log(3.14159265358979323846264338327950288419716939937510F64)), log_pi, 2);
CHECK_ULP_CLOSE(static_cast<std::float64_t>(sqrt(3.14159265358979323846264338327950288419716939937510F64)), root_pi, 2);
CHECK_ULP_CLOSE(static_cast<std::float64_t>(sqrt(3.14159265358979323846264338327950288419716939937510F64 / 2)), root_half_pi, 2);
CHECK_ULP_CLOSE(static_cast<std::float64_t>(sqrt(2 * 3.14159265358979323846264338327950288419716939937510F64)), root_two_pi, 2);
@@ -724,6 +734,7 @@ void test_long_double_spots()
BOOST_MATH_STD_USING
CHECK_ULP_CLOSE(static_cast<long double>(3.14159265358979323846264338327950288419716939937510L), pi, 2);
CHECK_ULP_CLOSE(static_cast<long double>(log(3.14159265358979323846264338327950288419716939937510L)), log_pi, 2);
CHECK_ULP_CLOSE(static_cast<long double>(sqrt(3.14159265358979323846264338327950288419716939937510L)), root_pi, 2);
CHECK_ULP_CLOSE(static_cast<long double>(sqrt(3.14159265358979323846264338327950288419716939937510L/2)), root_half_pi, 2);
CHECK_ULP_CLOSE(static_cast<long double>(sqrt(3.14159265358979323846264338327950288419716939937510L * 2)), root_two_pi, 2);
@@ -775,6 +786,7 @@ void test_long_double_spots()
CHECK_ULP_CLOSE(static_cast<long double>(1 / (3.14159265358979323846264338327950288419716939937510L)), one_div_pi, 2);
CHECK_ULP_CLOSE(static_cast<long double>(2 / (3.14159265358979323846264338327950288419716939937510L)), two_div_pi, 2);
CHECK_ULP_CLOSE(static_cast<long double>(1 / (2 * 3.14159265358979323846264338327950288419716939937510L)), one_div_two_pi, 2);
CHECK_ULP_CLOSE(static_cast<long double>(log(3.14159265358979323846264338327950288419716939937510L)), log_pi, 2);
CHECK_ULP_CLOSE(static_cast<long double>(sqrt(3.14159265358979323846264338327950288419716939937510L)), root_pi, 2);
CHECK_ULP_CLOSE(static_cast<long double>(sqrt(3.14159265358979323846264338327950288419716939937510L / 2)), root_half_pi, 2);
CHECK_ULP_CLOSE(static_cast<long double>(sqrt(2 * 3.14159265358979323846264338327950288419716939937510L)), root_two_pi, 2);
@@ -851,6 +863,7 @@ void test_real_concept_policy(const Policy&)
BOOST_MATH_STD_USING
CHECK_ULP_CLOSE(3.14159265358979323846264338327950288419716939937510L, (pi<real_concept, Policy>)(), 2);
CHECK_ULP_CLOSE(log(3.14159265358979323846264338327950288419716939937510L), (log_pi<real_concept, Policy>)(), 2);
CHECK_ULP_CLOSE(sqrt(3.14159265358979323846264338327950288419716939937510L), (root_pi<real_concept, Policy>)(), 2);
CHECK_ULP_CLOSE(sqrt(3.14159265358979323846264338327950288419716939937510L/2), (root_half_pi<real_concept, Policy>)(), 2);
CHECK_ULP_CLOSE(sqrt(3.14159265358979323846264338327950288419716939937510L * 2), (root_two_pi<real_concept, Policy>)(), 2);
@@ -902,6 +915,7 @@ void test_real_concept_policy(const Policy&)
CHECK_ULP_CLOSE(1 / (3.14159265358979323846264338327950288419716939937510L), (one_div_pi<real_concept, Policy>)(), 2);
CHECK_ULP_CLOSE(2 / (3.14159265358979323846264338327950288419716939937510L), (two_div_pi<real_concept, Policy>)(), 2);
CHECK_ULP_CLOSE(1 / (2 * 3.14159265358979323846264338327950288419716939937510L), (one_div_two_pi<real_concept, Policy>)(), 2);
CHECK_ULP_CLOSE(log(3.14159265358979323846264338327950288419716939937510L), (log_pi<real_concept, Policy>)(), 2);
CHECK_ULP_CLOSE(sqrt(3.14159265358979323846264338327950288419716939937510L), (root_pi<real_concept, Policy>)(), 2);
CHECK_ULP_CLOSE(sqrt(3.14159265358979323846264338327950288419716939937510L / 2), (root_half_pi<real_concept, Policy>)(), 2);
CHECK_ULP_CLOSE(sqrt(2 * 3.14159265358979323846264338327950288419716939937510L), (root_two_pi<real_concept, Policy>)(), 2);

91
test/test_logistic_sigmoid.cpp Normal file
View File

@@ -0,0 +1,91 @@
// Copyright Matt Borland 2025.
// Use, modification and distribution are subject to 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)
#include <boost/math/special_functions/logistic_sigmoid.hpp>
#include <boost/multiprecision/cpp_bin_float.hpp>
#include <boost/multiprecision/cpp_dec_float.hpp>
#include "math_unit_test.hpp"
#include <array>
#include <cfloat>
#include <cfenv>
#pragma STDC FENV_ACCESS ON
template <typename RealType>
void test()
{
const std::array<RealType, 5> x_values = {
0,
1,
1000,
0.5,
0.75
};
const std::array<RealType, 5> y_values = {
static_cast<RealType>(1) / 2,
static_cast<RealType>(0.73105857863000487925115924182183627436514464016505651927636590791904045307),
static_cast<RealType>(1),
static_cast<RealType>(0.62245933120185456463890056574550847875327936530891016305943716265854500),
static_cast<RealType>(0.6791786991753929731596801157765790212342212482195760219829517436805)
};
for (std::size_t i = 0; i < x_values.size(); ++i)
{
const RealType test_value {boost::math::logistic_sigmoid(x_values[i])};
BOOST_MATH_IF_CONSTEXPR (std::is_same<RealType, float>::value || std::is_same<RealType, double>::value)
{
CHECK_ULP_CLOSE(test_value, y_values[i], 1);
}
else
{
RealType comparison_value = y_values[i];
CHECK_MOLLIFIED_CLOSE(test_value, comparison_value, static_cast<RealType>(1e-15));
}
bool fe {false};
if (std::fetestexcept(FE_OVERFLOW))
{
fe = true; // LCOV_EXCL_LINE
std::cerr << "FE_OVERFLOW" << std::endl; // LCOV_EXCL_LINE
}
if (std::fetestexcept(FE_UNDERFLOW))
{
fe = true; // LCOV_EXCL_LINE
std::cerr << "FE_UNDERFLOW" << std::endl; // LCOV_EXCL_LINE
}
if (std::fetestexcept(FE_DIVBYZERO))
{
fe = true; // LCOV_EXCL_LINE
std::cerr << "FE_DIVBYZERO" << std::endl; // LCOV_EXCL_LINE
}
if (std::fetestexcept(FE_INVALID))
{
fe = true; // LCOV_EXCL_LINE
std::cerr << "FE_INVALID" << std::endl; // LCOV_EXCL_LINE
}
CHECK_EQUAL(fe, false);
}
}
int main()
{
std::feclearexcept(FE_ALL_EXCEPT);
test<float>();
std::feclearexcept(FE_ALL_EXCEPT);
test<double>();
std::feclearexcept(FE_ALL_EXCEPT);
test<long double>();
std::feclearexcept(FE_ALL_EXCEPT);
test<boost::multiprecision::cpp_bin_float_quad>();
test<boost::multiprecision::cpp_dec_float_50>();
return boost::math::test::report_errors();
}

107
test/test_logit.cpp Normal file
View File

@@ -0,0 +1,107 @@
// Copyright Matt Borland 2025.
// Use, modification and distribution are subject to 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)
#include <boost/math/special_functions/logit.hpp>
#include <boost/multiprecision/cpp_bin_float.hpp>
#include "math_unit_test.hpp"
#include <array>
#include <limits>
#include <cfenv>
#pragma STDC FENV_ACCESS ON
template <typename RealType>
void test()
{
const std::array<RealType, 5> x_values = {
0.01,
0.24,
0.5,
0.75,
0.995,
};
const std::array<RealType, 5> y_values = {
RealType{-4.595119850134589926852434051810180709116687969582916078687956376405},
RealType{-1.15267950993838545919655007350715126451438856911612411268258589327840479},
RealType{0},
RealType{1.09861228866810969139524523692252570464749055782274945173469433363749429},
RealType{5.2933048247244923954101212918685372018911052805694724989064609879440992}
};
for (std::size_t i = 0; i < x_values.size(); ++i)
{
const RealType test_value {boost::math::logit(x_values[i])};
BOOST_MATH_IF_CONSTEXPR (std::is_same<RealType, float>::value || std::is_same<RealType, double>::value)
{
CHECK_ULP_CLOSE(test_value, y_values[i], 5);
}
else
{
CHECK_MOLLIFIED_CLOSE(test_value, y_values[i], 1e-15);
}
bool fe {false};
if (std::fetestexcept(FE_OVERFLOW))
{
fe = true; // LCOV_EXCL_LINE
std::cerr << "FE_OVERFLOW" << std::endl; // LCOV_EXCL_LINE
}
if (std::fetestexcept(FE_UNDERFLOW))
{
fe = true; // LCOV_EXCL_LINE
std::cerr << "FE_UNDERFLOW" << std::endl; // LCOV_EXCL_LINE
}
if (std::fetestexcept(FE_DIVBYZERO))
{
fe = true; // LCOV_EXCL_LINE
std::cerr << "FE_DIVBYZERO" << std::endl; // LCOV_EXCL_LINE
}
if (std::fetestexcept(FE_INVALID))
{
fe = true; // LCOV_EXCL_LINE
std::cerr << "FE_INVALID" << std::endl; // LCOV_EXCL_LINE
}
CHECK_EQUAL(fe, false);
}
#if defined(_CPPUNWIND) || defined(__EXCEPTIONS)
BOOST_MATH_IF_CONSTEXPR (std::is_arithmetic<RealType>::value)
{
bool thrown {false};
try
{
boost::math::logit(std::numeric_limits<RealType>::denorm_min());
}
catch (...)
{
thrown = true;
}
CHECK_EQUAL(thrown, true);
}
#endif // Exceptional environments
}
int main()
{
std::feclearexcept(FE_ALL_EXCEPT);
test<float>();
std::feclearexcept(FE_ALL_EXCEPT);
test<double>();
std::feclearexcept(FE_ALL_EXCEPT);
test<long double>();
std::feclearexcept(FE_ALL_EXCEPT);
test<boost::multiprecision::cpp_bin_float_quad>();
return boost::math::test::report_errors();
}