libs/math/test/test_pFq_precision.cpp

   1 //  (C) Copyright John Maddock 2006.
   2 //  Use, modification and distribution are subject to the
   3 //  Boost Software License, Version 1.0. (See accompanying file
   4 //  LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
   5
   6 #define BOOST_TEST_MAIN
   7 #include <boost/test/unit_test.hpp>
   8 #include <boost/test/tools/floating_point_comparison.hpp>
   9 #include <boost/math/tools/stats.hpp>
  10 #include <boost/math/tools/test.hpp>
  11 #include <boost/math/tools/big_constant.hpp>
  12 #include <boost/math/constants/constants.hpp>
  13 #include <boost/type_traits/is_floating_point.hpp>
  14 #include <boost/array.hpp>
  15 #include "functor.hpp"
  16
  17 #include "handle_test_result.hpp"
  18 #include "table_type.hpp"
  19
  20 #include <boost/math/special_functions/hypergeometric_pFq.hpp>
  21 #include <boost/multiprecision/mpfr.hpp>
  22 #include <boost/math/special_functions/relative_difference.hpp>
  23
  24 #ifdef BOOST_MSVC
  25 #pragma warning(disable:4127)
  26 #endif
  27
  28 #ifndef SC_
  29 #define SC_(x) BOOST_MATH_BIG_CONSTANT(mp_type, 1000000, x)
  30 #endif
  31
  32 typedef boost::multiprecision::mpfr_float mp_type;
  33
  34 void test_spots_1F0()
  35 {
  36    using std::pow;
  37
  38    mp_type tolerance = 2e-20;
  39
  40    BOOST_CHECK_CLOSE_FRACTION(boost::math::hypergeometric_pFq_precision({ mp_type(-3) }, {}, mp_type(2), 20), mp_type(-1), tolerance);
  41    BOOST_CHECK_CLOSE_FRACTION(boost::math::hypergeometric_pFq_precision({ mp_type(-3) }, {}, mp_type(4), 20), mp_type(-27), tolerance);
  42    BOOST_CHECK_CLOSE_FRACTION(boost::math::hypergeometric_pFq_precision({ mp_type(-3) }, {}, mp_type(0.5), 20), mp_type(0.125), tolerance);
  43    BOOST_CHECK_CLOSE_FRACTION(boost::math::hypergeometric_pFq_precision({ mp_type(3) }, {}, mp_type(0.5), 20), mp_type(8), tolerance);
  44    BOOST_CHECK_CLOSE_FRACTION(boost::math::hypergeometric_pFq_precision({ mp_type(3) }, {}, mp_type(2), 20), mp_type(-1), tolerance);
  45    BOOST_CHECK_CLOSE_FRACTION(boost::math::hypergeometric_pFq_precision({ mp_type(3) }, {}, mp_type(4), 20), mp_type(mp_type(-1) / 27), tolerance);
  46    BOOST_CHECK_CLOSE_FRACTION(boost::math::hypergeometric_pFq_precision({ mp_type(3) }, {}, mp_type(-0.5), 20), pow(mp_type(1.5), -3), tolerance);
  47    BOOST_CHECK_CLOSE_FRACTION(boost::math::hypergeometric_pFq_precision({ mp_type(3) }, {}, mp_type(-2), 20), mp_type(1 / mp_type(27)), tolerance);
  48    BOOST_CHECK_CLOSE_FRACTION(boost::math::hypergeometric_pFq_precision({ mp_type(3) }, {}, mp_type(-4), 20), mp_type(mp_type(1) / 125), tolerance);
  49    BOOST_CHECK_CLOSE_FRACTION(boost::math::hypergeometric_pFq_precision({ mp_type(-3) }, {}, mp_type(-0.5), 20), pow(mp_type(1.5), 3), tolerance);
  50    BOOST_CHECK_CLOSE_FRACTION(boost::math::hypergeometric_pFq_precision({ mp_type(-3) }, {}, mp_type(-2), 20), mp_type(27), tolerance);
  51    BOOST_CHECK_CLOSE_FRACTION(boost::math::hypergeometric_pFq_precision({ mp_type(-3) }, {}, mp_type(-4), 20), mp_type(125), tolerance);
  52
  53    BOOST_CHECK_THROW(boost::math::hypergeometric_pFq_precision({ mp_type(3) }, {}, mp_type(1), 20), std::domain_error);
  54    BOOST_CHECK_THROW(boost::math::hypergeometric_pFq_precision({ mp_type(3.25) }, {}, mp_type(1), 20), std::domain_error);
  55    BOOST_CHECK_THROW(boost::math::hypergeometric_pFq_precision({ mp_type(3.25) }, {}, mp_type(2), 20), std::domain_error);
  56    BOOST_CHECK_THROW(boost::math::hypergeometric_pFq_precision({ mp_type(-3.25) }, {}, mp_type(2), 20), std::domain_error);
  57 }
  58
  59 void test_spots_0F1()
  60 {
  61    mp_type tolerance = 2e-20;
  62
  63    BOOST_CHECK_EQUAL(boost::math::hypergeometric_pFq_precision({}, { mp_type(3) }, mp_type(0), 20), 1);
  64    BOOST_CHECK_EQUAL(boost::math::hypergeometric_pFq_precision({}, { mp_type(-3) }, mp_type(0), 20), 1);
  65    //BOOST_CHECK_EQUAL(boost::math::hypergeometric_pFq_precision({}, { mp_type(0) }, mp_type(0), 20), 1);
  66
  67    BOOST_CHECK_THROW(boost::math::hypergeometric_pFq_precision({}, { mp_type(0) }, mp_type(-1), 20), std::domain_error);
  68    BOOST_CHECK_THROW(boost::math::hypergeometric_pFq_precision({}, { mp_type(-1) }, mp_type(-1), 20), std::domain_error);
  69    BOOST_CHECK_THROW(boost::math::hypergeometric_pFq_precision({}, { mp_type(-10) }, mp_type(-5), 20), std::domain_error);
  70
  71    static const boost::array<boost::array<mp_type, 3>, 35> hypergeometric_pFq_integer_data = { {
  72       { SC_(4.0), SC_(-20.0),  SC_(-0.012889714201783047561923257996127233830940165138385) },
  73       { SC_(8.0), SC_(-20.0),  SC_(0.046498609282365144223175012935939437508273248399881) },
  74       { SC_(12.0), SC_(-20.0),  SC_(0.16608847431869756642136191351311569335145459224622) },
  75       { SC_(16.0), SC_(-20.0),  SC_(0.27230484709157170329168048388841880599105216477631) },
  76       //{ SC_(20.0), SC_(-20.0),  SC_(0.35865872656868844615709101792040025253126126604266) },
  77       { SC_(4.0), SC_(-16.0),  SC_(-0.027293644412433023379286103818840667403690937153604) },
  78       { SC_(8.0), SC_(-16.0),  SC_(0.098618710511372349330666801041676087431136532039702) },
  79       { SC_(12.0), SC_(-16.0),  SC_(0.24360114226383905073379763460037817885919457531523) },
  80       //{ SC_(16.0), SC_(-16.0),  SC_(0.35635186318802906043824855864337727878754460163525) },
  81       //{ SC_(20.0), SC_(-16.0),  SC_(0.44218381382689101428948260613085371477815110358789) },
  82       { SC_(4.0), SC_(-12.0),  SC_(-0.021743572290699436419371120781513860006290363262907) },
  83       { SC_(8.0), SC_(-12.0),  SC_(0.19025625754362006866949730683824627505504067855043) },
  84       //{ SC_(12.0), SC_(-12.0),  SC_(0.35251228238278927379621049815222218665165551016489) },
  85       //{ SC_(16.0), SC_(-12.0),  SC_(0.46415411486674623230458980010115972932474705884865) },
  86       //{ SC_(20.0), SC_(-12.0),  SC_(0.54394918325286018927327004362535051310016558628741) },
  87       { SC_(4.0), SC_(-8.0),  SC_(0.056818744289274872033266550620647787396712125304880) },
  88       //{ SC_(8.0), SC_(-8.0),  SC_(0.34487371876996263249797701802458885718691612997456) },
  89       //{ SC_(12.0), SC_(-8.0),  SC_(0.50411654015891701804499796523449656998841355305043) },
  90       //{ SC_(16.0), SC_(-8.0),  SC_(0.60191459981670594041254437708158847428118361245442) },
  91       //{ SC_(20.0), SC_(-8.0),  SC_(0.66770752550930138035694866478078941681114294465418) },
  92       //{ SC_(4.0), SC_(-4.0),  SC_(0.32262860540671645526863760914000166725449779629143) },
  93       //{ SC_(8.0), SC_(-4.0),  SC_(0.59755773349355150397404772151441126513126998265958) },
  94       //{ SC_(12.0), SC_(-4.0),  SC_(0.71337465206009117934071859694314971137807212605147) },
  95       //{ SC_(16.0), SC_(-4.0),  SC_(0.77734333649378860739496954157535257278092349684783) },
  96       //{ SC_(20.0), SC_(-4.0),  SC_(0.81794177985447769150469288350369205683856312760890) },
  97
  98       { SC_(4.0), SC_(4.0),  SC_(2.5029568338152582758923890008139391395035041790831) },
  99       { SC_(8.0), SC_(4.0),  SC_(1.6273673128576761227855719910743734060605725722129) },
 100       { SC_(12.0), SC_(4.0),  SC_(1.3898419290864057799739567227851793491657442624207) },
 101       { SC_(16.0), SC_(4.0),  SC_(1.2817098157957427946677711269410726972209834860612) },
 102       { SC_(20.0), SC_(4.0),  SC_(1.2202539302152377230940386181201477276788392792437) },
 103       { SC_(4.0), SC_(8.0),  SC_(5.5616961007411965409200003309686924059253894118586) },
 104       { SC_(8.0), SC_(8.0),  SC_(2.5877053985451664722152913482683136948296873738479) },
 105       { SC_(12.0), SC_(8.0),  SC_(1.9166410733572697158003086323981583993970490592046) },
 106       { SC_(16.0), SC_(8.0),  SC_(1.6370675016890669952237854163997946987362497613701) },
 107       { SC_(20.0), SC_(8.0),  SC_(1.4862852701827990444915220582410007454379891584086) },
 108       { SC_(4.0), SC_(12.0),  SC_(11.419268276211177842169936131590385979116019595164) },
 109       { SC_(8.0), SC_(12.0),  SC_(4.0347215359576567066789638314925802225312840819037) },
 110       { SC_(12.0), SC_(12.0),  SC_(2.6242497527837800417573064942486918368886996538285) },
 111       { SC_(16.0), SC_(12.0),  SC_(2.0840468784170876805932772732753387258909164486511) },
 112       { SC_(20.0), SC_(12.0),  SC_(1.8071042457762091748544382847762106786633952487005) },
 113       { SC_(4.0), SC_(16.0),  SC_(22.132051970576036053853444648907108439504682530918) },
 114       { SC_(8.0), SC_(16.0),  SC_(6.1850485247748975008808779795786699492711191898792) },
 115       { SC_(12.0), SC_(16.0),  SC_(3.5694322843488018916484224923627864928705138154372) },
 116       { SC_(16.0), SC_(16.0),  SC_(2.6447371137201451261118187672029372265909501355722) },
 117       { SC_(20.0), SC_(16.0),  SC_(2.1934058398888071720297525592515838555602675797235) },
 118       { SC_(4.0), SC_(20.0),  SC_(41.021743268279206331672552645354782698296383424328) },
 119       { SC_(8.0), SC_(20.0),  SC_(9.3414225299809886395081381945971250426599939097753) },
 120       { SC_(12.0), SC_(20.0),  SC_(4.8253866205826406499959001774187695527272168375992) },
 121       { SC_(16.0), SC_(20.0),  SC_(3.3462305133519485784864062004430532216764447939942) },
 122       { SC_(20.0), SC_(20.0),  SC_(2.6578698872220394617444624241257799193518140676691) },
 123       } };
 124
 125    for (auto row = hypergeometric_pFq_integer_data.begin(); row != hypergeometric_pFq_integer_data.end(); ++row)
 126    {
 127       BOOST_CHECK_CLOSE_FRACTION(boost::math::hypergeometric_pFq_precision({}, { (*row)[0] }, (*row)[1], 20), (*row)[2], tolerance);
 128    }
 129 }
 130
 131 void test_spots_1F1()
 132 {
 133    typedef mp_type T;
 134 #include "hypergeometric_1F1.ipp"
 135
 136    mp_type tolerance = 2e-20;
 137
 138    for (auto row = hypergeometric_1F1.begin(); row != hypergeometric_1F1.end(); ++row)
 139    {
 140       try {
 141          mp_type norm;
 142          mp_type result = boost::math::hypergeometric_pFq_precision({ (*row)[0] }, { (*row)[1] }, (*row)[2], 20);
 143          BOOST_CHECK_CLOSE_FRACTION(result, (*row)[3], tolerance);
 144       }
 145       catch (const boost::math::evaluation_error&) {}
 146    }
 147 }
 148
 149 void test_spots_1F1_b()
 150 {
 151    typedef mp_type T;
 152 #include "hypergeometric_1F1_big.ipp"
 153
 154    mp_type tolerance = 2e-20;
 155
 156    for (auto row = hypergeometric_1F1_big.begin(); row != hypergeometric_1F1_big.end(); ++row)
 157    {
 158       try {
 159          mp_type result = boost::math::hypergeometric_pFq_precision({ (*row)[0] }, { (*row)[1] }, (*row)[2], 20);
 160          BOOST_CHECK_CLOSE_FRACTION(result, (*row)[3], tolerance);
 161       }
 162       catch (const boost::math::evaluation_error&) {}
 163    }
 164 }
 165
 166 void test_spots_2F1()
 167 {
 168    typedef mp_type T;
 169 #include "hypergeometric_2F1.ipp"
 170
 171    mp_type tolerance = 2e-20;
 172
 173    for (auto row = hypergeometric_2F1.begin(); row != hypergeometric_2F1.end(); ++row)
 174    {
 175       try {
 176          mp_type result = boost::math::hypergeometric_pFq_precision({ (*row)[0], (*row)[1] }, { (*row)[2] }, (*row)[3], 20);
 177          BOOST_CHECK_CLOSE_FRACTION(result, (*row)[4], tolerance);
 178       }
 179       catch (const boost::math::evaluation_error&) {}
 180    }
 181 }
 182
 183 void test_spots_0F2()
 184 {
 185    typedef mp_type T;
 186 #include "hypergeometric_0F2.ipp"
 187
 188    mp_type tolerance = 2e-20;
 189
 190    for (auto row = hypergeometric_0F2.begin(); row != hypergeometric_0F2.end(); ++row)
 191    {
 192       try {
 193          T result = boost::math::hypergeometric_pFq_precision({}, { (*row)[0], (*row)[1] }, (*row)[2], 20);
 194          BOOST_CHECK_CLOSE_FRACTION(result, (*row)[3], tolerance);
 195       }
 196       catch (const boost::math::evaluation_error&) {}
 197    }
 198 }
 199
 200 void test_spots_1F2()
 201 {
 202    typedef mp_type T;
 203 #include "hypergeometric_1F2.ipp"
 204
 205    mp_type tolerance = 2e-20;
 206
 207    for (auto row = hypergeometric_1F2.begin(); row != hypergeometric_1F2.end(); ++row)
 208    {
 209       try {
 210          mp_type result = boost::math::hypergeometric_pFq_precision({ (*row)[0] }, { (*row)[1], (*row)[2] }, (*row)[3], 20);
 211          BOOST_CHECK_CLOSE_FRACTION(result, (*row)[4], tolerance);
 212       }
 213       catch (const boost::math::evaluation_error&) {}
 214    }
 215 }
 216
 217 void test_spots_2F2()
 218 {
 219    typedef mp_type T;
 220 #include "hypergeometric_2F2.ipp"
 221
 222    mp_type tolerance = 2e-20;
 223
 224    for (auto row = hypergeometric_2F2.begin(); row != hypergeometric_2F2.end(); ++row)
 225    {
 226       try {
 227          mp_type result = boost::math::hypergeometric_pFq_precision({ (*row)[0], (*row)[1] }, { (*row)[2], (*row)[3] }, (*row)[4], 20);
 228          BOOST_CHECK_CLOSE_FRACTION(result, (*row)[5], tolerance);
 229       }
 230       catch (const boost::math::evaluation_error&) {}
 231    }
 232 }
 233
 234 BOOST_AUTO_TEST_CASE( test_main )
 235 {
 236    test_spots_1F0();
 237    test_spots_0F1();
 238    test_spots_1F1();
 239    test_spots_1F1_b();
 240    test_spots_2F1();
 241    test_spots_0F2();
 242    test_spots_1F2();
 243    test_spots_2F2();
 244
 245    mpfr_free_cache();
 246 }
 247