From c159693d7efe6d7c1e5ffbe47793e99bcbd430c4 Mon Sep 17 00:00:00 2001 From: jzmaddock Date: Mon, 5 Aug 2019 13:17:42 +0100 Subject: [PATCH] 1F1: tidy up docs and suppress/fix some compiler warnings. --- doc/sf/hypergeometric.qbk | 2 +- .../special_functions/detail/bessel_k0.hpp | 16 ++++++++-------- .../detail/hypergeometric_1F1_bessel.hpp | 18 +++++++++--------- .../detail/hypergeometric_1F1_large_abz.hpp | 2 +- .../hypergeometric_1F1_negative_b_regions.hpp | 10 +++++----- .../detail/hypergeometric_1F1_recurrence.hpp | 2 +- .../hypergeometric_pFq_checked_series.hpp | 2 +- .../detail/hypergeometric_pade.hpp | 2 +- .../detail/hypergeometric_series.hpp | 3 --- .../special_functions/hypergeometric_1F1.hpp | 6 +++--- include/boost/math/tools/recurrence.hpp | 4 ++-- include/boost/math/tools/roots.hpp | 5 ++--- 12 files changed, 34 insertions(+), 38 deletions(-) diff --git a/doc/sf/hypergeometric.qbk b/doc/sf/hypergeometric.qbk index 7dbe65d4e..77c727453 100644 --- a/doc/sf/hypergeometric.qbk +++ b/doc/sf/hypergeometric.qbk @@ -260,7 +260,7 @@ Note however, that the some of the methods used here fail to converge to much be [? __build_html ''''''] [?! __build_html [$../graphs/hypergeometric_1f1/negative_b.png]] -For both /a/ and /b/ less than zero, the errors the worst errors are clustered in a "difficult zone" with /b < a/ and /z/ [cong] /a/: +For both /a/ and /b/ less than zero, the errors the worst errors are clustered in a "difficult zone" with /b < a/ and /z/ [asymp] /a/: [$../graphs/hypergeometric_1f1/negative_ab_bins.svg] [? __build_html ''''''] diff --git a/include/boost/math/special_functions/detail/bessel_k0.hpp b/include/boost/math/special_functions/detail/bessel_k0.hpp index c24d89093..d3df18255 100644 --- a/include/boost/math/special_functions/detail/bessel_k0.hpp +++ b/include/boost/math/special_functions/detail/bessel_k0.hpp @@ -135,17 +135,17 @@ T bessel_k0_imp(const T& x, const mpl::int_<24>&) static const T P[] = { - 2.533141220e-01, - 5.221502603e-01, - 6.380180669e-02, - -5.934976547e-02 + 2.533141220e-01f, + 5.221502603e-01f, + 6.380180669e-02f, + -5.934976547e-02f }; static const T Q[] = { - 1.000000000e+00, - 2.679722431e+00, - 1.561635813e+00, - 1.573660661e-01 + 1.000000000e+00f, + 2.679722431e+00f, + 1.561635813e+00f, + 1.573660661e-01f }; if(x < tools::log_max_value()) return ((tools::evaluate_rational(P, Q, T(1 / x)) + 1) * exp(-x) / sqrt(x)); diff --git a/include/boost/math/special_functions/detail/hypergeometric_1F1_bessel.hpp b/include/boost/math/special_functions/detail/hypergeometric_1F1_bessel.hpp index 74504838c..0febf6f9a 100644 --- a/include/boost/math/special_functions/detail/hypergeometric_1F1_bessel.hpp +++ b/include/boost/math/special_functions/detail/hypergeometric_1F1_bessel.hpp @@ -184,11 +184,11 @@ // if ((j < cache_size - 2) && (tools::max_value() / fabs(64 * bessel_cache[j] / bessel_cache[j + 1]) < fabs(bessel_cache[j]))) { - T scale = pow(fabs(bessel_cache[j] / bessel_cache[j + 1]), j + 1) * 2; - if (!((boost::math::isfinite)(scale))) - scale = tools::max_value(); + T rescale = pow(fabs(bessel_cache[j] / bessel_cache[j + 1]), j + 1) * 2; + if (!((boost::math::isfinite)(rescale))) + rescale = tools::max_value(); for (int k = j; k < cache_size; ++k) - bessel_cache[k] /= scale; + bessel_cache[k] /= rescale; bessel_j_backwards_iterator ti(b_minus_1_plus_n + j, 2 * sqrt(bessel_arg), bessel_cache[j + 1], bessel_cache[j]); i = ti; } @@ -259,11 +259,11 @@ // if ((j < cache_size - 2) && (tools::max_value() / fabs(64 * bessel_cache[j] / bessel_cache[j + 1]) < fabs(bessel_cache[j]))) { - T scale = pow(fabs(bessel_cache[j] / bessel_cache[j + 1]), j + 1) * 2; - if (!((boost::math::isfinite)(scale))) - scale = tools::max_value(); + T rescale = pow(fabs(bessel_cache[j] / bessel_cache[j + 1]), j + 1) * 2; + if (!((boost::math::isfinite)(rescale))) + rescale = tools::max_value(); for (int k = j; k < cache_size; ++k) - bessel_cache[k] /= scale; + bessel_cache[k] /= rescale; i = bessel_i_backwards_iterator(b_minus_1_plus_n + j, 2 * sqrt(-bessel_arg), bessel_cache[j + 1], bessel_cache[j]); } } @@ -417,7 +417,7 @@ // Overflow: scale = itrunc(tools::log_max_value()) - 10; log_scale += scale; - result /= exp(scale); + result /= exp(T(scale)); } result *= prefix; } diff --git a/include/boost/math/special_functions/detail/hypergeometric_1F1_large_abz.hpp b/include/boost/math/special_functions/detail/hypergeometric_1F1_large_abz.hpp index bfdb4c7a3..5cfb89d20 100644 --- a/include/boost/math/special_functions/detail/hypergeometric_1F1_large_abz.hpp +++ b/include/boost/math/special_functions/detail/hypergeometric_1F1_large_abz.hpp @@ -264,7 +264,7 @@ { // Ooops, need to rescale h: int rescale = itrunc(log(fabs(h))); - T scale = exp(-rescale); + T scale = exp(T(-rescale)); h *= scale; log_scaling += rescale; } diff --git a/include/boost/math/special_functions/detail/hypergeometric_1F1_negative_b_regions.hpp b/include/boost/math/special_functions/detail/hypergeometric_1F1_negative_b_regions.hpp index 02086296c..f63de9c73 100644 --- a/include/boost/math/special_functions/detail/hypergeometric_1F1_negative_b_regions.hpp +++ b/include/boost/math/special_functions/detail/hypergeometric_1F1_negative_b_regions.hpp @@ -435,7 +435,7 @@ namespace boost { double b1 = domain[index][1]; double z0 = domain[index - 1][3]; double z1 = domain[index][3]; - T upper_z_limit = (z0 * (b1 - b) + z1 * (b - b0)) / (b1 - b0); + T upper_z_limit = static_cast((z0 * (b1 - b) + z1 * (b - b0)) / (b1 - b0)); if (z > upper_z_limit) return 1; return z < -b / (4 - 5 * sqrt(log(a)) * a / b) ? -1 : 0; @@ -483,10 +483,10 @@ namespace boost { // of making our lower bound smaller than it would otherwise be, and more so nearer // the centre of the bounding box. // - T effective_x = a + (std::min)(T(a - x1), T(x2 - a)) * 0.25; - T effective_y = b + (std::min)(T(b - y1), T(y2 - b)) * 0.25; + T effective_x = static_cast(a + (std::min)(T(a - x1), T(x2 - a)) * 0.25); + T effective_y = static_cast(b + (std::min)(T(b - y1), T(y2 - b)) * 0.25); - T lower_limit = 1 / ((x2 - x1) * (y2 - y1)) * (f11 * (x2 - effective_x) * (y2 - effective_y) + f21 * (effective_x - x1) * (y2 - effective_y) + f12 * (x2 - effective_x) * (effective_y - y1) + f22 * (effective_x - x1) * (effective_y - y1)); + T lower_limit = static_cast(1 / ((x2 - x1) * (y2 - y1)) * (f11 * (x2 - effective_x) * (y2 - effective_y) + f21 * (effective_x - x1) * (y2 - effective_y) + f12 * (x2 - effective_x) * (effective_y - y1) + f22 * (effective_x - x1) * (effective_y - y1))); // // If one f value was zero, take the whole box as zero: @@ -509,7 +509,7 @@ namespace boost { BOOST_ASSERT(f22 <= domain[index][3]); f22 = domain[index][3]; - T upper_limit = 1 / ((x2 - x1) * (y2 - y1)) * (f11 * (x2 - a) * (y2 - b) + f21 * (a - x1) * (y2 - b) + f12 * (x2 - a) * (b - y1) + f22 * (a - x1) * (b - y1)); + T upper_limit = static_cast(1 / ((x2 - x1) * (y2 - y1)) * (f11 * (x2 - a) * (y2 - b) + f21 * (a - x1) * (y2 - b) + f12 * (x2 - a) * (b - y1) + f22 * (a - x1) * (b - y1))); if (z > upper_limit) return 1; return 0; diff --git a/include/boost/math/special_functions/detail/hypergeometric_1F1_recurrence.hpp b/include/boost/math/special_functions/detail/hypergeometric_1F1_recurrence.hpp index fbfe4fd6c..f219f46f1 100644 --- a/include/boost/math/special_functions/detail/hypergeometric_1F1_recurrence.hpp +++ b/include/boost/math/special_functions/detail/hypergeometric_1F1_recurrence.hpp @@ -417,7 +417,7 @@ // second = boost::math::detail::hypergeometric_1F1_imp(T(a + a_shift - 1), T(b + b_shift), z, pol, scale2); if (scale1 != scale2) - second *= exp(scale2 - scale1); + second *= exp(T(scale2 - scale1)); log_scaling += scale1; // diff --git a/include/boost/math/special_functions/detail/hypergeometric_pFq_checked_series.hpp b/include/boost/math/special_functions/detail/hypergeometric_pFq_checked_series.hpp index 7eb8dcde7..da1410462 100644 --- a/include/boost/math/special_functions/detail/hypergeometric_pFq_checked_series.hpp +++ b/include/boost/math/special_functions/detail/hypergeometric_pFq_checked_series.hpp @@ -129,7 +129,7 @@ Real lower_limit(1 / upper_limit); int log_scaling_factor = itrunc(boost::math::tools::log_max_value()) - 2; Real scaling_factor = exp(Real(log_scaling_factor)); - Real term_m1 = 0; + Real term_m1; int local_scaling = 0; if ((aj.size() == 1) && (bj.size() == 0)) diff --git a/include/boost/math/special_functions/detail/hypergeometric_pade.hpp b/include/boost/math/special_functions/detail/hypergeometric_pade.hpp index 49bc18289..15f86a898 100644 --- a/include/boost/math/special_functions/detail/hypergeometric_pade.hpp +++ b/include/boost/math/special_functions/detail/hypergeometric_pade.hpp @@ -37,7 +37,7 @@ const unsigned max_iterations = boost::math::policies::get_max_series_iterations(); T b2 = T(0), a2 = T(0); - T result = T(0), prev_result = a1 / b1; + T result = T(0), prev_result; for (unsigned k = 1; k < max_iterations; ++k) { diff --git a/include/boost/math/special_functions/detail/hypergeometric_series.hpp b/include/boost/math/special_functions/detail/hypergeometric_series.hpp index 0bf3c6444..f788cdd77 100644 --- a/include/boost/math/special_functions/detail/hypergeometric_series.hpp +++ b/include/boost/math/special_functions/detail/hypergeometric_series.hpp @@ -337,7 +337,6 @@ // Backup state: // term = saved_term * exp(T(local_scaling - saved_scale)); - term_m1 = term; n = summit_location; term *= (b + (n - 1)) * n / ((a + (n - 1)) * z); --n; @@ -389,14 +388,12 @@ sum /= scaling_factor; term /= scaling_factor; log_scaling += log_scaling_factor; - local_scaling += log_scaling_factor; } if (fabs(sum) < lower_limit) { sum *= scaling_factor; term *= scaling_factor; log_scaling -= log_scaling_factor; - local_scaling -= log_scaling_factor; } //term_m1 = term; term *= (((a + n) / ((b + n) * (n + 1))) * z); diff --git a/include/boost/math/special_functions/hypergeometric_1F1.hpp b/include/boost/math/special_functions/hypergeometric_1F1.hpp index 5879ba2d0..6f994b5e3 100644 --- a/include/boost/math/special_functions/hypergeometric_1F1.hpp +++ b/include/boost/math/special_functions/hypergeometric_1F1.hpp @@ -634,8 +634,8 @@ namespace boost { namespace math { namespace detail { // // Actual result will be result * e^log_scaling / tgamma(b). // - int sign = 1; - T scale = log_scaling - boost::math::lgamma(b, &sign, pol); + int result_sign = 1; + T scale = log_scaling - boost::math::lgamma(b, &result_sign, pol); #ifndef BOOST_NO_CXX11_THREAD_LOCAL static const thread_local T max_scaling = boost::math::tools::log_max_value() - 2; static const thread_local T max_scale_factor = exp(max_scaling); @@ -656,7 +656,7 @@ namespace boost { namespace math { namespace detail { } if (scale != 0) result *= exp(scale); - return result * sign; + return result * result_sign; } } // namespace detail diff --git a/include/boost/math/tools/recurrence.hpp b/include/boost/math/tools/recurrence.hpp index 5e9275c13..b090ae7a0 100644 --- a/include/boost/math/tools/recurrence.hpp +++ b/include/boost/math/tools/recurrence.hpp @@ -145,7 +145,7 @@ namespace boost { using boost::math::tuple; using boost::math::get; - T third = 0; + T third; T a, b, c; for (unsigned k = 0; k < number_of_steps; ++k) @@ -200,7 +200,7 @@ namespace boost { using boost::math::tuple; using boost::math::get; - T next = 0; + T next; T a, b, c; for (unsigned k = 0; k < number_of_steps; ++k) diff --git a/include/boost/math/tools/roots.hpp b/include/boost/math/tools/roots.hpp index 552f7f7ae..50b204585 100644 --- a/include/boost/math/tools/roots.hpp +++ b/include/boost/math/tools/roots.hpp @@ -191,7 +191,6 @@ std::pair bisect(F f, T min, T max, Tol tol, boost::uintmax_t& max_iter, c else if(sign(fmid) * sign(fmin) < 0) { max = mid; - fmax = fmid; } else { @@ -280,7 +279,7 @@ T newton_raphson_iterate(F f, T guess, T min, T max, int digits, boost::uintmax_ T shift = (delta > 0) ? (result - min) / 2 : (result - max) / 2; if ((result != 0) && (fabs(shift) > fabs(result))) { - delta = sign(delta) * fabs(result) * 0.9; // Protect against huge jumps! + delta = sign(delta) * fabs(result) * 0.9f; // Protect against huge jumps! //delta = sign(delta) * result; // Protect against huge jumps! Failed for negative result. https://github.com/boostorg/math/issues/216 } else @@ -556,7 +555,7 @@ namespace detail{ // last two steps haven't converged. delta = (delta > 0) ? (result - min) / 2 : (result - max) / 2; if ((result != 0) && (fabs(delta) > result)) - delta = sign(delta) * fabs(result) * 0.9; // protect against huge jumps! + delta = sign(delta) * fabs(result) * 0.9f; // protect against huge jumps! // reset delta2 so that this branch will *not* be taken on the // next iteration: delta2 = delta * 3;