From b0096e649f74f4ce26d9233633f6772e699d7ff5 Mon Sep 17 00:00:00 2001 From: mckib2 Date: Mon, 18 Jul 2022 19:05:44 -0700 Subject: [PATCH] Avoid overflow in intermediate bessel_ik computation --- .../boost/math/special_functions/detail/bessel_ik.hpp | 10 +++++++++- 1 file changed, 9 insertions(+), 1 deletion(-) diff --git a/include/boost/math/special_functions/detail/bessel_ik.hpp b/include/boost/math/special_functions/detail/bessel_ik.hpp index 0d61835b2..f6d222db1 100644 --- a/include/boost/math/special_functions/detail/bessel_ik.hpp +++ b/include/boost/math/special_functions/detail/bessel_ik.hpp @@ -377,7 +377,15 @@ int bessel_ik(T v, T x, T* I, T* K, int kind, const Policy& pol) for (k = 1; k <= n; k++) // forward recurrence for K { T fact = 2 * (u + k) / x; - if((tools::max_value() - fabs(prev)) / fact < fabs(current)) + // Check for overflow: if (max - |prev|) / fact > max, then overflow + // (max - |prev|) / fact > max + // max * (1 - fact) > |prev| + // if fact < 1: safe to compute overflow check + // if fact >= 1: won't overflow + const bool will_overflow = (fact < 1) + ? tools::max_value() * (1 - fact) > fabs(prev) + : false; + if(!will_overflow && ((tools::max_value() - fabs(prev)) / fact < fabs(current))) { prev /= current; scale /= current;