Skip to content

Commit 5130435

Browse files
jzmaddockjzmaddock
authored
Merge pull request #1111 from boostorg/improve_coverage_3
Improve erf/expm1/expint coverage.
2 parents 01dbedc + 3e3ac03 commit 5130435

Some content is hidden

Large Commits have some content hidden by default. Use the searchbox below for content that may be hidden.

56 files changed

+1326
-1176
lines changed

‎include/boost/math/special_functions/detail/bessel_ik.hpp

+1-1
Original file line numberDiff line numberDiff line change
@@ -416,7 +416,7 @@ BOOST_MATH_GPU_ENABLED int bessel_ik(T v, T x, T* result_I, T* result_K, int kin
416416
else
417417
Iv = boost::math::numeric_limits<T>::quiet_NaN(); // any value will do
418418
}
419-
if (reflect)
419+
if (reflect && (kind & need_i))
420420
{
421421
BOOST_MATH_ASSERT(fabs(v - n - u) < tools::forth_root_epsilon<T>());
422422
T z = (u + n % 2);

‎include/boost/math/special_functions/detail/bessel_j1.hpp

-32
Original file line numberDiff line numberDiff line change
@@ -32,41 +32,9 @@
3232

3333
namespace boost { namespace math{ namespace detail{
3434

35-
template <typename T>
36-
BOOST_MATH_GPU_ENABLED T bessel_j1(T x);
37-
38-
template <class T>
39-
struct bessel_j1_initializer
40-
{
41-
struct init
42-
{
43-
BOOST_MATH_GPU_ENABLED init()
44-
{
45-
do_init();
46-
}
47-
BOOST_MATH_GPU_ENABLED static void do_init()
48-
{
49-
bessel_j1(T(1));
50-
}
51-
BOOST_MATH_GPU_ENABLED void force_instantiate()const{}
52-
};
53-
BOOST_MATH_STATIC const init initializer;
54-
BOOST_MATH_GPU_ENABLED static void force_instantiate()
55-
{
56-
#ifndef BOOST_MATH_HAS_GPU_SUPPORT
57-
initializer.force_instantiate();
58-
#endif
59-
}
60-
};
61-
62-
template <class T>
63-
const typename bessel_j1_initializer<T>::init bessel_j1_initializer<T>::initializer;
64-
6535
template <typename T>
6636
BOOST_MATH_GPU_ENABLED T bessel_j1(T x)
6737
{
68-
bessel_j1_initializer<T>::force_instantiate();
69-
7038
BOOST_MATH_STATIC const T P1[] = {
7139
static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, -1.4258509801366645672e+11)),
7240
static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 6.6781041261492395835e+09)),

‎include/boost/math/special_functions/detail/bessel_k0.hpp

-40
Original file line numberDiff line numberDiff line change
@@ -47,45 +47,6 @@
4747

4848
namespace boost { namespace math { namespace detail{
4949

50-
template <typename T>
51-
BOOST_MATH_GPU_ENABLED T bessel_k0(const T& x);
52-
53-
template <class T, class tag>
54-
struct bessel_k0_initializer
55-
{
56-
struct init
57-
{
58-
BOOST_MATH_GPU_ENABLED init()
59-
{
60-
do_init(tag());
61-
}
62-
BOOST_MATH_GPU_ENABLED static void do_init(const boost::math::integral_constant<int, 113>&)
63-
{
64-
bessel_k0(T(0.5));
65-
bessel_k0(T(1.5));
66-
}
67-
BOOST_MATH_GPU_ENABLED static void do_init(const boost::math::integral_constant<int, 64>&)
68-
{
69-
bessel_k0(T(0.5));
70-
bessel_k0(T(1.5));
71-
}
72-
template <class U>
73-
BOOST_MATH_GPU_ENABLED static void do_init(const U&){}
74-
BOOST_MATH_GPU_ENABLED void force_instantiate()const{}
75-
};
76-
BOOST_MATH_STATIC const init initializer;
77-
BOOST_MATH_GPU_ENABLED static void force_instantiate()
78-
{
79-
#ifndef BOOST_MATH_HAS_GPU_SUPPORT
80-
initializer.force_instantiate();
81-
#endif
82-
}
83-
};
84-
85-
template <class T, class tag>
86-
const typename bessel_k0_initializer<T, tag>::init bessel_k0_initializer<T, tag>::initializer;
87-
88-
8950
template <typename T, int N>
9051
BOOST_MATH_GPU_ENABLED T bessel_k0_imp(const T&, const boost::math::integral_constant<int, N>&)
9152
{
@@ -511,7 +472,6 @@ BOOST_MATH_GPU_ENABLED inline T bessel_k0(const T& x)
511472
113 : -1
512473
> tag_type;
513474

514-
bessel_k0_initializer<T, tag_type>::force_instantiate();
515475
return bessel_k0_imp(x, tag_type());
516476
}
517477

‎include/boost/math/special_functions/detail/bessel_k1.hpp

-41
Original file line numberDiff line numberDiff line change
@@ -47,46 +47,6 @@
4747

4848
namespace boost { namespace math { namespace detail{
4949

50-
template <typename T>
51-
BOOST_MATH_GPU_ENABLED T bessel_k1(const T&);
52-
53-
template <class T, class tag>
54-
struct bessel_k1_initializer
55-
{
56-
struct init
57-
{
58-
BOOST_MATH_GPU_ENABLED init()
59-
{
60-
do_init(tag());
61-
}
62-
BOOST_MATH_GPU_ENABLED static void do_init(const boost::math::integral_constant<int, 113>&)
63-
{
64-
bessel_k1(T(0.5));
65-
bessel_k1(T(2));
66-
bessel_k1(T(6));
67-
}
68-
BOOST_MATH_GPU_ENABLED static void do_init(const boost::math::integral_constant<int, 64>&)
69-
{
70-
bessel_k1(T(0.5));
71-
bessel_k1(T(6));
72-
}
73-
template <class U>
74-
BOOST_MATH_GPU_ENABLED static void do_init(const U&) {}
75-
BOOST_MATH_GPU_ENABLED void force_instantiate()const {}
76-
};
77-
BOOST_MATH_STATIC const init initializer;
78-
BOOST_MATH_GPU_ENABLED static void force_instantiate()
79-
{
80-
#ifndef BOOST_MATH_HAS_GPU_SUPPORT
81-
initializer.force_instantiate();
82-
#endif
83-
}
84-
};
85-
86-
template <class T, class tag>
87-
const typename bessel_k1_initializer<T, tag>::init bessel_k1_initializer<T, tag>::initializer;
88-
89-
9050
template <typename T, int N>
9151
inline BOOST_MATH_GPU_ENABLED T bessel_k1_imp(const T&, const boost::math::integral_constant<int, N>&)
9252
{
@@ -553,7 +513,6 @@ namespace boost { namespace math { namespace detail{
553513
113 : -1
554514
> tag_type;
555515

556-
bessel_k1_initializer<T, tag_type>::force_instantiate();
557516
return bessel_k1_imp(x, tag_type());
558517
}
559518

‎include/boost/math/special_functions/detail/hypergeometric_0F1_bessel.hpp

+4-5
Original file line numberDiff line numberDiff line change
@@ -20,12 +20,11 @@
2020
{
2121
BOOST_MATH_STD_USING
2222

23-
const bool is_z_nonpositive = z <= 0;
23+
//const bool is_z_nonpositive = z <= 0;
24+
BOOST_MATH_ASSERT(z < 0); // condition used at call site
2425

25-
const T sqrt_z = is_z_nonpositive ? T(sqrt(-z)) : T(sqrt(z));
26-
const T bessel_mult = is_z_nonpositive ?
27-
boost::math::cyl_bessel_j(b - 1, 2 * sqrt_z, pol) :
28-
boost::math::cyl_bessel_i(b - 1, 2 * sqrt_z, pol) ;
26+
const T sqrt_z = sqrt(-z);
27+
const T bessel_mult = boost::math::cyl_bessel_j(b - 1, 2 * sqrt_z, pol);
2928

3029
if (b > boost::math::max_factorial<T>::value)
3130
{

‎include/boost/math/special_functions/detail/hypergeometric_1F1_bessel.hpp

+24-3
Original file line numberDiff line numberDiff line change
@@ -74,8 +74,12 @@
7474
{
7575
// We get very limited precision due to rapid denormalisation/underflow of the Bessel values, raise an exception and try something else:
7676
policies::raise_evaluation_error("hypergeometric_1F1_AS_13_3_7_tricomi_series<%1%>", "Underflow in Bessel functions", bessel_cache[cache_size - 1], pol);
77+
// Exceptions are off if we get here, just fill the cache with NaN's and we'll let this method fail and fallback later:
78+
std::fill(bessel_cache.begin(), bessel_cache.end(), std::numeric_limits<T>::quiet_NaN());
79+
cache_offset = -cache_size;
80+
return;
7781
}
78-
if ((term * bessel_cache[cache_size - 1] < tools::min_value<T>() / (tools::epsilon<T>() * tools::epsilon<T>())) || !(boost::math::isfinite)(term) || (!std::numeric_limits<T>::has_infinity && (fabs(term) > tools::max_value<T>())))
82+
if ((fabs(term * bessel_cache[cache_size - 1]) < tools::min_value<T>() / (tools::epsilon<T>() * tools::epsilon<T>())) || !(boost::math::isfinite)(term) || (!std::numeric_limits<T>::has_infinity && (fabs(term) > tools::max_value<T>())))
7983
{
8084
term = -log(fabs(bessel_arg)) * b_minus_1_plus_n / 2;
8185
log_scale = lltrunc(term);
@@ -88,15 +92,27 @@
8892
if constexpr (std::numeric_limits<T>::has_infinity)
8993
{
9094
if (!(boost::math::isfinite)(bessel_cache[cache_size - 1]))
95+
{
9196
policies::raise_evaluation_error("hypergeometric_1F1_AS_13_3_7_tricomi_series<%1%>", "Expected finite Bessel function result but got %1%", bessel_cache[cache_size - 1], pol);
97+
// Exceptions are off if we get here, just fill the cache with NaN's and we'll let this method fail and fallback later:
98+
std::fill(bessel_cache.begin(), bessel_cache.end(), std::numeric_limits<T>::quiet_NaN());
99+
}
92100
}
93101
else
94102
if ((boost::math::isnan)(bessel_cache[cache_size - 1]) || (fabs(bessel_cache[cache_size - 1]) >= tools::max_value<T>()))
103+
{
95104
policies::raise_evaluation_error("hypergeometric_1F1_AS_13_3_7_tricomi_series<%1%>", "Expected finite Bessel function result but got %1%", bessel_cache[cache_size - 1], pol);
105+
// Exceptions are off if we get here, just fill the cache with NaN's and we'll let this method fail and fallback later:
106+
std::fill(bessel_cache.begin(), bessel_cache.end(), std::numeric_limits<T>::quiet_NaN());
107+
}
96108
#else
97109
if ((std::numeric_limits<T>::has_infinity && !(boost::math::isfinite)(bessel_cache[cache_size - 1]))
98110
|| (!std::numeric_limits<T>::has_infinity && ((boost::math::isnan)(bessel_cache[cache_size - 1]) || (fabs(bessel_cache[cache_size - 1]) >= tools::max_value<T>()))))
111+
{
99112
policies::raise_evaluation_error("hypergeometric_1F1_AS_13_3_7_tricomi_series<%1%>", "Expected finite Bessel function result but got %1%", bessel_cache[cache_size - 1], pol);
113+
// Exceptions are off if we get here, just fill the cache with NaN's and we'll let this method fail and fallback later:
114+
std::fill(bessel_cache.begin(), bessel_cache.end(), std::numeric_limits<T>::quiet_NaN());
115+
}
100116
#endif
101117
cache_offset = -cache_size;
102118
refill_cache();
@@ -108,8 +124,13 @@
108124
// very small (or zero) when b == 2a:
109125
//
110126
BOOST_MATH_STD_USING
127+
//
128+
// Except in the multiprecision case, we have probably illiminated anything
129+
// would need more than the default 64 Bessel Functions. Anything more
130+
// than that risks becoming a divergent series anyway...
131+
//
111132
if(n - 2 - cache_offset >= cache_size)
112-
refill_cache();
133+
refill_cache(); // LCOV_EXCL_LINE
113134
T result = A_minus_2 * term * bessel_cache[n - 2 - cache_offset];
114135
term *= mult;
115136
++n;
@@ -122,7 +143,7 @@
122143
if (A_minus_2 != 0)
123144
{
124145
if (n - 2 - cache_offset >= cache_size)
125-
refill_cache();
146+
refill_cache(); // LCOV_EXCL_LINE
126147
result += A_minus_2 * term * bessel_cache[n - 2 - cache_offset];
127148
}
128149
term *= mult;

‎include/boost/math/special_functions/detail/igamma_large.hpp

+9
Original file line numberDiff line numberDiff line change
@@ -95,6 +95,7 @@ BOOST_MATH_GPU_ENABLED T igamma_temme_large(T a, T x, const Policy& pol, const b
9595

9696
T workspace[13];
9797

98+
// LCOV_EXCL_START
9899
BOOST_MATH_STATIC const T C0[] = {
99100
BOOST_MATH_BIG_CONSTANT(T, 64, -0.333333333333333333333),
100101
BOOST_MATH_BIG_CONSTANT(T, 64, 0.0833333333333333333333),
@@ -272,6 +273,8 @@ BOOST_MATH_GPU_ENABLED T igamma_temme_large(T a, T x, const Policy& pol, const b
272273
BOOST_MATH_BIG_CONSTANT(T, 64, 0.00640336283380806979482),
273274
BOOST_MATH_BIG_CONSTANT(T, 64, -0.00404101610816766177474),
274275
};
276+
// LCOV_EXCL_END
277+
275278
workspace[12] = tools::evaluate_polynomial(C12, z);
276279

277280
T result = tools::evaluate_polynomial<13, T, T>(workspace, 1/a);
@@ -303,6 +306,7 @@ BOOST_MATH_GPU_ENABLED T igamma_temme_large(T a, T x, const Policy& pol, const b
303306

304307
T workspace[10];
305308

309+
// LCOV_EXCL_START
306310
BOOST_MATH_STATIC const T C0[] = {
307311
static_cast<T>(-0.33333333333333333L),
308312
static_cast<T>(0.083333333333333333L),
@@ -416,6 +420,7 @@ BOOST_MATH_GPU_ENABLED T igamma_temme_large(T a, T x, const Policy& pol, const b
416420
static_cast<T>(0.00083949872067208728L),
417421
static_cast<T>(-0.00043829709854172101L),
418422
};
423+
// LCOV_EXCL_END
419424
workspace[8] = tools::evaluate_polynomial(C8, z);
420425
workspace[9] = static_cast<T>(-0.00059676129019274625L);
421426

@@ -456,6 +461,7 @@ BOOST_MATH_GPU_ENABLED T igamma_temme_large(T a, T x, const Policy& pol, const b
456461

457462
T workspace[3];
458463

464+
// LCOV_EXCL_START
459465
BOOST_MATH_STATIC const T C0[] = {
460466
static_cast<T>(-0.333333333L),
461467
static_cast<T>(0.0833333333L),
@@ -482,6 +488,7 @@ BOOST_MATH_GPU_ENABLED T igamma_temme_large(T a, T x, const Policy& pol, const b
482488
static_cast<T>(0.000771604938L),
483489
};
484490
workspace[2] = tools::evaluate_polynomial(C2, z);
491+
// LCOV_EXCL_END
485492

486493
T result = tools::evaluate_polynomial(workspace, 1/a);
487494
result *= exp(-y) / sqrt(2 * constants::pi<T>() * a);
@@ -510,6 +517,7 @@ BOOST_MATH_GPU_ENABLED T igamma_temme_large(T a, T x, const Policy& pol, const b
510517
// It's use for a < 200 is not recommended, that would
511518
// require many more terms in the polynomials.
512519
//
520+
// LCOV_EXCL_START: 128-bit floats not deliberately tested in our coverage tests (takes too long)
513521
#ifndef BOOST_MATH_HAS_GPU_SUPPORT
514522

515523
template <class T, class Policy>
@@ -802,6 +810,7 @@ BOOST_MATH_GPU_ENABLED T igamma_temme_large(T a, T x, const Policy& pol, const b
802810

803811
return result;
804812
}
813+
// LCOV_EXCL_END
805814

806815
#endif
807816

0 commit comments

Comments
 (0)