2
0
mirror of https://github.com/boostorg/math.git synced 2026-02-10 23:42:23 +00:00

1F1: hypergeometric_1F1_small_a_negative_b_by_ratio doesn't work when b + i == a for some integer i.

This commit is contained in:
jzmaddock
2019-02-16 12:13:07 +00:00
parent 00dfd44470
commit bd9599fe3f
2 changed files with 12 additions and 3 deletions

View File

@@ -45,13 +45,14 @@
// We grab the ratio for M[a, b, z] / M[a, b+1, z] and use it to seed 2 initial values,
// then recurse until b > 0, compute a reference value and normalize (Millers method).
//
int iterations = boost::math::itrunc(-b, pol);
boost::uintmax_t max_iter = boost::math::policies::get_max_series_iterations<Policy>();
T ratio = boost::math::tools::function_ratio_from_forwards_recurrence(boost::math::detail::hypergeometric_1F1_recurrence_b_coefficients<T>(a, b, z), boost::math::tools::epsilon<T>(), max_iter);
boost::math::policies::check_series_iterations<T>("boost::math::hypergeometric_1F1_small_a_negative_b_by_ratio<%1%>(%1%,%1%,%1%)", max_iter, pol);
T first = 1;
T second = 1 / ratio;
int iterations = boost::math::itrunc(-b, pol);
int scaling1 = 0;
BOOST_ASSERT(b + iterations != a);
second = boost::math::tools::apply_recurrence_relation_forward(boost::math::detail::hypergeometric_1F1_recurrence_b_coefficients<T>(a, b + 1, z), iterations, first, second, &scaling1);
int scaling2 = 0;
first = hypergeometric_1F1_imp(a, b + iterations + 1, z, pol, scaling2);

View File

@@ -256,7 +256,16 @@ namespace boost { namespace math { namespace detail {
//
// This is a tricky area, potentially we have no good method at all:
//
if (max_b_for_1F1_small_a_negative_b_by_ratio(z) < b)
if (b - ceil(b) == a)
{
// Fractional parts of a and b are genuinely equal, we might as well
// apply Kummer's relation and get a truncated series:
int scaling = itrunc(z);
T r = exp(z - scaling) * detail::hypergeometric_1F1_imp<T>(b_minus_a, b, -z, pol, log_scaling);
log_scaling += scaling;
return r;
}
if ((b < -1) && (max_b_for_1F1_small_a_negative_b_by_ratio(z) < b))
return hypergeometric_1F1_small_a_negative_b_by_ratio(a, b, z, pol, log_scaling);
if ((b > -1) && (b < -0.5f))
{
@@ -408,7 +417,6 @@ namespace boost { namespace math { namespace detail {
// 13.3.6 appears to be the most efficient and often the most accurate method.
return boost::math::detail::hypergeometric_1F1_AS_13_3_6(a, b, z, b_minus_a, pol, log_scaling);
}
if ((a > 0) && (b > 0) && (a * z / b > 2))
{
//