diff --git a/include/boost/math/special_functions/acosh.hpp b/include/boost/math/special_functions/acosh.hpp
index c601960fc..63e8e2dff 100644
--- a/include/boost/math/special_functions/acosh.hpp
+++ b/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>())
             {
diff --git a/include/boost/math/special_functions/asinh.hpp b/include/boost/math/special_functions/asinh.hpp
index d66551e49..52362dd29 100644
--- a/include/boost/math/special_functions/asinh.hpp
+++ b/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>())
             {
diff --git a/include/boost/math/special_functions/bessel.hpp b/include/boost/math/special_functions/bessel.hpp
index 5fe518a36..12e5884f4 100644
--- a/include/boost/math/special_functions/bessel.hpp
+++ b/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:
    //
@@ -190,9 +191,7 @@ 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>("boost::math::cyl_bessel_i<%1%>(%1%,%1%)", "Got x = %1%, but we need x >= 0", x, pol);
    }
    if(x == 0)
    {
@@ -221,7 +220,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;
 }
@@ -240,11 +239,9 @@ inline T cyl_bessel_k_imp(T v, T x, const bessel_no_int_tag& /* t */, const Poli
    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;
 }
@@ -278,11 +275,9 @@ inline T cyl_neumann_imp(T v, T x, const bessel_no_int_tag&, const Policy& pol)
    {
       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_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,
@@ -290,7 +285,7 @@ inline T cyl_neumann_imp(T v, T x, const bessel_no_int_tag&, const Policy& pol)
    // 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,9 +324,7 @@ 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);
@@ -339,7 +332,7 @@ inline T sph_neumann_imp(unsigned v, T x, const Policy& 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 +410,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;
diff --git a/test/test_bernoulli_constants.cpp b/test/test_bernoulli_constants.cpp
index 03ddc8c05..6d73d82a4 100644
--- a/test/test_bernoulli_constants.cpp
+++ b/test/test_bernoulli_constants.cpp
@@ -199,6 +199,11 @@ void test(const char* name)
    BOOST_MATH_CHECK_THROW(boost::math::bernoulli_b2n<T>(overflow_index, boost::math::policies::make_policy(boost::math::policies::overflow_error<boost::math::policies::throw_on_error>())), std::overflow_error);
    BOOST_MATH_CHECK_THROW(boost::math::tangent_t2n<T>(overflow_index, boost::math::policies::make_policy(boost::math::policies::overflow_error<boost::math::policies::throw_on_error>())), std::overflow_error);
 #endif
+   BOOST_MATH_CHECK_THROW(boost::math::bernoulli_b2n<T>(-1), std::domain_error);
+   BOOST_MATH_CHECK_THROW(boost::math::tangent_t2n<T>(-1), std::domain_error);
+   std::vector<T> v;
+   BOOST_MATH_CHECK_THROW(boost::math::bernoulli_b2n<T>(-1, 5, std::back_inserter(v)), std::domain_error);
+   BOOST_MATH_CHECK_THROW(boost::math::tangent_t2n<T>(-1, 5, std::back_inserter(v)), std::domain_error);
 }
 
 void test_real_concept_extra()
diff --git a/test/test_bessel_airy_zeros.cpp b/test/test_bessel_airy_zeros.cpp
index 69132053d..2e5fbe937 100644
--- a/test/test_bessel_airy_zeros.cpp
+++ b/test/test_bessel_airy_zeros.cpp
@@ -116,8 +116,7 @@ void test_bessel_zeros(RealType)
    using boost::math::isnan;
 
   BOOST_MATH_CHECK_THROW(cyl_bessel_j_zero(static_cast<RealType>(0), 0), std::domain_error);
-  // BOOST_MATH_CHECK_THROW(cyl_bessel_j_zero(static_cast<RealType>(-1), 2), std::domain_error);
-  // From 83051 negative orders are supported.
+  BOOST_MATH_CHECK_THROW(cyl_bessel_j_zero(static_cast<RealType>(-1.5), 0), std::domain_error);
 
   // Abuse with infinity and max.
   if (std::numeric_limits<RealType>::has_infinity)
@@ -493,7 +492,23 @@ n |
     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-3), 4), static_cast<RealType>(14.623077742393873174076722507725200649352970569915L), tolerance);
     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-3), 5), static_cast<RealType>(17.818455232945520262553239064736739443380352162752L), tolerance);
 
