Improve test coverage.

.drone.star

@@ -15,8 +15,8 @@ windowsglobalimage="cppalliance/dronevs2019"
 
 def main(ctx):
 
-  things_to_test = [ "special_fun", "distribution_tests", "mp", "misc", "interpolators", "quadrature", "autodiff", "long-running-tests", "float128_tests" ]
-  gcc13_things_to_test = [ "special_fun", "distribution_tests", "mp", "misc", "interpolators", "quadrature", "autodiff", "long-running-tests", "float128_tests", "new_floats" ]
+  things_to_test = [ "special_fun", "distribution_tests", "mp", "misc", "interpolators", "quadrature", "autodiff", "long-running-tests", "float128_tests", "concepts" ]
+  gcc13_things_to_test = [ "special_fun", "distribution_tests", "mp", "misc", "interpolators", "quadrature", "autodiff", "long-running-tests", "float128_tests", "concepts", "new_floats" ]
   sanitizer_test = [ "special_fun", "distribution_tests", "misc", "interpolators", "quadrature", "float128_tests" ]
   gnu_5_stds = [ "gnu++14", "c++14" ]
   gnu_6_stds = [ "gnu++14", "c++14", "gnu++17", "c++17" ]

include/boost/math/quadrature/gauss.hpp

@@ -798,8 +798,8 @@ class gauss : public detail::gauss_detail<Real, N, detail::gauss_constant_catego
       Real L1 = abs(result);
       for (unsigned i = non_zero_start; i < base::abscissa().size(); ++i)
       {
-         K fp = f(base::abscissa()[i]);
-         K fm = f(-base::abscissa()[i]);
+         K fp = f(static_cast<Real>(base::abscissa()[i]));
+         K fm = f(static_cast<Real>(-base::abscissa()[i]));
          result += (fp + fm) * static_cast<Real>(base::weights()[i]);
          L1 += (abs(fp) + abs(fm)) * static_cast<Real>(base::weights()[i]);
       }

include/boost/math/quadrature/gauss_kronrod.hpp

@@ -1159,16 +1159,16 @@ class gauss_kronrod : public detail::gauss_kronrod_detail<Real, N, detail::gauss
       Real L1 = abs(kronrod_result);
       for (unsigned i = gauss_start; i < base::abscissa().size(); i += 2)
       {
-         fp = f(base::abscissa()[i]);
-         fm = f(-base::abscissa()[i]);
+         fp = f(static_cast<Real>(base::abscissa()[i]));
+         fm = f(static_cast<Real>(-base::abscissa()[i]));
          kronrod_result += (fp + fm) * static_cast<Real>(base::weights()[i]);
          L1 += (abs(fp) + abs(fm)) * static_cast<Real>(base::weights()[i]);
          gauss_result += (fp + fm) * static_cast<Real>(gauss<Real, (N - 1) / 2>::weights()[i / 2]);
       }
       for (unsigned i = kronrod_start; i < base::abscissa().size(); i += 2)
       {
-         fp = f(base::abscissa()[i]);
-         fm = f(-base::abscissa()[i]);
+         fp = f(static_cast<Real>(base::abscissa()[i]));
+         fm = f(static_cast<Real>(-base::abscissa()[i]));
          kronrod_result += (fp + fm) * static_cast<Real>(base::weights()[i]);
          L1 += (abs(fp) + abs(fm)) * static_cast<Real>(base::weights()[i]);
       }

include/boost/math/special_functions/acosh.hpp

