test_bernoulli_constants.cpp 11 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236
  1. // (C) Copyright John Maddock 2013.
  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. #define BOOST_TEST_MAIN
  6. #define BOOST_MATH_OVERFLOW_ERROR_POLICY ignore_error
  7. #include <boost/test/unit_test.hpp>
  8. #include <boost/math/concepts/real_concept.hpp>
  9. #include <boost/math/special_functions/bernoulli.hpp>
  10. #include <libs/math/test/table_type.hpp>
  11. #include <boost/math/tools/test.hpp>
  12. #include <iostream>
  13. #include <iomanip>
  14. #define SC_(x) static_cast<typename table_type<T>::type>(BOOST_JOIN(x, L))
  15. template <class T>
  16. void test(const char* name)
  17. {
  18. std::cout << "Testing type " << name << ":\n";
  19. static const typename table_type<T>::type data[] =
  20. {
  21. /* First 50 from 2 to 100 inclusive: */
  22. /* TABLE[N[BernoulliB[n], 200], {n,2,100,2}] */
  23. SC_(0.1666666666666666666666666666666666666666),
  24. SC_(-0.0333333333333333333333333333333333333333),
  25. SC_(0.0238095238095238095238095238095238095238),
  26. SC_(-0.0333333333333333333333333333333333333333),
  27. SC_(0.0757575757575757575757575757575757575757),
  28. SC_(-0.2531135531135531135531135531135531135531),
  29. SC_(1.1666666666666666666666666666666666666666),
  30. SC_(-7.0921568627450980392156862745098039215686),
  31. SC_(54.9711779448621553884711779448621553884711),
  32. SC_(-529.1242424242424242424242424242424242424242),
  33. SC_(6192.1231884057971014492753623188405797101449),
  34. SC_(-86580.2531135531135531135531135531135531135531),
  35. SC_(1.4255171666666666666666666666666666666666e6),
  36. SC_(-2.7298231067816091954022988505747126436781e7),
  37. SC_(6.0158087390064236838430386817483591677140e8),
  38. SC_(-1.5116315767092156862745098039215686274509e10),
  39. SC_(4.2961464306116666666666666666666666666666e11),
  40. SC_(-1.3711655205088332772159087948561632772159e13),
  41. SC_(4.8833231897359316666666666666666666666666e14),
  42. SC_(-1.9296579341940068148632668144863266814486e16),
  43. SC_(8.4169304757368261500055370985603543743078e17),
  44. SC_(-4.0338071854059455413076811594202898550724e19),
  45. SC_(2.1150748638081991605601453900709219858156e21),
  46. SC_(-1.2086626522296525934602731193708252531781e23),
  47. SC_(7.5008667460769643668557200757575757575757e24),
  48. SC_(-5.0387781014810689141378930305220125786163e26),
  49. SC_(3.6528776484818123335110430842971177944862e28),
  50. SC_(-2.8498769302450882226269146432910678160919e30),
  51. SC_(2.3865427499683627644645981919219214971751e32),
  52. SC_(-2.1399949257225333665810744765191097392674e34),
  53. SC_(2.0500975723478097569921733095672310251666e36),
  54. SC_(-2.0938005911346378409095185290027970184709e38),
  55. SC_(2.2752696488463515559649260352769264581469e40),
  56. SC_(-2.6257710286239576047303049736158202081449e42),
  57. SC_(3.2125082102718032518204792304264985243521e44),
  58. SC_(-4.1598278166794710913917074495262358936689e46),
  59. SC_(5.6920695482035280023883456219121058644480e48),
  60. SC_(-8.2183629419784575692290653468617333014550e50),
  61. SC_(1.2502904327166993016732339829702895524177e53),
  62. SC_(-2.0015583233248370274925329198813298768724e55),
  63. SC_(3.3674982915364374233396676903338753016219e57),
  64. SC_(-5.9470970503135447718660496844051540840579e59),
  65. SC_(1.1011910323627977559564130790437691604630e62),
  66. SC_(-2.1355259545253501188658385019041065678973e64),
  67. SC_(4.3328896986641192419616613059379206218451e66),
  68. SC_(-9.1885528241669328226200555215501897138960e68),
  69. SC_(2.0346896776329074493455027990220020065975e71),
  70. SC_(-4.7003833958035731078575255535006060654596e73),
  71. SC_(1.1318043445484249270675186257733934267890e76),
  72. SC_(-2.8382249570693706959264156336481764738284e78),
  73. /* next 50 from 102 to 200: */
  74. /* TABLE[N[BernoulliB[n], 200], {n,102,200,2}] */
  75. SC_(7.4064248979678850629750827140920984176879e80),
  76. SC_(-2.0096454802756604483465619672715363186867e83),
  77. SC_(5.6657170050805941445719346030519356961419e85),
  78. SC_(-1.6584511154136216915823713374319912301494e88),
  79. SC_(5.0368859950492377419289421915180154812442e90),
  80. SC_(-1.5861468237658186369363401572966438782740e93),
  81. SC_(5.1756743617545626984073240682507122561240e95),
  82. SC_(-1.7488921840217117339690025877618159145141e98),
  83. SC_(6.1160519994952185255824525264264167780767e100),
  84. SC_(-2.2122776912707834942288323456712932445573e103),
  85. SC_(8.2722776798770969854221062459984595731204e105),
  86. SC_(-3.1958925111415709583591634369180814873526e108),
  87. SC_(1.2750082223387792982310024302926679866957e111),
  88. SC_(-5.2500923086774133899402824624565175446919e113),
  89. SC_(2.2301817894241625209869298198838728143738e116),
  90. SC_(-9.7684521930955204438633513398980239301166e118),
  91. SC_(4.4098361978452954272272622874813169191875e121),
  92. SC_(-2.0508570886464088839729337727583015486456e124),
  93. SC_(9.8214433279791277107572969602097521041491e126),
  94. SC_(-4.8412600798208880508789196709963412761130e129),
  95. SC_(2.4553088801480982609783467404088690399673e132),
  96. SC_(-1.2806926804084747548782513278601785721811e135),
  97. SC_(6.8676167104668581192101888598464400436092e137),
  98. SC_(-3.7846468581969104694978995416379556814489e140),
  99. SC_(2.1426101250665291550871323135148272096660e143),
  100. SC_(-1.2456727137183695007019642961637607219458e146),
  101. SC_(7.4345787551000152543679668394052061311780e148),
  102. SC_(-4.5535795304641704894063333223321274876772e151),
  103. SC_(2.8612112816858868345363847251017232522918e154),
  104. SC_(-1.8437723552033869727688202653628785487541e157),
  105. SC_(1.2181154536221046699501316506599521355817e160),
  106. SC_(-8.2482187185314121548481845729689344730141e162),
  107. SC_(5.7225877937832943329651649814297861591868e165),
  108. SC_(-4.0668530525059104726767969383115865560219e168),
  109. SC_(2.9596092064642050062875269581585187042637e171),
  110. SC_(-2.2049522565189457509031175227344598483637e174),
  111. SC_(1.6812597072889599805831152515136066575446e177),
  112. SC_(-1.3116736213556957648645280635581715300443e180),
  113. SC_(1.0467894009478038082183285392982308964382e183),
  114. SC_(-8.5432893578833707718598254629908277459327e185),
  115. SC_(7.1287821322486542352288406677143822472124e188),
  116. SC_(-6.0802931455535899300084711868647745846198e191),
  117. SC_(5.2996776424849923930094291004324726622848e194),
  118. SC_(-4.7194259168745862644364622901337991110376e197),
  119. SC_(4.2928413791402981089416829654107466904552e200),
  120. SC_(-3.9876744968232207443447765554293879510665e203),
  121. SC_(3.7819780419358882713894418116139332789822e206),
  122. SC_(-3.6614233683681191243685808215119734875519e209),
  123. SC_(3.6176090272372862348855460929891408947754e212),
  124. SC_(-3.6470772645191354362138308865549944904868e215),
  125. };
  126. T tol = boost::math::tools::epsilon<T>() * 20;
  127. for(unsigned i = 1; i <= 100; ++i)
  128. {
  129. T b2n = boost::math::bernoulli_b2n<T>(i);
  130. T t2n = boost::math::tangent_t2n<T>(i);
  131. if((boost::math::isinf)(b2n))
  132. {
  133. if(!(boost::math::isinf)(data[i - 1]))
  134. {
  135. std::cout << "When calculating B2N(" << i << ")\n";
  136. BOOST_ERROR("Got an infinity when one wasn't expected");
  137. }
  138. else if((b2n > 0) != (data[i - 1] > 0))
  139. {
  140. std::cout << "When calculating B2N(" << i << ")\n";
  141. BOOST_ERROR("Sign of infinity was incorrect");
  142. }
  143. }
  144. else if((boost::math::isnan)(b2n))
  145. {
  146. std::cout << "When calculating B2N(" << i << ")\n";
  147. BOOST_ERROR("Result of B2n was a NaN, and that should never happen!");
  148. }
  149. else
  150. {
  151. BOOST_CHECK_CLOSE_FRACTION(b2n, data[i - 1], tol);
  152. }
  153. if(i <= boost::math::max_bernoulli_b2n<T>::value)
  154. {
  155. BOOST_CHECK_EQUAL(b2n, boost::math::unchecked_bernoulli_b2n<T>(i));
  156. }
  157. if((boost::math::isfinite)(t2n) && (t2n < boost::math::tools::max_value<T>()))
  158. {
  159. T p = ldexp(T(1), 2 * i);
  160. int s = i & 1 ? 1 : -1;
  161. p = t2n / (s * p * (p - 1));
  162. p *= 2 * i;
  163. BOOST_CHECK_CLOSE_FRACTION(p, b2n, tol);
  164. }
  165. }
  166. //
  167. // Test consistency of array interface:
  168. //
  169. T bn[boost::math::max_bernoulli_b2n<T>::value + 20];
  170. boost::math::bernoulli_b2n<T>(0, boost::math::max_bernoulli_b2n<T>::value + 20, bn);
  171. for(unsigned i = 0; i < boost::math::max_bernoulli_b2n<T>::value + 20; ++i)
  172. {
  173. BOOST_CHECK_EQUAL(bn[i], boost::math::bernoulli_b2n<T>(i));
  174. }
  175. boost::math::tangent_t2n<T>(0, boost::math::max_bernoulli_b2n<T>::value + 20, bn);
  176. for(unsigned i = 0; i < boost::math::max_bernoulli_b2n<T>::value + 20; ++i)
  177. {
  178. BOOST_CHECK_EQUAL(bn[i], boost::math::tangent_t2n<T>(i));
  179. }
  180. //
  181. // Some spot tests for things that should throw exceptions:
  182. //
  183. static unsigned overflow_index = boost::is_same<T, boost::math::concepts::real_concept>::value ?
  184. boost::math::max_bernoulli_b2n<long double>::value + 5 : boost::math::max_bernoulli_b2n<T>::value + 5;
  185. #ifndef BOOST_NO_EXCEPTIONS
  186. 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);
  187. 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);
  188. #endif
  189. }
  190. void test_real_concept_extra()
  191. {
  192. boost::math::concepts::real_concept tol = boost::math::tools::epsilon<boost::math::concepts::real_concept>() * 20;
  193. for(unsigned i = 0; i <= boost::math::max_bernoulli_b2n<long double>::value; ++i)
  194. {
  195. BOOST_CHECK_CLOSE_FRACTION(static_cast<boost::math::concepts::real_concept>(boost::math::bernoulli_b2n<long double>(i)), boost::math::bernoulli_b2n<boost::math::concepts::real_concept>(i), tol);
  196. }
  197. for(unsigned i = 1; i < 500; ++i)
  198. {
  199. boost::math::concepts::real_concept r = boost::math::bernoulli_b2n<long double>(i + boost::math::max_bernoulli_b2n<long double>::value);
  200. if((i + boost::math::max_bernoulli_b2n<long double>::value) & 1)
  201. {
  202. BOOST_CHECK(r >= boost::math::tools::max_value<boost::math::concepts::real_concept>());
  203. }
  204. else
  205. {
  206. BOOST_CHECK(r <= -boost::math::tools::max_value<boost::math::concepts::real_concept>());
  207. }
  208. }
  209. }
  210. BOOST_AUTO_TEST_CASE( test_main )
  211. {
  212. test<float>("float");
  213. test<double>("double");
  214. #ifndef BOOST_MATH_NO_LONG_DOUBLE_MATH_FUNCTIONS
  215. test<long double>("long double");
  216. test<boost::math::concepts::real_concept>("real_concept");
  217. test_real_concept_extra();
  218. #endif
  219. }