-  { // Repeat rest using multiple zeros version.
+    /*
+      Table[N[BesselYZero[-5/2, n], 50], {n, 1, 5, 1}]
+      n | y_(-2.5000000000000000000000000000000000000000000000000, n)
+      1 | 5.7634591968945497914064666539527350764090876841674
+      2 | 9.0950113304763551563376983279896952524009293663831
+      3 | 12.322940970566582051969567925329726061189423834915
+      4 | 5.7634591968945497914064666539527350764090876841674
+      5 | 9.0950113304763551563376983279896952524009293663831
+    */
+
+    BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-2.5), 1), static_cast<RealType>(5.7634591968945497914064666539527350764090876841674L), tolerance);
+    BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-2.5), 2), static_cast<RealType>(9.0950113304763551563376983279896952524009293663831L), tolerance);
+    BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-2.5), 3), static_cast<RealType>(12.322940970566582051969567925329726061189423834915L), tolerance);
+    //BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-2.5), 4), static_cast<RealType>(5.7634591968945497914064666539527350764090876841674L), tolerance);
+    //BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-2.5), 5), static_cast<RealType>(9.0950113304763551563376983279896952524009293663831L), tolerance);
+    
+    { // Repeat rest using multiple zeros version.
     std::vector<RealType> zeros;
     cyl_neumann_zero(static_cast<RealType>(0.0), 1, 3, std::back_inserter(zeros) );
     BOOST_CHECK_CLOSE_FRACTION(zeros[0], static_cast<RealType>(0.89357696627916752158488710205833824122514686193001L), tolerance);
diff --git a/test/test_bessel_i.hpp b/test/test_bessel_i.hpp
index 94df9ea73..ed565a27b 100644
--- a/test/test_bessel_i.hpp
+++ b/test/test_bessel_i.hpp
@@ -176,5 +176,15 @@ void test_bessel(T, const char* name)
 
     if(0 != static_cast<T>(ldexp(0.5, -700)))
       do_test_cyl_bessel_i<T>(iv_large_data, name, "Bessel Iv: Mathworld Data (large values)");
+
+    //
+    // Special cases for full coverage:
+    //
+    BOOST_CHECK_THROW(boost::math::cyl_bessel_i(T(-2.5), T(-2.5)), std::domain_error);
+    T tolerance = boost::math::tools::epsilon<T>() * 100;
+    if ((boost::math::tools::digits<T>() <= std::numeric_limits<double>::digits) && (std::numeric_limits<T>::max_exponent > 1000))
+    {
+       BOOST_CHECK_CLOSE_FRACTION(boost::math::cyl_bessel_i(T(0.5), T(710)), static_cast<T>(3.3447452278080108123142599104927325061327359278058601201179e306L), tolerance);
+    }
 }
 
diff --git a/test/test_bessel_j.hpp b/test/test_bessel_j.hpp
index d2b1a937d..2ed1d56c4 100644
--- a/test/test_bessel_j.hpp
+++ b/test/test_bessel_j.hpp
@@ -270,5 +270,11 @@ void test_bessel(T, const char* name)
     BOOST_MATH_CHECK_THROW(boost::math::cyl_bessel_j(T(-2.5), T(0)), std::domain_error);
     BOOST_MATH_CHECK_THROW(boost::math::cyl_bessel_j(T(-2.5), T(-2)), std::domain_error);
     BOOST_MATH_CHECK_THROW(boost::math::cyl_bessel_j(T(2.5), T(-2)), std::domain_error);
+
+    //
+    // special cases for code coverage:
+    //
+    BOOST_CHECK_EQUAL(boost::math::sph_bessel(200, T(0.5)), T(3.070403008048099934928128420285169174541102108657574230431e-497L));
+    BOOST_MATH_CHECK_THROW(boost::math::sph_bessel(2, T(-2.0)), std::domain_error);
 }
 
diff --git a/test/test_bessel_k.hpp b/test/test_bessel_k.hpp
index 09c631a00..a4d8fb3d3 100644
--- a/test/test_bessel_k.hpp
+++ b/test/test_bessel_k.hpp
@@ -171,5 +171,15 @@ void test_bessel(T, const char* name)
     do_test_cyl_bessel_k<T>(bessel_k_int_data, name, "Bessel Kn: Random Data");
 #include "bessel_k_data.ipp"
     do_test_cyl_bessel_k<T>(bessel_k_data, name, "Bessel Kv: Random Data");
+
+    //
+    // Extra test coverage:
+    //
+    BOOST_CHECK_THROW(boost::math::cyl_bessel_k(T(2), T(-1)), std::domain_error);
+    if (std::numeric_limits<T>::has_infinity)
+    {
+       BOOST_CHECK_EQUAL(boost::math::cyl_bessel_k(T(0), T(0)), std::numeric_limits<T>::infinity());
+    }
+    BOOST_CHECK_THROW(boost::math::cyl_bessel_k(T(1.25), T(0)), std::domain_error);
 }
 
diff --git a/test/test_bessel_y.hpp b/test/test_bessel_y.hpp
index 1618171e7..a2d111a0f 100644
--- a/test/test_bessel_y.hpp
+++ b/test/test_bessel_y.hpp
@@ -214,5 +214,22 @@ void test_bessel(T, const char* name)
 
 #include "sph_neumann_data.ipp"
     do_test_sph_neumann_y<T>(sph_neumann_data, name, "y: Random Data");
+
+    //
+    // Additional test coverage:
+    //
+    if (std::numeric_limits<T>::has_infinity)
+    {
+       BOOST_CHECK_EQUAL(boost::math::cyl_neumann(T(0), T(0)), -std::numeric_limits<T>::infinity());
+       BOOST_CHECK_EQUAL(boost::math::sph_neumann(2, (std::numeric_limits<T>::min)() * 1.5f), -std::numeric_limits<T>::infinity());
+       T small = 5.69289e-1645L;
+       if (small != 0)
+       {
+          BOOST_CHECK_EQUAL(boost::math::sph_neumann(2, small), -std::numeric_limits<T>::infinity());
+       }
+    }
+    BOOST_CHECK_THROW(boost::math::cyl_neumann(T(0), T(-1)), std::domain_error);
+    BOOST_CHECK_THROW(boost::math::cyl_neumann(T(2), T(0)), std::domain_error);
+    BOOST_CHECK_THROW(boost::math::sph_neumann(2, T(-2)), std::domain_error);
 }