@@ -38,9 +38,7 @@ namespace boost
 
             if((x < 1) || (boost::math::isnan)(x))
             {
-               return policies::raise_domain_error<T>(
-                  "boost::math::acosh<%1%>(%1%)",
-                  "acosh requires x >= 1, but got x = %1%.", x, pol);
+               return policies::raise_domain_error<T>("boost::math::acosh<%1%>(%1%)", "acosh requires x >= 1, but got x = %1%.", x, pol);
             }
             else if    ((x - 1) >= tools::root_epsilon<T>())
             {

include/boost/math/special_functions/airy.hpp

@@ -186,7 +186,7 @@ T airy_ai_zero_imp(int m, const Policy& pol)
    std::uintmax_t iterations_used = iterations_allowed;
 
    // Use a dynamic tolerance because the roots get closer the higher m gets.
-   T tolerance;
+   T tolerance;  // LCOV_EXCL_LINE
 
    if     (m <=   10) { tolerance = T(0.3F); }
    else if(m <=  100) { tolerance = T(0.1F); }
@@ -238,7 +238,7 @@ T airy_bi_zero_imp(int m, const Policy& pol)
    std::uintmax_t iterations_used = iterations_allowed;
 
    // Use a dynamic tolerance because the roots get closer the higher m gets.
-   T tolerance;
+   T tolerance; // LCOV_EXCL_LINE
 
    if     (m <=   10) { tolerance = T(0.3F); }
    else if(m <=  100) { tolerance = T(0.1F); }

include/boost/math/special_functions/asinh.hpp

@@ -38,9 +38,7 @@ namespace boost
 
             if((boost::math::isnan)(x))
             {
-               return policies::raise_domain_error<T>(
-                  "boost::math::asinh<%1%>(%1%)",
-                  "asinh requires a finite argument, but got x = %1%.", x, pol);
+               return policies::raise_domain_error<T>("boost::math::asinh<%1%>(%1%)", "asinh requires a finite argument, but got x = %1%.", x, pol);
             }
             if        (x >= tools::forth_root_epsilon<T>())
             {

include/boost/math/special_functions/atanh.hpp

@@ -40,21 +40,15 @@ namespace boost
 
             if(x < -1)
             {
-               return policies::raise_domain_error<T>(
-                  function,
-                  "atanh requires x >= -1, but got x = %1%.", x, pol);
+               return policies::raise_domain_error<T>(function, "atanh requires x >= -1, but got x = %1%.", x, pol);
             }
             else if(x > 1)
             {
-               return policies::raise_domain_error<T>(
-                  function,
-                  "atanh requires x <= 1, but got x = %1%.", x, pol);
+               return policies::raise_domain_error<T>(function, "atanh requires x <= 1, but got x = %1%.", x, pol);
             }
             else if((boost::math::isnan)(x))
             {
-               return policies::raise_domain_error<T>(
-                  function,
-                  "atanh requires -1 <= x <= 1, but got x = %1%.", x, pol);
+               return policies::raise_domain_error<T>(function, "atanh requires -1 <= x <= 1, but got x = %1%.", x, pol);
             }
             else if(x < -1 + tools::epsilon<T>())
             {

include/boost/math/special_functions/bernoulli.hpp

@@ -88,7 +88,7 @@ inline OutputIterator bernoulli_b2n(const int start_index,
    if(start_index < 0)
    {
       *out_it = policies::raise_domain_error<T>("boost::math::bernoulli_b2n<%1%>", "Index should be >= 0 but got %1%", T(start_index), pol);
-      return ++out_it;
+      return ++out_it; // LCOV_EXCL_LINE we don't reach here, previous line throws.
    }
 
    return boost::math::detail::bernoulli_number_imp<T>(out_it, start_index, number_of_bernoullis_b2n, pol, tag_type());
@@ -130,7 +130,7 @@ inline OutputIterator tangent_t2n(const int start_index,
    if(start_index < 0)
    {
       *out_it = policies::raise_domain_error<T>("boost::math::tangent_t2n<%1%>", "Index should be >= 0 but got %1%", T(start_index), pol);
-      return ++out_it;
+      return ++out_it; // LCOV_EXCL_LINE we don't reach here, previous line throws.
    }
 
    return boost::math::detail::get_bernoulli_numbers_cache<T, Policy>().copy_tangent_numbers(out_it, start_index, number_of_tangent_t2n, pol);

include/boost/math/special_functions/bessel.hpp

@@ -99,20 +99,23 @@ T cyl_bessel_j_imp(T v, T x, const bessel_no_int_tag& t, const Policy& pol)
    if(x < 0)
    {
       // better have integer v:
-      if(floor(v) == v)
+      if (floor(v) == v)
       {
+         // LCOV_EXCL_START
+         // This branch is hit by multiprecision types only, and is
+         // tested by our real_concept tests, but thee are excluded from coverage
+         // due to time constraints.
          T r = cyl_bessel_j_imp(v, T(-x), t, pol);
-         if(iround(v, pol) & 1)
+         if (iround(v, pol) & 1)
             r = -r;
          return r;
+         // LCOV_EXCL_STOP
       }
       else
-         return policies::raise_domain_error<T>(
-            function,
-            "Got x = %1%, but we need x >= 0", x, pol);
+         return policies::raise_domain_error<T>(function, "Got x = %1%, but we need x >= 0", x, pol);
    }
 
-   T result_J, y;
+   T result_J, y; // LCOV_EXCL_LINE
    bessel_jy(v, x, &result_J, &y, need_j, pol);
    return result_J;
 }
@@ -143,9 +146,7 @@ inline T sph_bessel_j_imp(unsigned n, T x, const Policy& pol)
 {
    BOOST_MATH_STD_USING // ADL of std names
    if(x < 0)
-      return policies::raise_domain_error<T>(
-         "boost::math::sph_bessel_j<%1%>(%1%,%1%)",
-         "Got x = %1%, but function requires x > 0.", x, pol);
+      return policies::raise_domain_error<T>("boost::math::sph_bessel_j<%1%>(%1%,%1%)", "Got x = %1%, but function requires x > 0.", x, pol);
    //
    // Special case, n == 0 resolves down to the sinus cardinal of x:
    //
@@ -179,6 +180,7 @@ T cyl_bessel_i_imp(T v, T x, const Policy& pol)
    // case has better error handling too).
    //
    BOOST_MATH_STD_USING
+   static const char* function = "boost::math::cyl_bessel_i<%1%>(%1%,%1%)";
    if(x < 0)
    {
       // better have integer v:
@@ -190,12 +192,12 @@ T cyl_bessel_i_imp(T v, T x, const Policy& pol)
          return r;
       }
       else
-         return policies::raise_domain_error<T>(
-         "boost::math::cyl_bessel_i<%1%>(%1%,%1%)",
-            "Got x = %1%, but we need x >= 0", x, pol);
+         return policies::raise_domain_error<T>(function, "Got x = %1%, but we need x >= 0", x, pol);
    }
    if(x == 0)
    {
+      if(v < 0) 
+         return floor(v) == v ? static_cast<T>(0) : policies::raise_overflow_error<T>(function, nullptr, pol);
       return (v == 0) ? static_cast<T>(1) : static_cast<T>(0);
    }
    if(v == 0.5f)
@@ -221,7 +223,7 @@ T cyl_bessel_i_imp(T v, T x, const Policy& pol)
    }
    if((v > 0) && (x / v < 0.25))
       return bessel_i_small_z_series(v, x, pol);
-   T result_I, result_K;
+   T result_I, result_K; // LCOV_EXCL_LINE
    bessel_ik(v, x, &result_I, &result_K, need_i, pol);
    return result_I;
 }
@@ -233,18 +235,14 @@ inline T cyl_bessel_k_imp(T v, T x, const bessel_no_int_tag& /* t */, const Poli
    BOOST_MATH_STD_USING
    if(x < 0)
    {
-      return policies::raise_domain_error<T>(
-         function,
-         "Got x = %1%, but we need x > 0", x, pol);
+      return policies::raise_domain_error<T>(function, "Got x = %1%, but we need x > 0", x, pol);
    }
    if(x == 0)
    {
       return (v == 0) ? policies::raise_overflow_error<T>(function, nullptr, pol)
-         : policies::raise_domain_error<T>(
-         function,
-         "Got x = %1%, but we need x > 0", x, pol);
+         : policies::raise_domain_error<T>(function, "Got x = %1%, but we need x > 0", x, pol);
    }
-   T result_I, result_K;
+   T result_I, result_K; // LCOV_EXCL_LINE
    bessel_ik(v, x, &result_I, &result_K, need_k, pol);
    return result_K;
 }
@@ -277,20 +275,18 @@ inline T cyl_neumann_imp(T v, T x, const bessel_no_int_tag&, const Policy& pol)
    if(x <= 0)
    {
       return (v == 0) && (x == 0) ?
-         policies::raise_overflow_error<T>(function, nullptr, pol)
-         : policies::raise_domain_error<T>(
-               function,
-               "Got x = %1%, but result is complex for x <= 0", x, pol);
+         -policies::raise_overflow_error<T>(function, nullptr, pol) // LCOV_EXCL_LINE MP case only here, not tested in code coverage as it takes too long.
+         : policies::raise_domain_error<T>(function, "Got x = %1%, but result is complex for x <= 0", x, pol);
    }
-   T result_J, y;
+   T result_J, y; // LCOV_EXCL_LINE
    bessel_jy(v, x, &result_J, &y, need_y, pol);
    //
    // Post evaluation check for internal overflow during evaluation,
    // can occur when x is small and v is large, in which case the result
    // is -INF:
    //
    if(!(boost::math::isfinite)(y))
-      return -policies::raise_overflow_error<T>(function, nullptr, pol);
+      return -policies::raise_overflow_error<T>(function, nullptr, pol); // LCOV_EXCL_LINE defensive programming?  Might be dead code now...
    return y;
 }
 
@@ -329,17 +325,15 @@ inline T sph_neumann_imp(unsigned v, T x, const Policy& pol)
    // evaluate the function's definition directly:
    //
    if(x < 0)
-      return policies::raise_domain_error<T>(
-         function,
-         "Got x = %1%, but function requires x > 0.", x, pol);
+      return policies::raise_domain_error<T>(function, "Got x = %1%, but function requires x > 0.", x, pol);
 
    if(x < 2 * tools::min_value<T>())
       return -policies::raise_overflow_error<T>(function, nullptr, pol);
 
    T result = cyl_neumann_imp(T(T(v)+0.5f), x, bessel_no_int_tag(), pol);
    T tx = sqrt(constants::pi<T>() / (2 * x));
 
-   if((tx > 1) && (tools::max_value<T>() / tx < result))
+   if((tx > 1) && (tools::max_value<T>() / tx < fabs(result)))
       return -policies::raise_overflow_error<T>(function, nullptr, pol);
 
    return result * tx;
@@ -417,8 +411,7 @@ inline T cyl_bessel_j_zero_imp(T v, int m, const Policy& pol)
 
    if(number_of_iterations >= policies::get_max_root_iterations<Policy>())
    {
-      return policies::raise_evaluation_error<T>(function, "Unable to locate root in a reasonable time:"
-         "  Current best guess is %1%", jvm, Policy());
+      return policies::raise_evaluation_error<T>(function, "Unable to locate root in a reasonable time: Current best guess is %1%", jvm, Policy()); // LCOV_EXCL_LINE
    }
 
    return jvm;
@@ -496,8 +489,7 @@ inline T cyl_neumann_zero_imp(T v, int m, const Policy& pol)
 
    if(number_of_iterations >= policies::get_max_root_iterations<Policy>())
    {
-      return policies::raise_evaluation_error<T>(function, "Unable to locate root in a reasonable time:"
-         "  Current best guess is %1%", yvm, Policy());
+      return policies::raise_evaluation_error<T>(function, "Unable to locate root in a reasonable time: Current best guess is %1%", yvm, Policy()); //LCOV_EXCL_LINE
    }
 
    return yvm;

include/boost/math/special_functions/bessel_iterators.hpp

@@ -9,6 +9,7 @@
 #define BOOST_MATH_BESSEL_ITERATORS_HPP
 
 #include <boost/math/tools/recurrence.hpp>
+#include <boost/math/special_functions/bessel.hpp>
 
 namespace boost {
    namespace math {
@@ -154,8 +155,8 @@ namespace boost {
             if (v > 1)
                boost::math::policies::raise_domain_error("bessel_i_forwards_iterator<%1%>", "Order must be < 0 stable forwards recurrence but got %1%", v, Policy());
          }
-         bessel_i_forwards_iterator(const T& v, const T& x, const T& I_v_plus_1, const T& I_v)
-            : it(detail::bessel_ik_recurrence<T>(v, x), I_v_plus_1, I_v)
+         bessel_i_forwards_iterator(const T& v, const T& x, const T& I_v_minus_1, const T& I_v)
+            : it(detail::bessel_ik_recurrence<T>(v, x), I_v_minus_1, I_v)
          {
             if (v > 1)
                boost::math::policies::raise_domain_error("bessel_i_forwards_iterator<%1%>", "Order must be < 0 stable forwards recurrence but got %1%", v, Policy());