libs/math/test/test_bessel_airy_zeros.cpp

   1 //  Copyright John Maddock 2013
   2 //  Copyright Christopher Kormanyos 2013.
   3 //  Copyright Paul A. Bristow 2013.
   4
   5 //  Use, modification and distribution are subject to the
   6 //  Boost Software License, Version 1.0. (See accompanying file
   7 //  LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
   8
   9 #ifdef _MSC_VER
  10 #  pragma warning(disable : 4127) // conditional expression is constant.
  11 #  pragma warning(disable : 4512) // assignment operator could not be generated.
  12 #  pragma warning(disable : 4996) // use -D_SCL_SECURE_NO_WARNINGS.
  13 #endif
  14
  15 //#include <pch_light.hpp> // commented out during testing.
  16
  17 //#include <boost/math/special_functions/math_fwd.hpp>
  18 #include <boost/cstdint.hpp>
  19 #include <boost/math/special_functions/bessel.hpp>
  20 #include <boost/math/special_functions/airy.hpp>
  21 #include <boost/math/tools/test.hpp>
  22
  23 #include <boost/math/concepts/real_concept.hpp> // for real_concept
  24
  25 #define BOOST_TEST_MAIN
  26 #include <boost/test/unit_test.hpp> // Boost.Test
  27 #include <boost/test/tools/floating_point_comparison.hpp>
  28
  29 #include <typeinfo>
  30 #include <iostream>
  31 #include <iomanip>
  32
  33 // #include <boost/math/tools/
  34 //
  35 // DESCRIPTION:
  36 // ~~~~~~~~~~~~
  37 //
  38 // This file tests the functions that evaluate zeros (or roots) of Bessel, Neumann and Airy functions.
  39
  40 // Spot tests which compare our results with selected values computed
  41 // using the online special function calculator at functions.wolfram.com,
  42 // and values generated with Boost.Multiprecision at about 1000-bit or 100 decimal digits precision.
  43
  44 // We are most grateful for the invaluable
  45 // Weisstein, Eric W. "Bessel Function Zeros." From MathWorld--A Wolfram Web Resource.
  46 // http://mathworld.wolfram.com/BesselFunctionZeros.html
  47 // and the newer http://www.wolframalpha.com/
  48
  49 // See also NIST Handbook of Mathematical Function http://dlmf.nist.gov/10.21
  50 /*
  51 Tests of cyl Bessel and cyl Neumann zeros.
  52 ==========================================
  53
  54 The algorithms for estimating the roots of both cyl. Bessel
  55 as well as cyl. Neumann have the same cross-over points,
  56 and also use expansions that have the same order of approximation.
  57
  58 Therefore, tests will be equally effective for both functions in the regions of order.
  59
  60 I have recently changed a critical cross-over in the algorithms
  61 from a value of order of 1.2 to a value of order of 2.2.
  62 In addition, there is a critical cross-over in the rank of the
  63 zero from rank 1 to rank 2 and above. The first zero is
  64 treated differently than the remaining ones.
  65
  66 The test cover various regions of order,
  67 each one tested with several zeros:
  68   * Order 219/100: This checks a region just below a critical cutof.
  69   * Order 221/100: This checks a region just above a critical cutoff.
  70   * Order 0: Something always tends to go wrong at zero.
  71   * Order 1/1000: A small order.
  72   * Order 71/19: Merely an intermediate order.
  73   * Order 7001/19: A medium-large order, small enough to retain moderate efficiency of calculation.
  74
  75 There are also a few selected high zeros
  76 such as the 1000th zero for a few modest orders such as 71/19, etc.
  77
  78 Tests of Airy zeros.
  79 ====================
  80
  81 The Airy zeros algorithms use tabulated values for the first 10 zeros,
  82 whereby algorithms are used for rank 11 and higher.
  83 So testing the zeros of Ai and Bi from 1 through 20 handles
  84 this cross-over nicely.
  85
  86 In addition, the algorithms for the estimates of the zeros
  87 become increasingly accurate for larger, negative argument.
  88
  89 On the other hand, the zeros become increasingly close
  90 for large, negative argument. So another nice test
  91 involves testing pairs of zeros for different orders of
  92 magnitude of the zeros, to insure that the program
  93 properly resolves very closely spaced zeros.
  94 */
  95
  96
  97 template <class RealType>
  98 void test_bessel_zeros(RealType)
  99 {
 100   // Basic sanity checks for finding zeros of Bessel and Airy function.
 101   // where template parameter RealType can be float, double, long double,
 102   // or real_concept, a prototype for user-defined floating-point types.
 103
 104   // Parameter RealType is only used to communicate the RealType, float, double...
 105   // and is an arbitrary zero for all tests.
 106    RealType tolerance = 5 * (std::max)(
 107      static_cast<RealType>(boost::math::tools::epsilon<long double>()),
 108      boost::math::tools::epsilon<RealType>());
 109    std::cout << "Tolerance for type " << typeid(RealType).name()  << " is " << tolerance << "." << std::endl;
 110    //
 111    // An extra fudge factor for real_concept which has a less accurate tgamma:
 112    RealType tolerance_tgamma_extra = std::numeric_limits<RealType>::is_specialized ? 1 : 15;
 113
 114    // http://www.wolframalpha.com/
 115    using boost::math::cyl_bessel_j_zero; // (nu, j)
 116    using boost::math::isnan;
 117
 118   BOOST_MATH_CHECK_THROW(cyl_bessel_j_zero(static_cast<RealType>(0), 0), std::domain_error);
 119   // BOOST_MATH_CHECK_THROW(cyl_bessel_j_zero(static_cast<RealType>(-1), 2), std::domain_error);
 120   // From 83051 negative orders are supported.
 121
 122   // Abuse with infinity and max.
 123   if (std::numeric_limits<RealType>::has_infinity)
 124   {
 125     //BOOST_CHECK_EQUAL(cyl_bessel_j_zero(static_cast<RealType>(std::numeric_limits<RealType>::infinity()), 1),
 126     //  static_cast<RealType>(std::numeric_limits<RealType>::infinity()) );
 127     // unknown location(0): fatal error in "test_main":
 128     // class boost::exception_detail::clone_impl<struct boost::exception_detail::error_info_injector<class std::domain_error> >:
 129     // Error in function boost::math::cyl_bessel_j_zero<double>(double, int): Order argument is 1.#INF, but must be finite >= 0 !
 130     // Note that the reported type long double is not the type of the original call RealType,
 131     // but the promoted value, here long double, if applicable.
 132     BOOST_MATH_CHECK_THROW(cyl_bessel_j_zero(static_cast<RealType>(std::numeric_limits<RealType>::infinity()), 1),
 133      std::domain_error);
 134     BOOST_MATH_CHECK_THROW(cyl_bessel_j_zero(static_cast<RealType>(-std::numeric_limits<RealType>::infinity()), 1),
 135      std::domain_error);
 136
 137   }
 138   // Test with maximum value of v that will cause evaluation error
 139   //BOOST_MATH_CHECK_THROW(cyl_bessel_j_zero(boost::math::tools::max_value<RealType>(), 1), std::domain_error);
 140   // unknown location(0): fatal error in "test_main":
 141   // class boost::exception_detail::clone_impl<struct boost::exception_detail::error_info_injector<class boost::math::evaluation_error> >:
 142   // Error in function boost::math::bessel_jy<double>(double,double): Order of Bessel function is too large to evaluate: got 3.4028234663852886e+038
 143
 144   BOOST_MATH_CHECK_THROW(cyl_bessel_j_zero(boost::math::tools::max_value<RealType>(), 1), boost::math::evaluation_error);
 145
 146   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(boost::math::tools::min_value<RealType>(), 1),
 147     static_cast<RealType>(2.4048255576957727686216318793264546431242449091460L), tolerance);
 148
 149   BOOST_CHECK_CLOSE_FRACTION(-cyl_bessel_j_zero(boost::math::tools::min_value<RealType>(), 1),
 150     static_cast<RealType>(-2.4048255576957727686216318793264546431242449091460L), tolerance);
 151
 152   // Checks on some spot values.
 153
 154   // http://mathworld.wolfram.com/BesselFunctionZeros.html provides some spot values,
 155   // evaluation at 50 deciaml digits using WoldramAlpha.
 156
 157   /* Table[N[BesselJZero[0, n], 50], {n, 1, 5, 1}]
 158   n |
 159   1 | 2.4048255576957727686216318793264546431242449091460
 160   2 | 5.5200781102863106495966041128130274252218654787829
 161   3 | 8.6537279129110122169541987126609466855657952312754
 162   4 | 11.791534439014281613743044911925458922022924699695
 163   5 | 14.930917708487785947762593997388682207915850115633
 164   */
 165   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(0), 1), static_cast<RealType>(2.4048255576957727686216318793264546431242449091460L), tolerance);
 166   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(0), 2), static_cast<RealType>(5.5200781102863106495966041128130274252218654787829L), tolerance);
 167   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(0), 3), static_cast<RealType>(8.6537279129110122169541987126609466855657952312754L), tolerance);
 168   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(0), 4), static_cast<RealType>(11.791534439014281613743044911925458922022924699695L), tolerance);
 169   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(0), 5), static_cast<RealType>(14.930917708487785947762593997388682207915850115633L), tolerance);
 170
 171   { // Same test using the multiple zeros version.
 172     std::vector<RealType> zeros;
 173     cyl_bessel_j_zero(static_cast<RealType>(0.0), 1, 3, std::back_inserter(zeros) );
 174     BOOST_CHECK_CLOSE_FRACTION(zeros[0], static_cast<RealType>(2.4048255576957727686216318793264546431242449091460L), tolerance);
 175     BOOST_CHECK_CLOSE_FRACTION(zeros[1], static_cast<RealType>(5.5200781102863106495966041128130274252218654787829L), tolerance);
 176     BOOST_CHECK_CLOSE_FRACTION(zeros[2], static_cast<RealType>(8.6537279129110122169541987126609466855657952312754L), tolerance);
 177   }
 178   // 1/1000 a small order.
 179   /* Table[N[BesselJZero[1/1000, n], 50], {n, 1, 4, 1}]
 180   n |
 181 1 | 2.4063682720422009275161970278295108254321633626292
 182 2 | 5.5216426858401848664019464270992222126391378706092
 183 3 | 8.6552960859298799453893840513333150237193779482071
 184 4 | 11.793103797689738596231262077785930962647860975357
 185
 186 Table[N[BesselJZero[1/1000, n], 50], {n, 10, 20, 1}]
 187 n |
 188 10 | 30.636177039613574749066837922778438992469950755736
 189 11 | 33.777390823252864715296422192027816488172667994611
 190 12 | 36.918668992567585467000743488690258054442556198147
 191 13 | 40.059996426251227493370316149043896483196561190610
 192 14 | 43.201362392820317233698309483240359167380135262681
 193 15 | 46.342759065846108737848449985452774243376260538634
 194 16 | 49.484180603489984324820981438067325210499739716337
 195 17 | 52.625622557085775090390071484188995092211215108718
 196 18 | 55.767081479279692992978326069855684800673801918763
 197 19 | 58.908554657366270044071505013449016741804538135905
 198 20 | 62.050039927521244984641179233170843941940575857282
 199
 200 */
 201
 202   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(1)/1000, 1), static_cast<RealType>(2.4063682720422009275161970278295108254321633626292L), tolerance);
 203   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(1)/1000, 4), static_cast<RealType>(11.793103797689738596231262077785930962647860975357L), tolerance);
 204   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(1)/1000, 10), static_cast<RealType>(30.636177039613574749066837922778438992469950755736L), tolerance);
 205   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(1)/1000, 20), static_cast<RealType>(62.050039927521244984641179233170843941940575857282L), tolerance);
 206
 207     /*
 208   Table[N[BesselJZero[1, n], 50], {n, 1, 4, 1}]
 209   n |
 210   1 | 3.8317059702075123156144358863081607665645452742878
 211   2 | 7.0155866698156187535370499814765247432763115029142
 212   3 | 10.173468135062722077185711776775844069819512500192
 213   4 | 13.323691936314223032393684126947876751216644731358
 214   */
 215
 216   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(1), 1), static_cast<RealType>(3.8317059702075123156144358863081607665645452742878L), tolerance);
 217   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(1), 2), static_cast<RealType>(7.0155866698156187535370499814765247432763115029142L), tolerance);
 218   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(1), 3), static_cast<RealType>(10.173468135062722077185711776775844069819512500192L), tolerance);
 219   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(1), 4), static_cast<RealType>(13.323691936314223032393684126947876751216644731358L), tolerance);
 220
 221   /*
 222   Table[N[BesselJZero[5, n], 50], {n, 1, 5, 1}]
 223   n |
 224   1 | 8.7714838159599540191228671334095605629810770148974
 225   2 | 12.338604197466943986082097644459004412683491122239
 226   3 | 15.700174079711671037587715595026422501346662246893
 227   4 | 18.980133875179921120770736748466932306588828411497
 228   5 | 22.217799896561267868824764947529187163096116704354
 229 */
 230
 231   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(5), 1), static_cast<RealType>(8.7714838159599540191228671334095605629810770148974L), tolerance);
 232   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(5), 2), static_cast<RealType>(12.338604197466943986082097644459004412683491122239L), tolerance);
 233   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(5), 3), static_cast<RealType>(15.700174079711671037587715595026422501346662246893L), tolerance);
 234   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(5), 4), static_cast<RealType>(18.980133875179921120770736748466932306588828411497L), tolerance);
 235   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(5), 5), static_cast<RealType>(22.217799896561267868824764947529187163096116704354L), tolerance);
 236
 237   // An intermediate order
 238   /*
 239   Table[N[BesselJZero[71/19, n], 50], {n, 1, 20, 1}]
 240
 241   7.27317519383164895031856942622907655889631967016227,
 242   10.7248583088831417325361727458514166471107495990849,
 243   14.0185045994523881061204595580426602824274719315813,
 244   17.2524984591704171821624871665497773491959038386104,
 245   20.4566788740445175951802340838942858854605020778141,
 246   23.6436308971423452249455142271473195998540517250404,
 247   26.8196711402550877454213114709650192615223905192969,
 248   29.9883431174236747426791417966614320438788681941419,
 249   33.1517968976905208712508624699734452654447919661140,
 250   36.3114160002162074157243540350393860813165201842005,
 251   39.4681324675052365879451978080833378877659670320292,
 252   42.6225978013912364748550348312979540188444334802274,
 253   45.7752814645368477533902062078067265814959500124386,
 254   48.9265304891735661983677668174785539924717398947994,
 255   52.0766070453430027942797460418789248768734780634716,
 256   55.2257129449125713935942243278172656890590028901917,
 257   58.3740061015388864367751881504390252017351514189321,
 258   61.5216118730009652737267426593531362663909441035715,
 259   64.6686310537909303683464822148736607945659662871596,
 260   67.8151456196962909255567913755559511651114605854579
 261   */
 262   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(71)/19, 1), static_cast<RealType>(7.27317519383164895031856942622907655889631967016227L), tolerance);
 263   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(71)/19, 4), static_cast<RealType>(17.2524984591704171821624871665497773491959038386104L), tolerance);
 264   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(71)/19, 10), static_cast<RealType>(36.3114160002162074157243540350393860813165201842005L), tolerance);
 265   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(71)/19, 20), static_cast<RealType>(67.8151456196962909255567913755559511651114605854579L), tolerance);
 266   /*
 267
 268   Table[N[BesselJZero[7001/19, n], 50], {n, 1, 2, 1}]
 269
 270 1 | 381.92201523024489386917204470434842699154031135348
 271 2 | 392.17508657648737502651299853099852567001239217724
 272
 273 Table[N[BesselJZero[7001/19, n], 50], {n, 19, 20, 1}]
 274
 275 19 | 491.67809669154347398205298745712766193052308172472
 276 20 | 496.39435037938252557535375498577989720272298310802
 277
 278   */
 279   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(7001)/19, 1), static_cast<RealType>(381.92201523024489386917204470434842699154031135348L), tolerance);
 280   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(7001)/19, 2), static_cast<RealType>(392.17508657648737502651299853099852567001239217724L), tolerance);
 281   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(7001)/19, 20), static_cast<RealType>(496.39435037938252557535375498577989720272298310802L), tolerance);
 282
 283   // Some non-integral tests.
 284   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(3.73684210526315789473684210526315789473684210526315789L), 1), static_cast<RealType>(7.273175193831648950318569426229076558896319670162279791988152000556091140599946365217211157877052381L), tolerance);
 285   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(3.73684210526315789473684210526315789473684210526315789L), 20), static_cast<RealType>(67.81514561969629092555679137555595116511146058545787883557679231060644931096494584364894743334132014L), tolerance);
 286
 287   // Some non-integral tests in 'tough' regions.
 288   // Order 219/100: This checks a region just below a critical cutoff.
 289   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(219)/100, 1), static_cast<RealType>(5.37568854370623186731066365697341253761466705063679L), tolerance);
 290   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(219)/100, 2), static_cast<RealType>(8.67632060963888122764226633146460596009874991130394L), tolerance);
 291   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(219)/100, 20), static_cast<RealType>(65.4517712237598926858973399895944886397152223643028L), tolerance);
 292   // Order 221/100: This checks a region just above a critical cutoff.
 293   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(221)/100, 1), static_cast<RealType>(5.40084731984998184087380740054933778965260387203942L), tolerance);
 294   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(221)/100, 2), static_cast<RealType>(8.70347906513509618445695740167369153761310106851599L), tolerance);
 295   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(221)/100, 20), static_cast<RealType>(65.4825314862621271716158606625527548818843845600782L), tolerance);
 296
 297   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(7001)/19, 1), static_cast<RealType>(381.922015230244893869172044704348426991540311353476L), tolerance);
 298   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(7001)/19, 2), static_cast<RealType>(392.175086576487375026512998530998525670012392177242L), tolerance);
 299   BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(7001)/19, 20), static_cast<RealType>(496.394350379382525575353754985779897202722983108025L), tolerance);
 300
 301   // Zero'th cases.
 302   BOOST_MATH_CHECK_THROW(boost::math::cyl_bessel_j_zero(static_cast<RealType>(0), 0), std::domain_error); // Zero'th zero of J0(x).
 303   BOOST_CHECK(boost::math::cyl_bessel_j_zero(static_cast<RealType>(1), 0) == 0); // Zero'th zero of J1(x).
 304   BOOST_CHECK(boost::math::cyl_bessel_j_zero(static_cast<RealType>(2), 0) == 0); // Zero'th zero of J2(x).
 305
 306
 307   // Negative order cases.
 308   // Table[N[BesselJZero[-39, n], 51], {n, 1, 20, 1}]
 309
 310   //  45.597624026432090522996531982029164361723758769649
 311
 312     BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(-39), 1), static_cast<RealType>(45.597624026432090522996531982029164361723758769649L), tolerance);
 313     BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(-39), 2), static_cast<RealType>(50.930599960211455519691708196247756810739999585797L), tolerance);
 314     BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(-39), 4), static_cast<RealType>(59.810708207036942166964205243063534405954475825070L), tolerance);
 315     BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(-39), 10), static_cast<RealType>(82.490310026657839398140015188318580114553721419436L), tolerance);
 316     BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(-39), 15), static_cast<RealType>(99.886172950858129702511715161572827825877395517083L), tolerance);
 317     BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(-39), 20), static_cast<RealType>(116.73117751356457774415638043701531989536641098359L), tolerance);
 318
 319    // Table[N[BesselJZero[-39 - (1/3), n], 51], {n, 1, 20, 1}]
 320
 321    // 43.803165820025277290601047312311146608776920513241
 322    // 49.624678304306778749502719837270544976331123155017
 323
 324     RealType v = static_cast<RealType>(-39);
 325
 326     v -= boost::math::constants::third<RealType>();
 327
 328    // BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(v, 1), static_cast<RealType>(43.803165820025277290601047312311146608776920513241L), tolerance);
 329   //  BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(-39) - static_cast<RealType>(1)/3, 1), static_cast<RealType>(43.803165820025277290601047312311146608776920513241L), tolerance);
 330     BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(v, 2), static_cast<RealType>(49.624678304306778749502719837270544976331123155017L), tolerance * 4);
 331     BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(-39) - static_cast<RealType>(0.333333333333333333333333333333333333333333333L), 5), static_cast<RealType>(62.911281619408963609400485687996804820400102193455L), tolerance * 4);
 332     BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(-39) - static_cast<RealType>(0.333333333333333333333333333333333333333333333L), 10), static_cast<RealType>(81.705998611506506523381866527389118594062841737382L), tolerance * 4);
 333     BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(static_cast<RealType>(-39) - static_cast<RealType>(0.333333333333333333333333333333333333333333333L), 20), static_cast<RealType>(116.05368337161392034833932554892349580959931408963L), tolerance * 4);
 334
 335
 336   // Table[N[BesselJZero[-1/3, n], 51], {n, 1, 20, 1}]
 337     // 1.86635085887389517154698498025466055044627209492336
 338     // 4.98785323143515872689263163814239463653891121063534
 339     v = - boost::math::constants::third<RealType>();
 340     BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(v, 1), static_cast<RealType>(1.86635085887389517154698498025466055044627209492336L), tolerance);
 341     BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(v, 2), static_cast<RealType>(4.98785323143515872689263163814239463653891121063534L), tolerance);
 342     BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(v, 5), static_cast<RealType>(14.4037758801360172217813556328092353168458341692115L), tolerance);
 343     BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(v, 20), static_cast<RealType>(61.5239847181314647255554392599009248210564008120358L), tolerance);
 344
 345   // Table[N[BesselJZero[-3 - (999999/1000000), n], 51], {n, 1, 20, 1}]
 346     // 0.666908567552422764702292353801313970109968787260547
 347     //7.58834489983121936102504707121493271448122800440112
 348
 349     std::cout.precision(2 + std::numeric_limits<RealType>::digits * 3010/10000);
 350     v = -static_cast<RealType>(3);
 351     //std::cout << "v = " << v << std::endl;
 352     RealType d = static_cast<RealType>(999999)/1000000; // Value very near to unity.
 353     //std::cout << "d = " << d << std::endl;
 354     v -= d;
 355    // std::cout << "v = " << v << std::endl; //   v = -3.9999989999999999
 356
 357     // 1st is much less accurate.
 358     BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(v, 1), static_cast<RealType>(0.666908567552422764702292353801313970109968787260547L), tolerance * 500000);
 359     BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(v, 2), static_cast<RealType>(7.58834489983121936102504707121493271448122800440112L), tolerance);
 360     BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(v, 5), static_cast<RealType>(17.6159678964372778134202353240221384945968807948928L), tolerance);
 361     BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(v, 20), static_cast<RealType>(65.0669968910414433560468307554730940098734494938136L), tolerance);
 362
 363
 364     v = -static_cast<RealType>(1)/81799; // Largish prime, so small value.
 365     // std::cout << "v = " << v << std::endl; // v = -1.22251e-005
 366
 367    // Table[N[BesselJZero[-1/81799, n], 51], {n, 1, 20, 1}]
 368
 369     // 2.40480669570616362235270726259606288441474232101937
 370     //5.52005898213436490056801834487410496538653938730884
 371     BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(v, 1), static_cast<RealType>(2.40480669570616362235270726259606288441474232101937L), tolerance);
 372     BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(v, 2), static_cast<RealType>(5.52005898213436490056801834487410496538653938730884L), tolerance);
 373     BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(v, 5), static_cast<RealType>(14.9308985160466385806685583210609848822943295303368L), tolerance);
 374     BOOST_CHECK_CLOSE_FRACTION(cyl_bessel_j_zero(v, 20), static_cast<RealType>(62.0484499877253314338528593349200129641402661038743L), tolerance);
 375
 376   // Confirm that negative m throws domain_error.
 377   BOOST_MATH_CHECK_THROW(boost::math::cyl_bessel_j_zero(static_cast<RealType>(0), -1), std::domain_error);
 378   // unknown location(0): fatal error in "test_main":
 379   // class boost::exception_detail::clone_impl<struct boost::exception_detail::error_info_injector<class std::domain_error> >:
 380   // Error in function boost::math::cyl_bessel_j_zero<double>(double, int): Requested the -1'th zero, but must be > 0 !
 381
 382   // Confirm that a C-style ignore_all policy returns NaN for bad input.
 383   typedef boost::math::policies::policy<
 384     boost::math::policies::domain_error<boost::math::policies::ignore_error>,
 385     boost::math::policies::overflow_error<boost::math::policies::ignore_error>,
 386     boost::math::policies::underflow_error<boost::math::policies::ignore_error>,
 387     boost::math::policies::denorm_error<boost::math::policies::ignore_error>,
 388     boost::math::policies::pole_error<boost::math::policies::ignore_error>,
 389     boost::math::policies::evaluation_error<boost::math::policies::ignore_error>
 390               > ignore_all_policy;
 391
 392   if (std::numeric_limits<RealType>::has_quiet_NaN)
 393   {
 394     BOOST_MATH_CHECK_THROW(cyl_bessel_j_zero(static_cast<RealType>(std::numeric_limits<RealType>::quiet_NaN()), 1), std::domain_error);
 395     // Check that bad m returns NaN if policy is no throws.
 396     BOOST_CHECK((boost::math::isnan<RealType>)(cyl_bessel_j_zero(std::numeric_limits<RealType>::quiet_NaN(), 1, ignore_all_policy())) );
 397     BOOST_MATH_CHECK_THROW(boost::math::cyl_bessel_j_zero(static_cast<RealType>(std::numeric_limits<RealType>::quiet_NaN()), -1), std::domain_error);
 398   }
 399   else
 400   { // real_concept bad m returns zero.
 401     //std::cout << boost::math::cyl_bessel_j_zero(static_cast<RealType>(0), -1, ignore_all_policy()) << std::endl; // 0 for real_concept.
 402       BOOST_CHECK_EQUAL(boost::math::cyl_bessel_j_zero(static_cast<RealType>(0), -1, ignore_all_policy() ), 0);
 403   }
 404
 405   if (std::numeric_limits<RealType>::has_infinity)
 406   {
 407     BOOST_MATH_CHECK_THROW(cyl_bessel_j_zero(std::numeric_limits<RealType>::infinity(), 0), std::domain_error);
 408     BOOST_MATH_CHECK_THROW(cyl_bessel_j_zero(std::numeric_limits<RealType>::infinity(), 1), std::domain_error);
 409     // Check that NaN is returned if error ignored.
 410     BOOST_CHECK((boost::math::isnan<RealType>)(cyl_bessel_j_zero(std::numeric_limits<RealType>::infinity(), 1, ignore_all_policy())) );
 411   }
 412
 413   // Tests of cyc_neumann zero function (BesselYZero in Wolfram) for spot values.
 414   /*
 415   Table[N[BesselYZero[0, n], 50], {n, 1, 5, 1}]
 416 n |
 417 1 | 0.89357696627916752158488710205833824122514686193001
 418 2 | 3.9576784193148578683756771869174012814186037655636
 419 3 | 7.0860510603017726976236245968203524689715103811778
 420 4 | 10.222345043496417018992042276342187125994059613181
 421 5 | 13.361097473872763478267694585713786426579135174880
 422
 423 Table[N[BesselYZero[0, n], 50], {n, 1, 5, 1}]
 424
 425 n |
 426 1 | 0.89357696627916752158488710205833824122514686193001
 427 2 | 3.9576784193148578683756771869174012814186037655636
 428 3 | 7.0860510603017726976236245968203524689715103811778
 429 4 | 10.222345043496417018992042276342187125994059613181
 430 5 | 13.361097473872763478267694585713786426579135174880
 431
 432 So K == Y
 433
 434 Table[N[BesselYZero[1, n], 50], {n, 1, 5, 1}]
 435 n |
 436 1 | 2.1971413260310170351490335626989662730530183315003
 437 2 | 5.4296810407941351327720051908525841965837574760291
 438 3 | 8.5960058683311689264296061801639678511029215669749
 439 4 | 11.749154830839881243399421939922350714301165983279
 440 5 | 14.897442128336725378844819156429870879807150630875
 441
 442 Table[N[BesselYZero[2, n], 50], {n, 1, 5, 1}]
 443 n |
 444 1 | 3.3842417671495934727014260185379031127323883259329
 445 2 | 6.7938075132682675382911671098369487124493222183854
 446 3 | 10.023477979360037978505391792081418280789658279097
 447 4 | 13.209986710206416382780863125329852185107588501072
 448 5 | 16.378966558947456561726714466123708444627678549687
 449
 450 */
 451   // Some simple integer values.
 452
 453   using boost::math::cyl_neumann_zero;
 454   // Bad rank m.
 455   BOOST_MATH_CHECK_THROW(cyl_neumann_zero(static_cast<RealType>(0), 0), std::domain_error); //
 456   BOOST_MATH_CHECK_THROW(cyl_neumann_zero(static_cast<RealType>(0), -1), std::domain_error);
 457
 458   if (std::numeric_limits<RealType>::has_quiet_NaN)
 459   {
 460     BOOST_MATH_CHECK_THROW(cyl_neumann_zero(std::numeric_limits<RealType>::quiet_NaN(), 1), std::domain_error);
 461     BOOST_MATH_CHECK_THROW(cyl_neumann_zero(static_cast<RealType>(0), -1), std::domain_error);
 462   }
 463   if (std::numeric_limits<RealType>::has_infinity)
 464   {
 465     BOOST_MATH_CHECK_THROW(cyl_neumann_zero(std::numeric_limits<RealType>::infinity(), 2), std::domain_error);
 466     BOOST_MATH_CHECK_THROW(cyl_neumann_zero(static_cast<RealType>(0), -1), std::domain_error);
 467   }
 468   // else no infinity tests.
 469
 470   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(0), 1), static_cast<RealType>(0.89357696627916752158488710205833824122514686193001L), tolerance);
 471   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(1), 2), static_cast<RealType>(5.4296810407941351327720051908525841965837574760291L), tolerance);
 472   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(2), 3), static_cast<RealType>(10.023477979360037978505391792081418280789658279097L), tolerance);
 473
 474   /*
 475   Table[N[BesselYZero[3, n], 50], {n, 1, 5, 1}]
 476   1 | 4.5270246611496438503700268671036276386651555486109
 477   2 | 8.0975537628604907044022139901128042290432231369075
 478   3 | 11.396466739595866739252048190629504945984969192535
 479   4 | 14.623077742393873174076722507725200649352970569915
 480   5 | 17.818455232945520262553239064736739443380352162752
 481   */
 482
 483     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(3), 1), static_cast<RealType>(4.5270246611496438503700268671036276386651555486109L), tolerance * 2);
 484     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(3), 2), static_cast<RealType>(8.0975537628604907044022139901128042290432231369075L), tolerance * 2);
 485     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(3), 3), static_cast<RealType>(11.396466739595866739252048190629504945984969192535L), tolerance);
 486     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(3), 4), static_cast<RealType>(14.623077742393873174076722507725200649352970569915L), tolerance);
 487     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(3), 5), static_cast<RealType>(17.818455232945520262553239064736739443380352162752L), tolerance);
 488
 489
 490     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-3), 1), static_cast<RealType>(4.5270246611496438503700268671036276386651555486109L), tolerance);
 491     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-3), 2), static_cast<RealType>(8.0975537628604907044022139901128042290432231369075L), tolerance);
 492     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-3), 3), static_cast<RealType>(11.396466739595866739252048190629504945984969192535L), tolerance);
 493     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-3), 4), static_cast<RealType>(14.623077742393873174076722507725200649352970569915L), tolerance);
 494     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-3), 5), static_cast<RealType>(17.818455232945520262553239064736739443380352162752L), tolerance);
 495
 496   { // Repeat rest using multiple zeros version.
 497     std::vector<RealType> zeros;
 498     cyl_neumann_zero(static_cast<RealType>(0.0), 1, 3, std::back_inserter(zeros) );
 499     BOOST_CHECK_CLOSE_FRACTION(zeros[0], static_cast<RealType>(0.89357696627916752158488710205833824122514686193001L), tolerance);
 500     BOOST_CHECK_CLOSE_FRACTION(zeros[1], static_cast<RealType>(3.9576784193148578683756771869174012814186037655636L), tolerance);
 501     BOOST_CHECK_CLOSE_FRACTION(zeros[2], static_cast<RealType>(7.0860510603017726976236245968203524689715103811778L), tolerance);
 502   }
 503   // Order 0: Something always tends to go wrong at zero.
 504
 505   /* Order 219/100: This checks accuracy in a region just below a critical cutoff.
 506
 507   Table[N[BesselKZero[219/100, n], 50], {n, 1, 20, 4}]
 508 1 | 3.6039149425338727979151181355741147312162055042157
 509 5 | 16.655399111666833825247894251535326778980614938275
 510 9 | 29.280564448169163756478439692311605757712873534942
 511 13 | 41.870269811145814760551599481942750124112093564643
 512 17 | 54.449180021209532654553613813754733514317929678038
 513   */
 514
 515   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(219)/100, 1), static_cast<RealType>(3.6039149425338727979151181355741147312162055042157L), tolerance);
 516   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(219)/100, 5), static_cast<RealType>(16.655399111666833825247894251535326778980614938275L), tolerance);
 517   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(219)/100, 17), static_cast<RealType>(54.449180021209532654553613813754733514317929678038L), tolerance);
 518
 519   /* Order 221/100: This checks a region just above a critical cutoff.
 520   Table[N[BesselYZero[220/100, n], 50], {n, 1, 20, 5}]
 521   1 | 3.6154383428745996706772556069431792744372398748425
 522   6 | 19.833435100254138641131431268153987585842088078470
 523   11 | 35.592602956438811360473753622212346081080817891225
 524   16 | 51.320322762482062633162699745957897178885350674038
 525   */
 526
 527   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(220)/100, 1), static_cast<RealType>(3.6154383428745996706772556069431792744372398748425L), 2 * tolerance);
 528   // Note * 2 tolerance needed - using cpp_dec_float_50 it computes exactly, probably because of extra guard digits in multiprecision decimal version.
 529   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(220)/100, 6), static_cast<RealType>(19.833435100254138641131431268153987585842088078470L), tolerance);
 530   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(220)/100, 11), static_cast<RealType>(35.592602956438811360473753622212346081080817891225L), tolerance);
 531   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(220)/100, 16), static_cast<RealType>(51.320322762482062633162699745957897178885350674038L), tolerance);
 532
 533   /*  Order 1/1000: A small order.
 534   Table[N[BesselYZero[1/1000, n], 50], {n, 1, 20, 5}]
 535   1 | 0.89502371604431360670577815537297733265776195646969
 536   6 | 16.502492490954716850993456703662137628148182892787
 537   11 | 32.206774708309182755790609144739319753463907110990
 538   16 | 47.913467031941494147962476920863688176374357572509
 539   */
 540
 541   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(1)/1000, 1), static_cast<RealType>(0.89502371604431360670577815537297733265776195646969L), 2 * tolerance);
 542   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(1)/1000, 6), static_cast<RealType>(16.5024924909547168509934567036621376281481828927870L), tolerance);
 543   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(1)/1000, 11), static_cast<RealType>(32.206774708309182755790609144739319753463907110990L), tolerance);
 544   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(1)/1000, 16), static_cast<RealType>(47.913467031941494147962476920863688176374357572509L), tolerance);
 545
 546   /* Order 71/19: Merely an intermediate order.
 547   Table[N[BesselYZero[71/19, n], 50], {n, 1, 20, 5}]
 548   1 | 5.3527167881149432911848659069476821793319749146616
 549   6 | 22.051823727778538215953091664153117627848857279151
 550   11 | 37.890091170552491176745048499809370107665221628364
 551   16 | 53.651270581421816017744203789836444968181687858095
 552   */
 553     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(71)/19, 1), static_cast<RealType>(5.3527167881149432911848659069476821793319749146616L), tolerance);
 554     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(71)/19, 6), static_cast<RealType>(22.051823727778538215953091664153117627848857279151L), tolerance);
 555     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(71)/19, 11), static_cast<RealType>(37.890091170552491176745048499809370107665221628364L), tolerance);
 556     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(71)/19, 16), static_cast<RealType>(53.651270581421816017744203789836444968181687858095L), tolerance);
 557
 558   /* Order 7001/19: A medium-large order, small enough to retain moderate efficiency of calculation.
 559
 560   Table[N[BesselYZero[7001/19, n], 50], {n, 1}]
 561    1 | 375.18866334770357669101711932706658671250621098115
 562
 563   Table[N[BesselYZero[7001/19, n], 50], {n, 2}]
 564   Standard computation time exceeded :-(
 565   */
 566   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(7001)/19, 1), static_cast<RealType>(375.18866334770357669101711932706658671250621098115L), tolerance);
 567
 568  /* A high zero such as the 1000th zero for a modest order such as 71/19.
 569    Table[N[BesselYZero[71/19, n], 50], {n, 1000}]
 570    Standard computation time exceeded :-(
 571  */
 572
 573   /*
 574   Test Negative orders cyl_neumann.
 575
 576   Table[N[BesselYZero[-1, n], 50], {n, 1, 10, 1}]
 577   1 | 2.1971413260310170351490335626989662730530183315003
 578   2 | 5.4296810407941351327720051908525841965837574760291
 579   3 | 8.5960058683311689264296061801639678511029215669749
 580   4 | 11.749154830839881243399421939922350714301165983279
 581   5 | 14.897442128336725378844819156429870879807150630875
 582   6 | 18.043402276727855564304555507889508902163088324834
 583   7 | 21.188068934142213016142481528685423196935024604904
 584   8 | 24.331942571356912035992944051850129651414333340303
 585   9 | 27.475294980449223512212285525410668235700897307021
 586   10 | 30.618286491641114715761625696447448310277939570868
 587   11 | 33.761017796109325692471759911249650993879821495802
 588   16 | 49.472505679924095824128003887609267273294894411716
 589 */
 590
 591   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-1), 1), static_cast<RealType>(2.1971413260310170351490335626989662730530183315003L), tolerance * 3);
 592   // Note this test passes at tolerance for float, double and long double, but fails for real_concept if tolerance <= 2.
 593   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-1), 6), static_cast<RealType>(18.043402276727855564304555507889508902163088324834L), tolerance * 3);
 594   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-1), 11), static_cast<RealType>(33.761017796109325692471759911249650993879821495802L), tolerance * 3);
 595   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-1), 16), static_cast<RealType>(49.472505679924095824128003887609267273294894411716L), tolerance * 3);
 596
 597  /*
 598   Table[N[BesselYZero[-2, n], 50], {n, 1, 20, 5}]
 599   1 | 3.3842417671495934727014260185379031127323883259329
 600   6 | 19.539039990286384411511740683423888947393156497603
 601   11 | 35.289793869635804143323234828826075805683602368473
 602   16 | 51.014128749483902310217774804582826908060740157564
 603   */
 604
 605   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-2), 1), static_cast<RealType>(3.3842417671495934727014260185379031127323883259329L), tolerance * 3);
 606   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-2), 6), static_cast<RealType>(19.539039990286384411511740683423888947393156497603L), tolerance * 3);
 607   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-2), 11), static_cast<RealType>(35.289793869635804143323234828826075805683602368473L), tolerance * 3);
 608   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-2), 16), static_cast<RealType>(51.014128749483902310217774804582826908060740157564L), tolerance * 3);
 609
 610 /*
 611   Table[N[BesselYZero[-3, n], 51], {n, 1, 7, 1}]
 612   1 | 4.52702466114964385037002686710362763866515554861094
 613   2 | 8.09755376286049070440221399011280422904322313690750
 614   3 | 11.3964667395958667392520481906295049459849691925349
 615   4 | 14.6230777423938731740767225077252006493529705699150
 616   5 | 17.8184552329455202625532390647367394433803521627517
 617   6 | 20.9972847541877606834525058939528641630713169437070
 618   7 | 24.1662357585818282287385597668220226288453739040042
 619 */
 620
 621   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-3), 1), static_cast<RealType>(4.52702466114964385037002686710362763866515554861094L), tolerance);
 622   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-3), 2), static_cast<RealType>(8.09755376286049070440221399011280422904322313690750L), tolerance);
 623
 624   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-3), 4), static_cast<RealType>(14.6230777423938731740767225077252006493529705699150L), tolerance);
 625
 626 /* Table[N[BesselKZero[-39, n], 51], {n, 1, 20, 5}]
 627   1 | 42.2362394762664681287397356668342141701037684436723
 628   6 | 65.8250353430045981408288669790173009159561533403819
 629   11 | 84.2674082411341814641248554679382420802125973458922
 630   16 | 101.589776978258493441843447810649346266014624868410
 631 */
 632   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-39), 1), static_cast<RealType>(42.2362394762664681287397356668342141701037684436723L), tolerance );
 633   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-39), 6), static_cast<RealType>(65.8250353430045981408288669790173009159561533403819L), tolerance);
 634   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-39), 11), static_cast<RealType>(84.2674082411341814641248554679382420802125973458922L), tolerance);
 635   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-39), 16), static_cast<RealType>(101.589776978258493441843447810649346266014624868410L), tolerance);
 636
 637 /* Table[N[BesselKZero[-39 -(1/3), n], 51], {n, 1, 20, 5}]
 638   1 | 39.3336965099558453809241429692683050137281997313679
 639   6 | 64.9038181444904768984884565999608291433823953030822
 640   11 | 83.4922341795560713832607574604255239776551554961143
 641   16 | 100.878386349724826125265571457142254077564666532665
 642 */
 643
 644   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-39) - static_cast<RealType>(1)/3, 1), static_cast<RealType>(39.3336965099558453809241429692683050137281997313679L), tolerance * 4);
 645   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-39) - static_cast<RealType>(1)/3, 6), static_cast<RealType>(64.9038181444904768984884565999608291433823953030822L), tolerance * 4);
 646   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-39) - static_cast<RealType>(1)/3, 11), static_cast<RealType>(83.4922341795560713832607574604255239776551554961143L), tolerance * 4);
 647   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-39) - static_cast<RealType>(1)/3, 16), static_cast<RealType>(100.878386349724826125265571457142254077564666532665L), tolerance * 4);
 648 /* Table[N[BesselKZero[-(1/3), n], 51], {n, 1, 20, 5}]
 649 n |
 650 1 | 0.364442931311036254896373762996743259918847602789703
 651 6 | 15.9741013584105984633772025789145590038676373673203
 652 11 | 31.6799168750213003020847708007848147516190373648194
 653 16 | 47.3871543280673235432396563497681616285970326011211
 654 */
 655
 656     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(-static_cast<RealType>(1)/3, 1), static_cast<RealType>(0.364442931311036254896373762996743259918847602789703L), tolerance * 10);
 657     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(-static_cast<RealType>(1)/3, 6), static_cast<RealType>(15.9741013584105984633772025789145590038676373673203L), tolerance * 10);
 658     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(-static_cast<RealType>(1)/3, 11), static_cast<RealType>(31.6799168750213003020847708007848147516190373648194L), tolerance * 4);
 659     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(-static_cast<RealType>(1)/3, 16), static_cast<RealType>(47.3871543280673235432396563497681616285970326011211L), tolerance * 4);
 660
 661 /* Table[N[BesselKZero[-3 -(9999/10000), n], 51], {n, 1, 20, 5}]
 662   1 | 5.64546089250283694562642537496601708928630550185069
 663   2 | 9.36184180108088288881787970896747209376324330610979
 664   3 | 12.7303431758275183078115963473808796340618061355885
 665   4 | 15.9998152121877557837972245675029531998475502716021
 666   6 | 9.36184180108088288881787970896747209376324330610979
 667   11 | 25.6104419106589739931633042959774157385787405502820
 668   16 | 41.4361281441868132581487460354904567452973524446193
 669 */
 670
 671     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-3) -static_cast<RealType>(9999)/10000, 1), static_cast<RealType>(5.64546089250283694562642537496601708928630550185069L), tolerance * 4);
 672     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-3) -static_cast<RealType>(9999)/10000, 2), static_cast<RealType>(9.36184180108088288881787970896747209376324330610979L), tolerance * 4);
 673     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-3) -static_cast<RealType>(9999)/10000, 3), static_cast<RealType>(12.7303431758275183078115963473808796340618061355885L), tolerance * 4);
 674     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-3) -static_cast<RealType>(9999)/10000, 4), static_cast<RealType>(15.9998152121877557837972245675029531998475502716021L), tolerance * 4);
 675
 676 /*  Table[N[BesselYZero[-3 -(9999/10000), n], 51], {n, 1, 7, 1}]
 677 1 | 5.64546089250283694562642537496601708928630550185069
 678 2 | 9.36184180108088288881787970896747209376324330610979
 679 3 | 12.7303431758275183078115963473808796340618061355885
 680 4 | 15.9998152121877557837972245675029531998475502716021
 681
 682 // but 5 is same as 1!!  Acknowledged as fault Wolfram [TS 6475] 26 Feb 13.
 683
 684 5 | 5.64546089250283694562642537496601708928630550184982
 685 6 | 9.36184180108088288881787970896747209376324330610979
 686 7 | 12.7303431758275183078115963473808796340618061355885
 687
 688 In[26]:= FindRoot[BesselY[-3 -9999/10000, r] == 0, {r, 3}]  for r = 2,3, 4, 5 = {r->5.64546}
 689
 690 In[26]:= FindRoot[BesselY[-3 -9999/10000, r] == 0, {r, 19}] = 19.2246
 691
 692 So no very accurate reference value for these.
 693
 694 Calculated using cpp_dec_float_50
 695
 696   5.6454608925028369456264253749660170892863055018498
 697   9.3618418010808828888178797089674720937632433061099
 698   12.730343175827518307811596347380879634061806135589
 699   15.999815212187755783797224567502953199847550271602
 700
 701   19.224610865671563344572152795434688888375602299773
 702   22.424988389021059116212186912990863561607855849204
 703   25.610441910658973993163304295977415738578740550282
 704   28.786066313968546073981640755202085944374967166411
 705   31.954857624676521867923579695253822854717613513587
 706     */
 707   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-3) -static_cast<RealType>(9999)/10000, 1), static_cast<RealType>(5.64546089250283694562642537496601708928630550185069L), tolerance * 4);
 708   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-3) -static_cast<RealType>(9999)/10000, 2), static_cast<RealType>(9.36184180108088288881787970896747209376324330610979L), tolerance * 4);
 709   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-3) -static_cast<RealType>(9999)/10000, 3), static_cast<RealType>(12.7303431758275183078115963473808796340618061355885L), tolerance * 4);
 710   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-3) -static_cast<RealType>(9999)/10000, 4), static_cast<RealType>(15.9998152121877557837972245675029531998475502716021L), tolerance * 4);
 711 //
 712   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-3) -static_cast<RealType>(9999)/10000, 5), static_cast<RealType>(19.224610865671563344572152795434688888375602299773L), tolerance * 4);
 713   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-3) -static_cast<RealType>(9999)/10000, 6), static_cast<RealType>(22.424988389021059116212186912990863561607855849204L), tolerance * 4);
 714   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-3) -static_cast<RealType>(9999)/10000, 7), static_cast<RealType>(25.610441910658973993163304295977415738578740550282L), tolerance * 4);
 715   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-3) -static_cast<RealType>(9999)/10000, 8), static_cast<RealType>(28.786066313968546073981640755202085944374967166411L), tolerance * 4);
 716   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-3) -static_cast<RealType>(9999)/10000, 9), static_cast<RealType>(31.954857624676521867923579695253822854717613513587L), tolerance * 4);
 717
 718
 719 // Plot[BesselYZero[-7 - v, 1], {v, 0, 1}]  shows discontinuity at the mid-point between integers.
 720
 721 /* Table[N[BesselYZero[-7 - (4999/10000), n], 51], {n, 1, 4, 1}]
 722   1 | 3.59209698655443348407622952525352410710983745802573
 723   2 | 11.6573245781899449398248761667833391837824916603434
 724   3 | 15.4315262542144355217979771618575628291362029097236
 725   4 | 18.9232143766706670333395285892576635207736306576135
 726 */
 727
 728 /* Table[N[BesselYZero[-7 - (5001/10000), n], 51], {n, 1, 4, 1}]
 729   1 | 11.6567397956147934678808863468662427054245897492445
 730   2 | 15.4310521624769624067699131497395566368341140531722
 731   3 | 18.9227840182910629037411848072684247564491740961847
 732   4 | 22.2951449444372591060253508661432751300205474374696
 733 */
 734
 735     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-7) -static_cast<RealType>(4999)/10000, 1), static_cast<RealType>(3.59209698655443348407622952525352410710983745802573L), tolerance * 2000);
 736     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-7) -static_cast<RealType>(4999)/10000, 2), static_cast<RealType>(11.6573245781899449398248761667833391837824916603434L), tolerance * 100);
 737     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-7) -static_cast<RealType>(4999)/10000, 3), static_cast<RealType>(15.4315262542144355217979771618575628291362029097236L), tolerance * 100);
 738     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-7) -static_cast<RealType>(4999)/10000, 4), static_cast<RealType>(18.9232143766706670333395285892576635207736306576135L), tolerance * 100);
 739
 740     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-7) -static_cast<RealType>(5001)/10000, 1), static_cast<RealType>(11.6567397956147934678808863468662427054245897492445L), tolerance * 100);
 741     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-7) -static_cast<RealType>(5001)/10000, 2), static_cast<RealType>(15.4310521624769624067699131497395566368341140531722L), tolerance * 100);
 742     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-7) -static_cast<RealType>(5001)/10000, 3), static_cast<RealType>(18.9227840182910629037411848072684247564491740961847L), tolerance * 100);
 743     BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-7) -static_cast<RealType>(5001)/10000, 4), static_cast<RealType>(22.2951449444372591060253508661432751300205474374696L), tolerance * 100);
 744
 745   //BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(-(static_cast<RealType>(-3)-static_cast<RealType>(99)/100), 5),
 746   //  cyl_neumann_zero(+(static_cast<RealType>(-3)-static_cast<RealType>(99)/100), 5), tolerance * 100);
 747   {
 748     long double x = 1.L;
 749     BOOST_CHECK_CLOSE_FRACTION(
 750       cyl_neumann_zero(-(static_cast<RealType>(x)), 5),
 751       cyl_neumann_zero(+(static_cast<RealType>(x)), 5), tolerance * 100);
 752   }
 753   {
 754     long double x = 2.L;
 755     BOOST_CHECK_CLOSE_FRACTION(
 756       cyl_neumann_zero(-(static_cast<RealType>(x)), 5),
 757       cyl_neumann_zero(+(static_cast<RealType>(x)), 5), tolerance * 100);
 758   }
 759   {
 760     long double x = 3.L;
 761     BOOST_CHECK_CLOSE_FRACTION(
 762       cyl_neumann_zero(-(static_cast<RealType>(x)), 5),
 763       cyl_neumann_zero(+(static_cast<RealType>(x)), 5), tolerance * 100);
 764   }
 765   // These are very close but not exactly same.
 766   //{
 767   //  RealType x = static_cast<RealType>(-3) -static_cast<RealType>(9999)/10000;
 768   //  BOOST_CHECK_CLOSE_FRACTION(
 769   //    cyl_neumann_zero(-(static_cast<RealType>(x)), 5),
 770   //    cyl_neumann_zero(+(static_cast<RealType>(x)), 5), tolerance * 100);
 771   //  // 19.2242889 and 19.2246113
 772   //}
 773   //{
 774
 775   //  RealType x = static_cast<RealType>(-3) -static_cast<RealType>(9999)/10000;
 776   //  BOOST_CHECK_CLOSE_FRACTION(
 777   //    cyl_neumann_zero(-(static_cast<RealType>(x)), 6),
 778   //    cyl_neumann_zero(+(static_cast<RealType>(x)), 6), tolerance * 100);
 779   //  // 22.4246693 and 22.4249878
 780   //}
 781
 782
 783
 784   //  2.5  18.6890354  17.1033592
 785
 786
 787 /*Table[N[BesselYZero[-1/81799, n], 51], {n, 1, 10, 5}]
 788
 789 1 | 0.893559276290122922836047849416713592133322804889757
 790 2 | 3.95765935645507004204986415533750122885237402118726
 791 3 | 7.08603190350579828577279552434514387474680226004173
 792 4 | 10.2223258629823064789904339889550588869985272176335
 793 5 | 13.3610782840659145864973521693322670264135672594988
 794 3 | 7.08603190350579828577279552434514387474680226004173
 795 5 | 13.3610782840659145864973521693322670264135672594988
 796 6 | 16.5009032471619898684110089652474861084220781491575
 797 7 | 19.6412905039556082160052482410981245043314155416354
 798 9 | 25.9229384536173175152381652048590136247796591153244
 799 */
 800    BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(-static_cast<RealType>(1)/81799, 1), static_cast<RealType>(0.893559276290122922836047849416713592133322804889757L), tolerance * 4);
 801   // Doesn't converge!
 802   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(-static_cast<RealType>(1)/81799, 2), static_cast<RealType>(3.95765935645507004204986415533750122885237402118726L), tolerance * 4);
 803   BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(-static_cast<RealType>(1)/81799, 3), static_cast<RealType>(7.08603190350579828577279552434514387474680226004173L), tolerance * 4);
 804   /* try positive x
 805   Table[N[BesselYZero[1/81799, n], 51], {n, 1, 5, 1}]
 806
 807   1 | 0.893594656187326273432267210617481926490785928764963
 808   2 | 3.95769748213950546166537901626409026826595687994956
 809   3 | 7.08607021707716361104064671367526817399129653285580
 810   4 | 10.2223642239960815612515914411615233651316361060338
 811   5 | 13.3611166636685056799674772287389749065996094266976
 812 */
 813    BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(1)/81799, 2), static_cast<RealType>(3.95769748213950546166537901626409026826595687994956L), tolerance * 4);
 814    BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(1)/81799, 3), static_cast<RealType>(7.08607021707716361104064671367526817399129653285580L), tolerance * 4);
 815    BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(-static_cast<RealType>(1)/81799, 4), static_cast<RealType>(10.2223258629823064789904339889550588869985272176335L), tolerance * 4);
 816    BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(-static_cast<RealType>(1)/81799, 5), static_cast<RealType>(13.3610782840659145864973521693322670264135672594988L), tolerance * 4);
 817    BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(-static_cast<RealType>(1)/81799, 6), static_cast<RealType>(16.5009032471619898684110089652474861084220781491575L), tolerance * 4);
 818    BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(-static_cast<RealType>(1)/81799, 9), static_cast<RealType>(25.9229384536173175152381652048590136247796591153244L), tolerance * 4);
 819    BOOST_CHECK_CLOSE_FRACTION(cyl_neumann_zero(static_cast<RealType>(-7) - static_cast<RealType>(1)/3, 1), static_cast<RealType>(7.3352783956690540155848592759652828459644819344081L), tolerance * 1000);
 820
 821   // Test Data for airy_ai_zero and airy_bi_zero functions.
 822
 823    using boost::math::airy_ai_zero; //
 824
 825    using boost::math::isnan;
 826
 827    BOOST_MATH_CHECK_THROW(airy_ai_zero<RealType>(0), std::domain_error);
 828
 829   if (std::numeric_limits<RealType>::has_quiet_NaN)
 830   { // If ignore errors, return NaN.
 831     BOOST_CHECK((boost::math::isnan)(airy_ai_zero<RealType>(0, ignore_all_policy())));
 832     BOOST_CHECK((boost::math::isnan)(airy_ai_zero<RealType>((std::numeric_limits<unsigned>::min)() , ignore_all_policy())));
 833     // Can't abuse with NaN as won't compile.
 834     //BOOST_MATH_CHECK_THROW(airy_ai_zero<RealType>(std::numeric_limits<RealType>::quiet_NaN()), std::domain_error);
 835   }
 836   else
 837   { // real_concept NaN not available, so return zero.
 838     BOOST_CHECK_EQUAL(airy_ai_zero<RealType>(0, ignore_all_policy()), 0);
 839     // BOOST_CHECK_EQUAL(airy_ai_zero<RealType>(-1), 0); //  warning C4245: 'argument' : conversion from 'int' to 'unsigned int', signed/unsigned mismatch
 840   }
 841
 842   BOOST_MATH_CHECK_THROW(airy_ai_zero<RealType>(-1), std::domain_error);
 843   BOOST_CHECK_CLOSE_FRACTION(airy_ai_zero<RealType>((std::numeric_limits<boost::int32_t>::max)()), -static_cast<RealType>(4678579.33301973093739L), tolerance);
 844
 845   // Can't abuse with infinity because won't compile - no conversion.
 846   //if (std::numeric_limits<RealType>::has_infinity)
 847   //{
 848   //  BOOST_CHECK(isnan(airy_bi_zero<RealType>(-1)) );
 849   //}
 850
 851   // WolframAlpha  Table[N[AiryAiZero[n], 51], {n, 1, 20, 1}]
 852
 853   BOOST_CHECK_CLOSE_FRACTION(airy_ai_zero<RealType>(1), static_cast<RealType>(-2.33810741045976703848919725244673544063854014567239L), tolerance * 2 * tolerance_tgamma_extra);
 854   BOOST_CHECK_CLOSE_FRACTION(airy_ai_zero<RealType>(2), static_cast<RealType>(-4.08794944413097061663698870145739106022476469910853L), tolerance);
 855   BOOST_CHECK_CLOSE_FRACTION(airy_ai_zero<RealType>(3), static_cast<RealType>(-5.52055982809555105912985551293129357379721428061753L), tolerance);
 856   BOOST_CHECK_CLOSE_FRACTION(airy_ai_zero<RealType>(4), static_cast<RealType>(-6.78670809007175899878024638449617696605388247739349L), tolerance);
 857   BOOST_CHECK_CLOSE_FRACTION(airy_ai_zero<RealType>(5), static_cast<RealType>(-7.94413358712085312313828055579826853214067439697221L), tolerance);
 858   BOOST_CHECK_CLOSE_FRACTION(airy_ai_zero<RealType>(6), static_cast<RealType>(-9.02265085334098038015819083988008925652467753515608L), tolerance);
 859   BOOST_CHECK_CLOSE_FRACTION(airy_ai_zero<RealType>(7), static_cast<RealType>(-10.0401743415580859305945567373625180940429025691058L), tolerance);
 860   BOOST_CHECK_CLOSE_FRACTION(airy_ai_zero<RealType>(8), static_cast<RealType>(-11.0085243037332628932354396495901510167308253815040L), tolerance);
 861   BOOST_CHECK_CLOSE_FRACTION(airy_ai_zero<RealType>(9), static_cast<RealType>(-11.9360155632362625170063649029305843155778862321198L), tolerance);
 862   BOOST_CHECK_CLOSE_FRACTION(airy_ai_zero<RealType>(10), static_cast<RealType>(-12.8287767528657572004067294072418244773864155995734L), tolerance);
 863   BOOST_CHECK_CLOSE_FRACTION(airy_ai_zero<RealType>(11), static_cast<RealType>(-13.6914890352107179282956967794669205416653698092008L), tolerance);
 864   BOOST_CHECK_CLOSE_FRACTION(airy_ai_zero<RealType>(12), static_cast<RealType>(-14.5278299517753349820739814429958933787141648698348L), tolerance);
 865   BOOST_CHECK_CLOSE_FRACTION(airy_ai_zero<RealType>(13), static_cast<RealType>(-15.3407551359779968571462085134814867051175833202480L), tolerance);
 866   BOOST_CHECK_CLOSE_FRACTION(airy_ai_zero<RealType>(14), static_cast<RealType>(-16.1326851569457714393459804472025217905182723970763L), tolerance);
 867   BOOST_CHECK_CLOSE_FRACTION(airy_ai_zero<RealType>(15), static_cast<RealType>(-16.9056339974299426270352387706114765990900510950317L), tolerance);
 868   BOOST_CHECK_CLOSE_FRACTION(airy_ai_zero<RealType>(16), static_cast<RealType>(-17.6613001056970575092536503040180559521532186681200L), tolerance);
 869   BOOST_CHECK_CLOSE_FRACTION(airy_ai_zero<RealType>(17), static_cast<RealType>(-18.4011325992071154158613979295043367545938146060201L), tolerance);
 870   BOOST_CHECK_CLOSE_FRACTION(airy_ai_zero<RealType>(18), static_cast<RealType>(-19.1263804742469521441241486897324946890754583847531L), tolerance);
 871   BOOST_CHECK_CLOSE_FRACTION(airy_ai_zero<RealType>(19), static_cast<RealType>(-19.8381298917214997009475636160114041983356824945389L), tolerance);
 872   BOOST_CHECK_CLOSE_FRACTION(airy_ai_zero<RealType>(20), static_cast<RealType>(-20.5373329076775663599826814113081017453042180147375L), tolerance);
 873
 874   // Table[N[AiryAiZero[n], 51], {n, 1000, 1001, 1}]
 875
 876   BOOST_CHECK_CLOSE_FRACTION(airy_ai_zero<RealType>(1000), static_cast<RealType>(-281.031519612521552835336363963709689055717463965420L), tolerance);
 877   BOOST_CHECK_CLOSE_FRACTION(airy_ai_zero<RealType>(1001), static_cast<RealType>(-281.218889579130068414512015874511112547569713693446L), tolerance);
 878
 879   // Table[N[AiryAiZero[n], 51], {n, 1000000, 1000001, 1}]
 880   BOOST_CHECK_CLOSE_FRACTION(airy_ai_zero<RealType>(1000000), static_cast<RealType>(-28107.8319793795834876064419863203282898723750036048L), tolerance);
 881   BOOST_CHECK_CLOSE_FRACTION(airy_ai_zero<RealType>(1000001), static_cast<RealType>(-28107.8507179357979542838020057465277368471496446555L), tolerance);
 882
 883
 884   // Table[N[AiryAiZero[n], 51], {n, 1000000000, 1000000001, 1}]
 885   BOOST_CHECK_CLOSE_FRACTION(airy_ai_zero<RealType>(1000000000), static_cast<RealType>(-2.81078366593344513918947921096193426320298300481145E+6L), tolerance);
 886   BOOST_CHECK_CLOSE_FRACTION(airy_ai_zero<RealType>(1000000001), static_cast<RealType>(-2.81078366780730091663459728526906320267920607427246E+6L), tolerance);
 887
 888   // Test Data for airy_bi
 889   using boost::math::airy_bi_zero;
 890
 891   BOOST_MATH_CHECK_THROW(airy_bi_zero<RealType>(0), std::domain_error);
 892
 893   if (std::numeric_limits<RealType>::has_quiet_NaN)
 894   { // return NaN.
 895     BOOST_CHECK((boost::math::isnan)(airy_bi_zero<RealType>(0, ignore_all_policy())));
 896     BOOST_CHECK((boost::math::isnan)(airy_bi_zero<RealType>((std::numeric_limits<unsigned>::min)() , ignore_all_policy())));
 897     // Can't abuse with NaN as won't compile.
 898     // BOOST_MATH_CHECK_THROW(airy_bi_zero<RealType>(std::numeric_limits<RealType>::quiet_NaN()), std::domain_error);
 899     // cannot convert parameter 1 from 'boost::math::concepts::real_concept' to 'unsigned int'.
 900   }
 901   else
 902   { // real_concept NaN not available, so return zero.
 903     BOOST_CHECK_EQUAL(airy_bi_zero<RealType>(0, ignore_all_policy()), 0);
 904     // BOOST_CHECK_EQUAL(airy_bi_zero<RealType>(-1), 0);
 905     //  warning C4245: 'argument' : conversion from 'int' to 'unsigned int', signed/unsigned mismatch.
 906     // If ignore the warning, interpreted as max unsigned:
 907     // check airy_bi_zero<RealType>(-1) == 0 has failed [-7.42678e+006 != 0]
 908   }
 909
 910   BOOST_MATH_CHECK_THROW(airy_bi_zero<RealType>(-1), std::domain_error);
 911   BOOST_CHECK_CLOSE_FRACTION(airy_bi_zero<RealType>((std::numeric_limits<boost::int32_t>::max)()), -static_cast<RealType>(4678579.33229351984573L), tolerance * 300);
 912
 913   // Can't abuse with infinity because won't compile - no conversion.
 914   //if (std::numeric_limits<RealType>::has_infinity)
 915   //{
 916   //  BOOST_CHECK(isnan(airy_bi_zero<RealType>(std::numeric_limits<RealType>::infinity)) );
 917   //}
 918
 919   // Table[N[AiryBiZero[n], 51], {n, 1, 20, 1}]
 920   BOOST_CHECK_CLOSE_FRACTION(airy_bi_zero<RealType>(1), static_cast<RealType>(-1.17371322270912792491997996247390210454364638917570L), tolerance * 4 * tolerance_tgamma_extra);
 921   BOOST_CHECK_CLOSE_FRACTION(airy_bi_zero<RealType>(2), static_cast<RealType>(-3.27109330283635271568022824016641380630093596910028L), tolerance * tolerance_tgamma_extra);
 922   BOOST_CHECK_CLOSE_FRACTION(airy_bi_zero<RealType>(3), static_cast<RealType>(-4.83073784166201593266770933990517817696614261732301L), tolerance);
 923   BOOST_CHECK_CLOSE_FRACTION(airy_bi_zero<RealType>(4), static_cast<RealType>(-6.16985212831025125983336452055593667996554943427563L), tolerance);
 924   BOOST_CHECK_CLOSE_FRACTION(airy_bi_zero<RealType>(5), static_cast<RealType>(-7.37676207936776371359995933044254122209152229939710L), tolerance);
 925   BOOST_CHECK_CLOSE_FRACTION(airy_bi_zero<RealType>(6), static_cast<RealType>(-8.49194884650938801344803949280977672860508755505546L), tolerance);
 926   BOOST_CHECK_CLOSE_FRACTION(airy_bi_zero<RealType>(7), static_cast<RealType>(-9.53819437934623888663298854515601962083907207638247L), tolerance);
 927   BOOST_CHECK_CLOSE_FRACTION(airy_bi_zero<RealType>(8), static_cast<RealType>(-10.5299135067053579244005555984531479995295775946214L), tolerance);
 928   BOOST_CHECK_CLOSE_FRACTION(airy_bi_zero<RealType>(9), static_cast<RealType>(-11.4769535512787794379234649247328196719482538148877L), tolerance);
 929   BOOST_CHECK_CLOSE_FRACTION(airy_bi_zero<RealType>(10), static_cast<RealType>(-12.3864171385827387455619015028632809482597983846856L), tolerance);
 930   BOOST_CHECK_CLOSE_FRACTION(airy_bi_zero<RealType>(11), static_cast<RealType>(-13.2636395229418055541107433243954907752411519609813L), tolerance);
 931   BOOST_CHECK_CLOSE_FRACTION(airy_bi_zero<RealType>(12), static_cast<RealType>(-14.1127568090686577915873097822240184716840428285509L), tolerance);
 932   BOOST_CHECK_CLOSE_FRACTION(airy_bi_zero<RealType>(13),  static_cast<RealType>(-14.9370574121541640402032143104909046396121763517782L), tolerance);
 933   BOOST_CHECK_CLOSE_FRACTION(airy_bi_zero<RealType>(14), static_cast<RealType>(-15.7392103511904827708949784797481833807180162767841L), tolerance);
 934   BOOST_CHECK_CLOSE_FRACTION(airy_bi_zero<RealType>(15), static_cast<RealType>(-16.5214195506343790539179499652105457167110310370581L), tolerance);
 935   BOOST_CHECK_CLOSE_FRACTION(airy_bi_zero<RealType>(16), static_cast<RealType>(-17.2855316245812425329342366922535392425279753602710L), tolerance);
 936   BOOST_CHECK_CLOSE_FRACTION(airy_bi_zero<RealType>(17), static_cast<RealType>(-18.0331132872250015721711125433391920008087291416406L), tolerance);
 937   BOOST_CHECK_CLOSE_FRACTION(airy_bi_zero<RealType>(18), static_cast<RealType>(-18.7655082844800810413429789236105128440267189551421L), tolerance);
 938   BOOST_CHECK_CLOSE_FRACTION(airy_bi_zero<RealType>(19), static_cast<RealType>(-19.4838801329892340136659986592413575122062977793610L), tolerance);
 939   BOOST_CHECK_CLOSE_FRACTION(airy_bi_zero<RealType>(20), static_cast<RealType>(-20.1892447853962024202253232258275360764649783583934L), tolerance);
 940
 941  // Table[N[AiryBiZero[n], 51], {n, 1000, 1001, 1}]
 942   BOOST_CHECK_CLOSE_FRACTION(airy_bi_zero<RealType>(1000), static_cast<RealType>(-280.937811203415240157883427412260300146245056425646L), tolerance);
 943   BOOST_CHECK_CLOSE_FRACTION(airy_bi_zero<RealType>(1001), static_cast<RealType>(-281.125212400956392021977771104562061554648675044114L), tolerance);
 944
 945   // Table[N[AiryBiZero[n], 51], {n, 1000000, 1000001, 1}]
 946   BOOST_CHECK_CLOSE_FRACTION(airy_bi_zero<RealType>(1000000), static_cast<RealType>(-28107.8226100991339342855024130953986989636667226163L), tolerance);
 947   BOOST_CHECK_CLOSE_FRACTION(airy_bi_zero<RealType>(1000001), static_cast<RealType>(-28107.8413486584714939255315213519230566014624895515L), tolerance);
 948
 949   //Table[N[AiryBiZero[n], 51], {n, 1000000000, 1000000001, 1}]
 950   BOOST_CHECK_CLOSE_FRACTION(airy_bi_zero<RealType>(1000000000), static_cast<RealType>(-2.81078366499651725023268820158218492845371527054171E+6L), tolerance);
 951   BOOST_CHECK_CLOSE_FRACTION(airy_bi_zero<RealType>(1000000001), static_cast<RealType>(-2.81078366687037302799011557215619265502627118526716E+6L), tolerance);
 952
 953   // Check the multi-root versions.
 954   {
 955     unsigned int n_roots = 1U;
 956     std::vector<RealType> roots;
 957     boost::math::airy_ai_zero<RealType>(2U, n_roots, std::back_inserter(roots));
 958     BOOST_CHECK_CLOSE_FRACTION(roots[0], static_cast<RealType>(-4.08794944413097061663698870145739106022476469910853L), tolerance);
 959   }
 960   {
 961     unsigned int n_roots = 1U;
 962     std::vector<RealType> roots;
 963     boost::math::airy_bi_zero<RealType>(2U, n_roots, std::back_inserter(roots));
 964     BOOST_CHECK_CLOSE_FRACTION(roots[0], static_cast<RealType>(-3.27109330283635271568022824016641380630093596910028L), tolerance * tolerance_tgamma_extra);
 965   }
 966 } // template <class RealType> void test_spots(RealType)
 967
 968   #include <boost/multiprecision/cpp_dec_float.hpp>
 969
 970 BOOST_AUTO_TEST_CASE(test_main)
 971 {
 972    test_bessel_zeros(0.1F);
 973    test_bessel_zeros(0.1);
 974 #ifndef BOOST_MATH_NO_LONG_DOUBLE_MATH_FUNCTIONS
 975    test_bessel_zeros(0.1L);
 976 #ifndef BOOST_MATH_NO_REAL_CONCEPT_TESTS
 977    test_bessel_zeros(boost::math::concepts::real_concept(0.1));
 978 #endif
 979 #else
 980    std::cout << "<note>The long double tests have been disabled on this platform "
 981       "either because the long double overloads of the usual math functions are "
 982       "not available at all, or because they are too inaccurate for these tests "
 983       "to pass.</note>" << std::endl;
 984 #endif
 985 }