hyp_1f1_log_big_data.cpp 3.3 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132
  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. #define BOOST_MATH_MAX_SERIES_ITERATION_POLICY 10000000
  6. #include <boost/math/special_functions/hypergeometric_1f1.hpp>
  7. #include <boost/math/special_functions/hypergeometric_pfq.hpp>
  8. #include <boost/math/constants/constants.hpp>
  9. #include <boost/lexical_cast.hpp>
  10. #include <fstream>
  11. #include <map>
  12. #include <boost/math/tools/test_data.hpp>
  13. #include <boost/random.hpp>
  14. #define BOOST_MATH_USE_MPFR
  15. #include "mp_t.hpp"
  16. #include <boost/multiprecision/mpfr.hpp>
  17. using namespace boost::math::tools;
  18. using namespace boost::math;
  19. using namespace std;
  20. using namespace boost::multiprecision;
  21. struct hypergeometric_1f1_gen
  22. {
  23. mp_t operator()(mp_t a1, mp_t a2, mp_t z)
  24. {
  25. mp_t result;
  26. try {
  27. result = (mp_t)log(abs(boost::math::hypergeometric_pFq_precision({ mpfr_float(a1) }, { mpfr_float(a2) }, mpfr_float(z), 70, 25.0)));
  28. std::cout << a1 << " " << a2 << " " << z << " " << result << std::endl;
  29. }
  30. catch (...)
  31. {
  32. throw std::domain_error("");
  33. }
  34. if (fabs(result) > (std::numeric_limits<double>::max)())
  35. {
  36. std::cout << "Rejecting over large value\n";
  37. throw std::domain_error("");
  38. }
  39. if (fabs(result) < 1/1024.0)
  40. {
  41. std::cout << "Rejecting over small value\n";
  42. throw std::domain_error("");
  43. }
  44. return result;
  45. }
  46. };
  47. int main(int, char* [])
  48. {
  49. parameter_info<mp_t> arg1, arg2, arg3;
  50. test_data<mp_t> data;
  51. std::cout << "Welcome.\n"
  52. "This program will generate spot tests for 1F1 (Yeh!!):\n";
  53. std::string line;
  54. //bool cont;
  55. std::vector<mp_t> v;
  56. random_ns::mt19937 rnd;
  57. random_ns::uniform_real_distribution<float> ur_a(0, 1);
  58. mp_t p = ur_a(rnd);
  59. p *= 1e6;
  60. v.push_back(p);
  61. v.push_back(-p);
  62. p = ur_a(rnd);
  63. p *= 1e5;
  64. v.push_back(p);
  65. v.push_back(-p);
  66. p = ur_a(rnd);
  67. p *= 1e4;
  68. v.push_back(p);
  69. v.push_back(-p);
  70. p = ur_a(rnd);
  71. p *= 1e3;
  72. v.push_back(p);
  73. v.push_back(-p);
  74. p = ur_a(rnd);
  75. p *= 1e2;
  76. v.push_back(p);
  77. v.push_back(-p);
  78. p = ur_a(rnd);
  79. p *= 1e-5;
  80. v.push_back(p);
  81. v.push_back(-p);
  82. p = ur_a(rnd);
  83. p *= 1e-12;
  84. v.push_back(p);
  85. v.push_back(-p);
  86. p = ur_a(rnd);
  87. p *= 1e-30;
  88. v.push_back(p);
  89. v.push_back(-p);
  90. for (unsigned i = 0; i < v.size(); ++i)
  91. {
  92. for (unsigned j = 0; j < v.size(); ++j)
  93. {
  94. for (unsigned k = 0; k < v.size(); ++k)
  95. {
  96. std::cout << i << " " << j << " " << k << std::endl;
  97. std::cout << v[i] << " " << (v[j] * 3) / 2 << " " << (v[k] * 5) / 4 << std::endl;
  98. arg1 = make_single_param(v[i]);
  99. arg2 = make_single_param(mp_t((v[j] * 3) / 2));
  100. arg3 = make_single_param(mp_t((v[k] * 5) / 4));
  101. data.insert(hypergeometric_1f1_gen(), arg1, arg2, arg3);
  102. }
  103. }
  104. }
  105. std::cout << "Enter name of test data file [default=hypergeometric_1f1.ipp]";
  106. std::getline(std::cin, line);
  107. boost::algorithm::trim(line);
  108. if(line == "")
  109. line = "hypergeometric_1f1.ipp";
  110. std::ofstream ofs(line.c_str());
  111. ofs << std::scientific << std::setprecision(40);
  112. write_code(ofs, data, line.c_str());
  113. return 0;
  114. }