1 // Copyright John Maddock 2015.
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 #ifdef _MSC_VER
7 # pragma warning (disable : 4224)
8 #endif
9
10 #define BOOST_MATH_OVERFLOW_ERROR_POLICY ignore_error
11 #define DISTRIBUTIONS_TEST
12
13 #include <boost/math/distributions.hpp>
14 #include <boost/array.hpp>
15 #include <boost/lexical_cast.hpp>
16 #include "../../test/table_type.hpp"
17 #include "table_helper.hpp"
18 #include "performance.hpp"
19 #include <iostream>
20
21 #ifdef TEST_GSL
22 #include <gsl/gsl_cdf.h>
23 #endif
24
25 class distribution_tester
26 {
27 std::string distro_name;
28 static const double quantiles[19];
29 double sum;
30
31 struct param_info
32 {
33 std::vector<double> params;
34 std::vector<double> x_values;
35 };
36 std::vector<param_info> tests;
sanitize_x(double x)37 double sanitize_x(double x)
38 {
39 if(x > boost::math::tools::max_value<float>() / 2)
40 return boost::math::tools::max_value<float>() / 2;
41 if(x < -boost::math::tools::max_value<float>() / 2)
42 return -boost::math::tools::max_value<float>() / 2;
43 return x;
44 }
45 public:
distribution_tester(const char * name)46 distribution_tester(const char* name) : distro_name(name), sum(0) {}
47
48 template <class F>
add_test_case(F f)49 void add_test_case(F f)
50 {
51 tests.push_back(param_info());
52 for(unsigned i = 0; i < sizeof(quantiles) / sizeof(quantiles[0]); ++i)
53 {
54 tests.back().x_values.push_back(sanitize_x(f(quantiles[i])));
55 }
56 }
57 template <class F>
add_test_case(double p1,F f)58 void add_test_case(double p1, F f)
59 {
60 tests.push_back(param_info());
61 tests.back().params.push_back(p1);
62 for(unsigned i = 0; i < sizeof(quantiles) / sizeof(quantiles[0]); ++i)
63 {
64 tests.back().x_values.push_back(sanitize_x(f(p1, quantiles[i])));
65 }
66 }
67 template <class F>
add_test_case(double p1,double p2,F f)68 void add_test_case(double p1, double p2, F f)
69 {
70 tests.push_back(param_info());
71 tests.back().params.push_back(p1);
72 tests.back().params.push_back(p2);
73 for(unsigned i = 0; i < sizeof(quantiles) / sizeof(quantiles[0]); ++i)
74 {
75 tests.back().x_values.push_back(sanitize_x(f(p1, p2, quantiles[i])));
76 }
77 }
78 template <class F>
add_test_case(double p1,double p2,double p3,F f)79 void add_test_case(double p1, double p2, double p3, F f)
80 {
81 tests.push_back(param_info());
82 tests.back().params.push_back(p1);
83 tests.back().params.push_back(p2);
84 tests.back().params.push_back(p3);
85 for(unsigned i = 0; i < sizeof(quantiles) / sizeof(quantiles[0]); ++i)
86 {
87 tests.back().x_values.push_back(sanitize_x(f(p1, p2, p3, quantiles[i])));
88 }
89 }
90
91 enum
92 {
93 main_table = 1,
94 boost_only_table = 2,
95 both_tables = 3
96 };
97
98 template <class F>
run_timed_tests(F f,std::string sub_name,std::string column,bool p_value=false,int where=main_table)99 void run_timed_tests(F f, std::string sub_name, std::string column, bool p_value = false, int where = main_table)
100 {
101 std::cout << "Testing " << distro_name + " (" + std::string(sub_name) + ")" << " with library " << column << std::endl;
102 try{
103 double t = 0;
104 unsigned repeats = 1;
105 unsigned data_size;
106 do{
107 data_size = 0;
108 stopwatch<boost::chrono::high_resolution_clock> w;
109
110 for(unsigned count = 0; count < repeats; ++count)
111 {
112 for(unsigned i = 0; i < tests.size(); ++i)
113 {
114 for(unsigned j = 0; j < tests[i].x_values.size(); ++j)
115 {
116 if((boost::math::isfinite)(tests[i].x_values[j]))
117 sum += f(tests[i].params, p_value ? quantiles[j] : tests[i].x_values[j]);
118 ++data_size;
119 }
120 }
121 }
122
123 t = boost::chrono::duration_cast<boost::chrono::duration<double>>(w.elapsed()).count();
124 if(t < 0.5)
125 repeats *= 2;
126 } while(t < 0.5);
127
128 static const std::string main_table_name = std::string("Distribution performance comparison with ") + compiler_name() + std::string(" on ") + platform_name();
129 static const std::string boost_table_name = std::string("Distribution performance comparison for different performance options with ") + compiler_name() + std::string(" on ") + platform_name();
130
131 if (where & 1)
132 {
133 report_execution_time(
134 t / data_size,
135 main_table_name,
136 distro_name + " (" + std::string(sub_name) + ")",
137 column);
138 }
139 if (where & 2)
140 {
141 report_execution_time(
142 t / data_size,
143 boost_table_name,
144 distro_name + " (" + std::string(sub_name) + ")",
145 column);
146 }
147 }
148 catch(const std::exception& e)
149 {
150 std::cerr << "Aborting due to exception: " << e.what() << std::endl;
151 std::cerr << "In " << distro_name + " (" + std::string(sub_name) + ")" << std::endl;
152 report_execution_time(
153 (std::numeric_limits<boost::uintmax_t>::max)(),
154 std::string("Distribution performance comparison with ") + compiler_name() + std::string(" on ") + platform_name(),
155 distro_name + " (" + std::string(sub_name) + ")",
156 column);
157 }
158 }
159 };
160
161 const double distribution_tester::quantiles[19] =
162 {
163 0.000001,
164 0.00001,
165 0.0001,
166 0.001,
167 0.01,
168 0.1,
169 0.2,
170 0.3,
171 0.4,
172 0.5,
173 0.6,
174 0.7,
175 0.8,
176 0.9,
177 0.99,
178 0.999,
179 0.9999,
180 0.99999,
181 0.999999
182 };
183
184 template <class D>
185 struct three_param_quantile
186 {
187 template <class T, class U, class V, class X>
operator ()three_param_quantile188 double operator()(T x, U y, V z, X q)const
189 {
190 return quantile(D(x, y, z), q);
191 }
192 };
193
194 template <class D>
195 struct two_param_quantile
196 {
197 template <class T, class U, class V>
operator ()two_param_quantile198 double operator()(T x, U y, V q)const
199 {
200 return quantile(D(x, y), q);
201 }
202 };
203
204 template <class D>
205 struct one_param_quantile
206 {
207 template <class T, class V>
operator ()one_param_quantile208 double operator()(T x, V q)const
209 {
210 return quantile(D(x), q);
211 }
212 };
213
214 template <template <class T, class U> class D>
test_boost_1_param(distribution_tester & tester)215 void test_boost_1_param(distribution_tester& tester)
216 {
217 //
218 // Define some custom policies to test:
219 //
220 typedef boost::math::policies::policy<> default_policy;
221 typedef boost::math::policies::policy<boost::math::policies::promote_double<false> > no_promote_double_policy;
222 typedef boost::math::policies::policy<boost::math::policies::promote_double<false>, boost::math::policies::digits10<10> > no_promote_double_10_digits_policy;
223 typedef boost::math::policies::policy<boost::math::policies::promote_float<false> > no_promote_float_policy;
224
225 tester.run_timed_tests([](const std::vector<double>& v, double x){ return pdf(D<double, default_policy>(v[0]), x); }, "PDF", boost_name(), false, distribution_tester::both_tables);
226 tester.run_timed_tests([](const std::vector<double>& v, double x){ return cdf(D<double, default_policy>(v[0]), x); }, "CDF", boost_name(), false, distribution_tester::both_tables);
227 tester.run_timed_tests([](const std::vector<double>& v, double x){ return quantile(D<double, default_policy>(v[0]), x); }, "quantile", boost_name(), true, distribution_tester::both_tables);
228 if(sizeof(double) != sizeof(long double))
229 {
230 tester.run_timed_tests([](const std::vector<double>& v, double x){ return pdf(D<double, no_promote_double_policy>(v[0]), x); }, "PDF", "Boost[br]promote_double<false>", false, distribution_tester::both_tables);
231 tester.run_timed_tests([](const std::vector<double>& v, double x){ return cdf(D<double, no_promote_double_policy>(v[0]), x); }, "CDF", "Boost[br]promote_double<false>", false, distribution_tester::both_tables);
232 tester.run_timed_tests([](const std::vector<double>& v, double x){ return quantile(D<double, no_promote_double_policy>(v[0]), x); }, "quantile", "Boost[br]promote_double<false>", true, distribution_tester::both_tables);
233 }
234 tester.run_timed_tests([](const std::vector<double>& v, double x){ return pdf(D<double, no_promote_double_10_digits_policy>(v[0]), x); }, "PDF", "Boost[br]promote_double<false>[br]digits10<10>", false, distribution_tester::boost_only_table);
235 tester.run_timed_tests([](const std::vector<double>& v, double x){ return cdf(D<double, no_promote_double_10_digits_policy>(v[0]), x); }, "CDF", "Boost[br]promote_double<false>[br]digits10<10>", false, distribution_tester::boost_only_table);
236 tester.run_timed_tests([](const std::vector<double>& v, double x){ return quantile(D<double, no_promote_double_10_digits_policy>(v[0]), x); }, "quantile", "Boost[br]promote_double<false>[br]digits10<10>", true, distribution_tester::boost_only_table);
237
238 tester.run_timed_tests([](const std::vector<double>& v, double x){ return pdf(D<float, no_promote_float_policy>(static_cast<float>(v[0])), static_cast<float>(x)); }, "PDF", "Boost[br]float[br]promote_float<false>", false, distribution_tester::boost_only_table);
239 tester.run_timed_tests([](const std::vector<double>& v, double x){ return cdf(D<float, no_promote_float_policy>(static_cast<float>(v[0])), static_cast<float>(x)); }, "CDF", "Boost[br]float[br]promote_float<false>", false, distribution_tester::boost_only_table);
240 tester.run_timed_tests([](const std::vector<double>& v, double x){ return quantile(D<float, no_promote_float_policy>(static_cast<float>(v[0])), static_cast<float>(x)); }, "quantile", "Boost[br]float[br]promote_float<false>", true, distribution_tester::boost_only_table);
241 }
242
243 template <template <class T, class U> class D>
test_boost_2_param(distribution_tester & tester)244 void test_boost_2_param(distribution_tester& tester)
245 {
246 //
247 // Define some custom policies to test:
248 //
249 typedef boost::math::policies::policy<> default_policy;
250 typedef boost::math::policies::policy<boost::math::policies::promote_double<false> > no_promote_double_policy;
251 typedef boost::math::policies::policy<boost::math::policies::promote_double<false>, boost::math::policies::digits10<10> > no_promote_double_10_digits_policy;
252 typedef boost::math::policies::policy<boost::math::policies::promote_float<false> > no_promote_float_policy;
253
254 tester.run_timed_tests([](const std::vector<double>& v, double x){ return pdf(D<double, default_policy>(v[0], v[1]), x); }, "PDF", boost_name(), false, distribution_tester::both_tables);
255 tester.run_timed_tests([](const std::vector<double>& v, double x){ return cdf(D<double, default_policy>(v[0], v[1]), x); }, "CDF", boost_name(), false, distribution_tester::both_tables);
256 tester.run_timed_tests([](const std::vector<double>& v, double x){ return quantile(D<double, default_policy>(v[0], v[1]), x); }, "quantile", boost_name(), true, distribution_tester::both_tables);
257 if(sizeof(double) != sizeof(long double))
258 {
259 tester.run_timed_tests([](const std::vector<double>& v, double x){ return pdf(D<double, no_promote_double_policy>(v[0], v[1]), x); }, "PDF", "Boost[br]promote_double<false>", false, distribution_tester::both_tables);
260 tester.run_timed_tests([](const std::vector<double>& v, double x){ return cdf(D<double, no_promote_double_policy>(v[0], v[1]), x); }, "CDF", "Boost[br]promote_double<false>", false, distribution_tester::both_tables);
261 tester.run_timed_tests([](const std::vector<double>& v, double x){ return quantile(D<double, no_promote_double_policy>(v[0], v[1]), x); }, "quantile", "Boost[br]promote_double<false>", true, distribution_tester::both_tables);
262 }
263 tester.run_timed_tests([](const std::vector<double>& v, double x){ return pdf(D<double, no_promote_double_10_digits_policy>(v[0], v[1]), x); }, "PDF", "Boost[br]promote_double<false>[br]digits10<10>", false, distribution_tester::boost_only_table);
264 tester.run_timed_tests([](const std::vector<double>& v, double x){ return cdf(D<double, no_promote_double_10_digits_policy>(v[0], v[1]), x); }, "CDF", "Boost[br]promote_double<false>[br]digits10<10>", false, distribution_tester::boost_only_table);
265 tester.run_timed_tests([](const std::vector<double>& v, double x){ return quantile(D<double, no_promote_double_10_digits_policy>(v[0], v[1]), x); }, "quantile", "Boost[br]promote_double<false>[br]digits10<10>", true, distribution_tester::boost_only_table);
266
267 tester.run_timed_tests([](const std::vector<double>& v, double x){ return pdf(D<float, no_promote_float_policy>(static_cast<float>(v[0]), static_cast<float>(v[1])), static_cast<float>(x)); }, "PDF", "Boost[br]float[br]promote_float<false>", false, distribution_tester::boost_only_table);
268 tester.run_timed_tests([](const std::vector<double>& v, double x){ return cdf(D<float, no_promote_float_policy>(static_cast<float>(v[0]), static_cast<float>(v[1])), static_cast<float>(x)); }, "CDF", "Boost[br]float[br]promote_float<false>", false, distribution_tester::boost_only_table);
269 tester.run_timed_tests([](const std::vector<double>& v, double x){ return quantile(D<float, no_promote_float_policy>(static_cast<float>(v[0]), static_cast<float>(v[1])), static_cast<float>(x)); }, "quantile", "Boost[br]float[br]promote_float<false>", true, distribution_tester::boost_only_table);
270 }
271
272 template <template <class T, class U> class D>
test_boost_3_param(distribution_tester & tester)273 void test_boost_3_param(distribution_tester& tester)
274 {
275 //
276 // Define some custom policies to test:
277 //
278 typedef boost::math::policies::policy<> default_policy;
279 typedef boost::math::policies::policy<boost::math::policies::promote_double<false> > no_promote_double_policy;
280 typedef boost::math::policies::policy<boost::math::policies::promote_double<false>, boost::math::policies::digits10<10> > no_promote_double_10_digits_policy;
281 typedef boost::math::policies::policy<boost::math::policies::promote_float<false> > no_promote_float_policy;
282
283 tester.run_timed_tests([](const std::vector<double>& v, double x){ return pdf(D<double, default_policy>(v[0], v[1], v[2]), x); }, "PDF", boost_name(), false, distribution_tester::both_tables);
284 tester.run_timed_tests([](const std::vector<double>& v, double x){ return cdf(D<double, default_policy>(v[0], v[1], v[2]), x); }, "CDF", boost_name(), false, distribution_tester::both_tables);
285 tester.run_timed_tests([](const std::vector<double>& v, double x){ return quantile(D<double, default_policy>(v[0], v[1], v[2]), x); }, "quantile", boost_name(), true, distribution_tester::both_tables);
286 if(sizeof(double) != sizeof(long double))
287 {
288 tester.run_timed_tests([](const std::vector<double>& v, double x){ return pdf(D<double, no_promote_double_policy>(v[0], v[1], v[2]), x); }, "PDF", "Boost[br]promote_double<false>", false, distribution_tester::both_tables);
289 tester.run_timed_tests([](const std::vector<double>& v, double x){ return cdf(D<double, no_promote_double_policy>(v[0], v[1], v[2]), x); }, "CDF", "Boost[br]promote_double<false>", false, distribution_tester::both_tables);
290 tester.run_timed_tests([](const std::vector<double>& v, double x){ return quantile(D<double, no_promote_double_policy>(v[0], v[1], v[2]), x); }, "quantile", "Boost[br]promote_double<false>", true, distribution_tester::both_tables);
291 }
292 tester.run_timed_tests([](const std::vector<double>& v, double x){ return pdf(D<double, no_promote_double_10_digits_policy>(v[0], v[1], v[2]), x); }, "PDF", "Boost[br]promote_double<false>[br]digits10<10>", false, distribution_tester::boost_only_table);
293 tester.run_timed_tests([](const std::vector<double>& v, double x){ return cdf(D<double, no_promote_double_10_digits_policy>(v[0], v[1], v[2]), x); }, "CDF", "Boost[br]promote_double<false>[br]digits10<10>", false, distribution_tester::boost_only_table);
294 tester.run_timed_tests([](const std::vector<double>& v, double x){ return quantile(D<double, no_promote_double_10_digits_policy>(v[0], v[1], v[2]), x); }, "quantile", "Boost[br]promote_double<false>[br]digits10<10>", true, distribution_tester::boost_only_table);
295
296 tester.run_timed_tests([](const std::vector<double>& v, double x){ return pdf(D<float, no_promote_float_policy>(static_cast<float>(v[0]), static_cast<float>(v[1]), static_cast<float>(v[2])), static_cast<float>(x)); }, "PDF", "Boost[br]float[br]promote_float<false>", false, distribution_tester::boost_only_table);
297 tester.run_timed_tests([](const std::vector<double>& v, double x){ return cdf(D<float, no_promote_float_policy>(static_cast<float>(v[0]), static_cast<float>(v[1]), static_cast<float>(v[2])), static_cast<float>(x)); }, "CDF", "Boost[br]float[br]promote_float<false>", false, distribution_tester::boost_only_table);
298 tester.run_timed_tests([](const std::vector<double>& v, double x){ return quantile(D<float, no_promote_float_policy>(static_cast<float>(v[0]), static_cast<float>(v[1]), static_cast<float>(v[2])), static_cast<float>(x)); }, "quantile", "Boost[br]float[br]promote_float<false>", true, distribution_tester::boost_only_table);
299 }
300
main()301 int main()
302 {
303 try {
304 //
305 // Normal:
306 //
307 distribution_tester n("Normal");
308 n.add_test_case(0, 1, two_param_quantile<boost::math::normal_distribution<> >());
309 n.add_test_case(20, 20, two_param_quantile<boost::math::normal_distribution<> >());
310 n.add_test_case(-20, 0.0125, two_param_quantile<boost::math::normal_distribution<> >());
311
312 test_boost_2_param<boost::math::normal_distribution>(n);
313
314 distribution_tester arcsine("ArcSine");
315 arcsine.add_test_case(0, 1, two_param_quantile<boost::math::arcsine_distribution<> >());
316 arcsine.add_test_case(20, 500, two_param_quantile<boost::math::arcsine_distribution<> >());
317 arcsine.add_test_case(-20, 100000, two_param_quantile<boost::math::arcsine_distribution<> >());
318
319 test_boost_2_param<boost::math::arcsine_distribution>(arcsine);
320
321 distribution_tester beta("Beta");
322 beta.add_test_case(1, 4, two_param_quantile<boost::math::beta_distribution<> >());
323 beta.add_test_case(20, 500, two_param_quantile<boost::math::beta_distribution<> >());
324 beta.add_test_case(0.1, 0.01, two_param_quantile<boost::math::beta_distribution<> >());
325
326 test_boost_2_param<boost::math::beta_distribution>(beta);
327
328 distribution_tester binomial("Binomial");
329 binomial.add_test_case(5, 0.125, two_param_quantile<boost::math::binomial_distribution<> >());
330 binomial.add_test_case(200, 0.75, two_param_quantile<boost::math::binomial_distribution<> >());
331 binomial.add_test_case(2000, 0.5, two_param_quantile<boost::math::binomial_distribution<> >());
332 binomial.add_test_case(20000, 0.001, two_param_quantile<boost::math::binomial_distribution<> >());
333 binomial.add_test_case(200000, 0.99, two_param_quantile<boost::math::binomial_distribution<> >());
334
335 test_boost_2_param<boost::math::binomial_distribution>(binomial);
336
337 distribution_tester cauchy("Cauchy");
338 cauchy.add_test_case(0, 1, two_param_quantile<boost::math::cauchy_distribution<> >());
339 cauchy.add_test_case(20, 20, two_param_quantile<boost::math::cauchy_distribution<> >());
340 cauchy.add_test_case(-20, 0.0125, two_param_quantile<boost::math::cauchy_distribution<> >());
341
342 test_boost_2_param<boost::math::cauchy_distribution>(cauchy);
343
344 distribution_tester chi_squared("ChiSquared");
345 chi_squared.add_test_case(3, one_param_quantile<boost::math::chi_squared_distribution<> >());
346 chi_squared.add_test_case(20, one_param_quantile<boost::math::chi_squared_distribution<> >());
347 chi_squared.add_test_case(200, one_param_quantile<boost::math::chi_squared_distribution<> >());
348 chi_squared.add_test_case(2000, one_param_quantile<boost::math::chi_squared_distribution<> >());
349 chi_squared.add_test_case(20000, one_param_quantile<boost::math::chi_squared_distribution<> >());
350 chi_squared.add_test_case(200000, one_param_quantile<boost::math::chi_squared_distribution<> >());
351
352 test_boost_1_param<boost::math::chi_squared_distribution>(chi_squared);
353
354 distribution_tester exponential("Exponential");
355 exponential.add_test_case(0.001, one_param_quantile<boost::math::exponential_distribution<> >());
356 exponential.add_test_case(0.01, one_param_quantile<boost::math::exponential_distribution<> >());
357 exponential.add_test_case(0.1, one_param_quantile<boost::math::exponential_distribution<> >());
358 exponential.add_test_case(1, one_param_quantile<boost::math::exponential_distribution<> >());
359 exponential.add_test_case(10, one_param_quantile<boost::math::exponential_distribution<> >());
360 exponential.add_test_case(100, one_param_quantile<boost::math::exponential_distribution<> >());
361 exponential.add_test_case(1000, one_param_quantile<boost::math::exponential_distribution<> >());
362
363 test_boost_1_param<boost::math::exponential_distribution>(exponential);
364
365 distribution_tester extreme_value("ExtremeValue");
366 extreme_value.add_test_case(0, 1, two_param_quantile<boost::math::extreme_value_distribution<> >());
367 extreme_value.add_test_case(20, 20, two_param_quantile<boost::math::extreme_value_distribution<> >());
368 extreme_value.add_test_case(-20, 0.0125, two_param_quantile<boost::math::extreme_value_distribution<> >());
369
370 test_boost_2_param<boost::math::extreme_value_distribution>(extreme_value);
371
372 distribution_tester fisher("F");
373 for (unsigned i = 2; i <= 200000; i *= 10)
374 {
375 for (unsigned j = 2; j <= 200000; j *= 10)
376 {
377 fisher.add_test_case(i, j, two_param_quantile<boost::math::fisher_f_distribution<> >());
378 }
379 }
380 test_boost_2_param<boost::math::fisher_f_distribution>(fisher);
381
382 distribution_tester gamma("Gamma");
383 gamma.add_test_case(0.1, 1, two_param_quantile<boost::math::gamma_distribution<> >());
384 gamma.add_test_case(20, 20, two_param_quantile<boost::math::gamma_distribution<> >());
385 gamma.add_test_case(200, 0.0125, two_param_quantile<boost::math::gamma_distribution<> >());
386 gamma.add_test_case(2000, 500, two_param_quantile<boost::math::gamma_distribution<> >());
387
388 test_boost_2_param<boost::math::gamma_distribution>(gamma);
389
390 distribution_tester geometric("Geometric");
391 geometric.add_test_case(0.001, one_param_quantile<boost::math::geometric_distribution<> >());
392 geometric.add_test_case(0.01, one_param_quantile<boost::math::geometric_distribution<> >());
393 geometric.add_test_case(0.1, one_param_quantile<boost::math::geometric_distribution<> >());
394 geometric.add_test_case(0.5, one_param_quantile<boost::math::geometric_distribution<> >());
395 geometric.add_test_case(0.9, one_param_quantile<boost::math::geometric_distribution<> >());
396 geometric.add_test_case(0.99, one_param_quantile<boost::math::geometric_distribution<> >());
397 geometric.add_test_case(0.999, one_param_quantile<boost::math::geometric_distribution<> >());
398
399 test_boost_1_param<boost::math::geometric_distribution>(geometric);
400
401 distribution_tester hypergeometric("Hypergeometric");
402 hypergeometric.add_test_case(10, 5, 100, three_param_quantile<boost::math::hypergeometric_distribution<> >());
403 hypergeometric.add_test_case(50, 75, 100, three_param_quantile<boost::math::hypergeometric_distribution<> >());
404 hypergeometric.add_test_case(30, 20, 100, three_param_quantile<boost::math::hypergeometric_distribution<> >());
405 hypergeometric.add_test_case(100, 50, 1000000, three_param_quantile<boost::math::hypergeometric_distribution<> >());
406 hypergeometric.add_test_case(500000, 3000, 1000000, three_param_quantile<boost::math::hypergeometric_distribution<> >());
407 hypergeometric.add_test_case(20000, 800000, 1000000, three_param_quantile<boost::math::hypergeometric_distribution<> >());
408 hypergeometric.add_test_case(100, 5, 1000, three_param_quantile<boost::math::hypergeometric_distribution<> >());
409 hypergeometric.add_test_case(500, 50, 1000, three_param_quantile<boost::math::hypergeometric_distribution<> >());
410 hypergeometric.add_test_case(2, 25, 1000, three_param_quantile<boost::math::hypergeometric_distribution<> >());
411 hypergeometric.add_test_case(1, 5, 1000, three_param_quantile<boost::math::hypergeometric_distribution<> >());
412 hypergeometric.add_test_case(100, 500, 1000, three_param_quantile<boost::math::hypergeometric_distribution<> >());
413
414 test_boost_3_param<boost::math::hypergeometric_distribution>(hypergeometric);
415
416 distribution_tester inverse_chi_squared("InverseChiSquared");
417 inverse_chi_squared.add_test_case(5, 0.125, two_param_quantile<boost::math::inverse_chi_squared_distribution<> >());
418 inverse_chi_squared.add_test_case(200, 0.75, two_param_quantile<boost::math::inverse_chi_squared_distribution<> >());
419 inverse_chi_squared.add_test_case(2000, 1, two_param_quantile<boost::math::inverse_chi_squared_distribution<> >());
420 inverse_chi_squared.add_test_case(20000, 10, two_param_quantile<boost::math::inverse_chi_squared_distribution<> >());
421 inverse_chi_squared.add_test_case(200000, 100, two_param_quantile<boost::math::inverse_chi_squared_distribution<> >());
422
423 test_boost_2_param<boost::math::inverse_chi_squared_distribution>(inverse_chi_squared);
424
425 distribution_tester inverse_gamma("InverseGamma");
426 inverse_gamma.add_test_case(0.1, 1, two_param_quantile<boost::math::inverse_gamma_distribution<> >());
427 inverse_gamma.add_test_case(20, 20, two_param_quantile<boost::math::inverse_gamma_distribution<> >());
428 inverse_gamma.add_test_case(200, 0.0125, two_param_quantile<boost::math::inverse_gamma_distribution<> >());
429 inverse_gamma.add_test_case(2000, 500, two_param_quantile<boost::math::inverse_gamma_distribution<> >());
430
431 test_boost_2_param<boost::math::inverse_gamma_distribution>(inverse_gamma);
432
433 distribution_tester inverse_gaussian("InverseGaussian");
434 inverse_gaussian.add_test_case(0.001, 1, two_param_quantile<boost::math::inverse_gaussian_distribution<> >());
435 inverse_gaussian.add_test_case(20, 20, two_param_quantile<boost::math::inverse_gaussian_distribution<> >());
436
437 test_boost_2_param<boost::math::inverse_gaussian_distribution>(inverse_gaussian);
438
439 distribution_tester laplace("Laplace");
440 laplace.add_test_case(0, 1, two_param_quantile<boost::math::laplace_distribution<> >());
441 laplace.add_test_case(20, 20, two_param_quantile<boost::math::laplace_distribution<> >());
442 laplace.add_test_case(-20, 0.0125, two_param_quantile<boost::math::laplace_distribution<> >());
443
444 test_boost_2_param<boost::math::laplace_distribution>(laplace);
445
446 distribution_tester logistic("Logistic");
447 logistic.add_test_case(0, 1, two_param_quantile<boost::math::logistic_distribution<> >());
448 logistic.add_test_case(20, 20, two_param_quantile<boost::math::logistic_distribution<> >());
449 logistic.add_test_case(-20, 0.0125, two_param_quantile<boost::math::logistic_distribution<> >());
450
451 test_boost_2_param<boost::math::logistic_distribution>(logistic);
452
453 distribution_tester lognormal("LogNormal");
454 lognormal.add_test_case(0, 1, two_param_quantile<boost::math::lognormal_distribution<> >());
455 lognormal.add_test_case(20, 20, two_param_quantile<boost::math::lognormal_distribution<> >());
456 lognormal.add_test_case(-20, 0.0125, two_param_quantile<boost::math::lognormal_distribution<> >());
457
458 test_boost_2_param<boost::math::lognormal_distribution>(lognormal);
459
460 distribution_tester negative_binomial("NegativeBinomial");
461 negative_binomial.add_test_case(5, 0.125, two_param_quantile<boost::math::negative_binomial_distribution<> >());
462 negative_binomial.add_test_case(200, 0.75, two_param_quantile<boost::math::negative_binomial_distribution<> >());
463 negative_binomial.add_test_case(2000, 0.001, two_param_quantile<boost::math::negative_binomial_distribution<> >());
464 negative_binomial.add_test_case(20000, 0.5, two_param_quantile<boost::math::negative_binomial_distribution<> >());
465 negative_binomial.add_test_case(200000, 0.99, two_param_quantile<boost::math::negative_binomial_distribution<> >());
466
467 test_boost_2_param<boost::math::negative_binomial_distribution>(negative_binomial);
468
469 distribution_tester non_central_beta("NonCentralBeta");
470 non_central_beta.add_test_case(2, 5, 2.1, three_param_quantile<boost::math::non_central_beta_distribution<> >());
471 non_central_beta.add_test_case(0.25, 0.01, 20, three_param_quantile<boost::math::non_central_beta_distribution<> >());
472 non_central_beta.add_test_case(20, 3, 30, three_param_quantile<boost::math::non_central_beta_distribution<> >());
473 non_central_beta.add_test_case(100, 200, 400, three_param_quantile<boost::math::non_central_beta_distribution<> >());
474 non_central_beta.add_test_case(100, 0.25, 20, three_param_quantile<boost::math::non_central_beta_distribution<> >());
475
476 test_boost_3_param<boost::math::non_central_beta_distribution>(non_central_beta);
477
478 distribution_tester non_central_chi_squared("NonCentralChiSquared");
479 non_central_chi_squared.add_test_case(5, 0.5, two_param_quantile<boost::math::non_central_chi_squared_distribution<> >());
480 non_central_chi_squared.add_test_case(200, 2, two_param_quantile<boost::math::non_central_chi_squared_distribution<> >());
481 non_central_chi_squared.add_test_case(2000, 20, two_param_quantile<boost::math::non_central_chi_squared_distribution<> >());
482 non_central_chi_squared.add_test_case(20000, 10, two_param_quantile<boost::math::non_central_chi_squared_distribution<> >());
483 non_central_chi_squared.add_test_case(200000, 50, two_param_quantile<boost::math::non_central_chi_squared_distribution<> >());
484
485 test_boost_2_param<boost::math::non_central_chi_squared_distribution>(non_central_chi_squared);
486
487 distribution_tester non_central_f("NonCentralF");
488 non_central_f.add_test_case(20, 20, 3, three_param_quantile<boost::math::non_central_f_distribution<> >());
489 non_central_f.add_test_case(20, 50, 20, three_param_quantile<boost::math::non_central_f_distribution<> >());
490 non_central_f.add_test_case(100, 20, 30, three_param_quantile<boost::math::non_central_f_distribution<> >());
491 non_central_f.add_test_case(100, 200, 100, three_param_quantile<boost::math::non_central_f_distribution<> >());
492 non_central_f.add_test_case(1000, 100000, 20, three_param_quantile<boost::math::non_central_f_distribution<> >());
493
494 test_boost_3_param<boost::math::non_central_f_distribution>(non_central_f);
495
496 distribution_tester non_central_t("NonCentralT");
497 non_central_t.add_test_case(5, 0.5, two_param_quantile<boost::math::non_central_t_distribution<> >());
498 non_central_t.add_test_case(200, 2, two_param_quantile<boost::math::non_central_t_distribution<> >());
499 non_central_t.add_test_case(2000, 20, two_param_quantile<boost::math::non_central_t_distribution<> >());
500 non_central_t.add_test_case(20000, 10, two_param_quantile<boost::math::non_central_t_distribution<> >());
501 non_central_t.add_test_case(200000, 50, two_param_quantile<boost::math::non_central_t_distribution<> >());
502
503 test_boost_2_param<boost::math::non_central_t_distribution>(non_central_t);
504
505 distribution_tester pareto("Pareto");
506 pareto.add_test_case(0.1, 1, two_param_quantile<boost::math::pareto_distribution<> >());
507 pareto.add_test_case(20, 20, two_param_quantile<boost::math::pareto_distribution<> >());
508 pareto.add_test_case(200, 0.0125, two_param_quantile<boost::math::pareto_distribution<> >());
509 pareto.add_test_case(2000, 500, two_param_quantile<boost::math::pareto_distribution<> >());
510
511 test_boost_2_param<boost::math::pareto_distribution>(pareto);
512
513 distribution_tester poisson("Poisson");
514 poisson.add_test_case(0.001, one_param_quantile<boost::math::poisson_distribution<> >());
515 poisson.add_test_case(0.01, one_param_quantile<boost::math::poisson_distribution<> >());
516 poisson.add_test_case(0.1, one_param_quantile<boost::math::poisson_distribution<> >());
517 poisson.add_test_case(1, one_param_quantile<boost::math::poisson_distribution<> >());
518 poisson.add_test_case(10, one_param_quantile<boost::math::poisson_distribution<> >());
519 poisson.add_test_case(100, one_param_quantile<boost::math::poisson_distribution<> >());
520 poisson.add_test_case(1000, one_param_quantile<boost::math::poisson_distribution<> >());
521
522 test_boost_1_param<boost::math::poisson_distribution>(poisson);
523
524 distribution_tester rayleigh("Rayleigh");
525 rayleigh.add_test_case(0.001, one_param_quantile<boost::math::rayleigh_distribution<> >());
526 rayleigh.add_test_case(0.01, one_param_quantile<boost::math::rayleigh_distribution<> >());
527 rayleigh.add_test_case(0.1, one_param_quantile<boost::math::rayleigh_distribution<> >());
528 rayleigh.add_test_case(1, one_param_quantile<boost::math::rayleigh_distribution<> >());
529 rayleigh.add_test_case(10, one_param_quantile<boost::math::rayleigh_distribution<> >());
530 rayleigh.add_test_case(100, one_param_quantile<boost::math::rayleigh_distribution<> >());
531 rayleigh.add_test_case(1000, one_param_quantile<boost::math::rayleigh_distribution<> >());
532
533 test_boost_1_param<boost::math::rayleigh_distribution>(rayleigh);
534
535 distribution_tester skew_norm("SkewNormal");
536 skew_norm.add_test_case(0, 1, 0.1, three_param_quantile<boost::math::skew_normal_distribution<> >());
537 skew_norm.add_test_case(20, 20, 30, three_param_quantile<boost::math::skew_normal_distribution<> >());
538 skew_norm.add_test_case(-20, 0.0125, 10, three_param_quantile<boost::math::skew_normal_distribution<> >());
539
540 test_boost_3_param<boost::math::skew_normal_distribution>(skew_norm);
541
542 distribution_tester students_t("StudentsT");
543 students_t.add_test_case(3, one_param_quantile<boost::math::students_t_distribution<> >());
544 students_t.add_test_case(20, one_param_quantile<boost::math::students_t_distribution<> >());
545 students_t.add_test_case(200, one_param_quantile<boost::math::students_t_distribution<> >());
546 students_t.add_test_case(2000, one_param_quantile<boost::math::students_t_distribution<> >());
547 students_t.add_test_case(20000, one_param_quantile<boost::math::students_t_distribution<> >());
548 students_t.add_test_case(200000, one_param_quantile<boost::math::students_t_distribution<> >());
549
550 test_boost_1_param<boost::math::students_t_distribution>(students_t);
551
552 distribution_tester weibull("Weibull");
553 weibull.add_test_case(0.1, 1, two_param_quantile<boost::math::weibull_distribution<> >());
554 weibull.add_test_case(20, 20, two_param_quantile<boost::math::weibull_distribution<> >());
555 weibull.add_test_case(200, 0.0125, two_param_quantile<boost::math::weibull_distribution<> >());
556 weibull.add_test_case(2000, 500, two_param_quantile<boost::math::weibull_distribution<> >());
557
558 test_boost_2_param<boost::math::weibull_distribution>(weibull);
559
560 #ifdef TEST_GSL
561 // normal, note no location param
562 n.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_gaussian_P(x, v[1]); }, "CDF", "GSL");
563 n.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_gaussian_Pinv(x, v[1]); }, "quantile", "GSL", true);
564 // exponential:
565 exponential.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_exponential_P(x, 1 / v[0]); }, "CDF", "GSL");
566 exponential.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_exponential_Pinv(x, 1 / v[0]); }, "quantile", "GSL", true);
567 // laplace, note no location param:
568 laplace.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_laplace_P(x, v[1]); }, "CDF", "GSL");
569 laplace.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_laplace_Pinv(x, v[1]); }, "quantile", "GSL", true);
570 // cauchy, note no location param:
571 cauchy.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_cauchy_P(x, v[1]); }, "CDF", "GSL");
572 cauchy.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_cauchy_Pinv(x, v[1]); }, "quantile", "GSL", true);
573 // rayleigh:
574 rayleigh.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_rayleigh_P(x, v[0]); }, "CDF", "GSL");
575 rayleigh.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_rayleigh_Pinv(x, v[0]); }, "quantile", "GSL", true);
576 // gamma:
577 gamma.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_gamma_P(x, v[0], v[1]); }, "CDF", "GSL");
578 gamma.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_gamma_Pinv(x, v[0], v[1]); }, "quantile", "GSL", true);
579 // lognormal:
580 lognormal.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_lognormal_P(x, v[0], v[1]); }, "CDF", "GSL");
581 lognormal.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_lognormal_Pinv(x, v[0], v[1]); }, "quantile", "GSL", true);
582 // chi squared:
583 chi_squared.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_chisq_P(x, v[0]); }, "CDF", "GSL");
584 chi_squared.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_chisq_Pinv(x, v[0]); }, "quantile", "GSL", true);
585 // F:
586 fisher.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_fdist_P(x, v[0], v[1]); }, "CDF", "GSL");
587 fisher.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_fdist_Pinv(x, v[0], v[1]); }, "quantile", "GSL", true);
588 // T:
589 students_t.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_tdist_P(x, v[0]); }, "CDF", "GSL");
590 students_t.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_tdist_Pinv(x, v[0]); }, "quantile", "GSL", true);
591 // beta:
592 beta.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_beta_P(x, v[0], v[1]); }, "CDF", "GSL");
593 beta.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_beta_Pinv(x, v[0], v[1]); }, "quantile", "GSL", true);
594 // logistic, note no location param
595 logistic.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_logistic_P(x, v[1]); }, "CDF", "GSL");
596 logistic.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_logistic_Pinv(x, v[1]); }, "quantile", "GSL", true);
597 // pareto:
598 pareto.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_pareto_P(x, v[1], v[0]); }, "CDF", "GSL");
599 pareto.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_pareto_Pinv(x, v[1], v[0]); }, "quantile", "GSL", true);
600 // weibull:
601 weibull.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_weibull_P(x, v[1], v[0]); }, "CDF", "GSL");
602 weibull.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_weibull_Pinv(x, v[1], v[0]); }, "quantile", "GSL", true);
603 // poisson:
604 poisson.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_poisson_P(x, v[0]); }, "CDF", "GSL");
605 //poisson.run_timed_tests([](const std::vector<double>& v, double x){ return gsl_cdf_poisson_Pinv(x, v[0]); }, "quantile", "GSL", true);
606 // binomial:
607 binomial.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_binomial_P(x, v[1], v[0]); }, "CDF", "GSL");
608 //binomial.run_timed_tests([](const std::vector<double>& v, double x){ return gsl_cdf_binomial_Pinv(x, v[1], v[0]); }, "quantile", "GSL", true);
609 // negative_binomial:
610 negative_binomial.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_negative_binomial_P(x, v[1], v[0]); }, "CDF", "GSL");
611 //negative_binomial.run_timed_tests([](const std::vector<double>& v, double x){ return gsl_cdf_negative_binomial_Pinv(x, v[1], v[0]); }, "quantile", "GSL", true);
612 // geometric:
613 geometric.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_geometric_P(x + 1, v[0]); }, "CDF", "GSL");
614 //geometric.run_timed_tests([](const std::vector<double>& v, double x){ return gsl_cdf_geometric_Pinv(x, v[0]) - 1; }, "quantile", "GSL", true);
615 // hypergeometric:
616 hypergeometric.run_timed_tests([](const std::vector<double>& v, double x) { return gsl_cdf_hypergeometric_P(x, v[0], v[2] - v[0], v[1]); }, "CDF", "GSL");
617 //hypergeometric.run_timed_tests([](const std::vector<double>& v, double x){ return gsl_cdf_hypergeometric_Pinv(x, v[0], v[2] - v[0], v[1]); }, "quantile", "GSL", true);
618 #endif
619
620 #ifdef TEST_RMATH
621 // beta
622 beta.run_timed_tests([](const std::vector<double>& v, double x) { return dbeta(x, v[0], v[1], 0); }, "PDF", "Rmath " R_VERSION_STRING);
623 beta.run_timed_tests([](const std::vector<double>& v, double x) { return pbeta(x, v[0], v[1], 1, 0); }, "CDF", "Rmath " R_VERSION_STRING);
624 beta.run_timed_tests([](const std::vector<double>& v, double x) { return qbeta(x, v[0], v[1], 1, 0); }, "quantile", "Rmath " R_VERSION_STRING, true);
625 // non-central beta
626 non_central_beta.run_timed_tests([](const std::vector<double>& v, double x) { return dnbeta(x, v[0], v[1], v[2], 0); }, "PDF", "Rmath " R_VERSION_STRING);
627 non_central_beta.run_timed_tests([](const std::vector<double>& v, double x) { return pnbeta(x, v[0], v[1], v[2], 1, 0); }, "CDF", "Rmath " R_VERSION_STRING);
628 non_central_beta.run_timed_tests([](const std::vector<double>& v, double x) { return qnbeta(x, v[0], v[1], v[2], 1, 0); }, "quantile", "Rmath " R_VERSION_STRING, true);
629 // binomial
630 binomial.run_timed_tests([](const std::vector<double>& v, double x) { return dbinom(x, v[0], v[1], 0); }, "PDF", "Rmath " R_VERSION_STRING);
631 binomial.run_timed_tests([](const std::vector<double>& v, double x) { return pbinom(x, v[0], v[1], 1, 0); }, "CDF", "Rmath " R_VERSION_STRING);
632 binomial.run_timed_tests([](const std::vector<double>& v, double x) { return qbinom(x, v[0], v[1], 1, 0); }, "quantile", "Rmath " R_VERSION_STRING, true);
633 // cauchy
634 cauchy.run_timed_tests([](const std::vector<double>& v, double x) { return dcauchy(x, v[0], v[1], 0); }, "PDF", "Rmath " R_VERSION_STRING);
635 cauchy.run_timed_tests([](const std::vector<double>& v, double x) { return pcauchy(x, v[0], v[1], 1, 0); }, "CDF", "Rmath " R_VERSION_STRING);
636 cauchy.run_timed_tests([](const std::vector<double>& v, double x) { return qcauchy(x, v[0], v[1], 1, 0); }, "quantile", "Rmath " R_VERSION_STRING, true);
637 // chi squared
638 chi_squared.run_timed_tests([](const std::vector<double>& v, double x) { return dchisq(x, v[0], 0); }, "PDF", "Rmath " R_VERSION_STRING);
639 chi_squared.run_timed_tests([](const std::vector<double>& v, double x) { return pchisq(x, v[0], 1, 0); }, "CDF", "Rmath " R_VERSION_STRING);
640 chi_squared.run_timed_tests([](const std::vector<double>& v, double x) { return qchisq(x, v[0], 1, 0); }, "quantile", "Rmath " R_VERSION_STRING, true);
641 // non central chi squared
642 non_central_chi_squared.run_timed_tests([](const std::vector<double>& v, double x) { return dnchisq(x, v[0], v[1], 0); }, "PDF", "Rmath " R_VERSION_STRING);
643 non_central_chi_squared.run_timed_tests([](const std::vector<double>& v, double x) { return pnchisq(x, v[0], v[1], 1, 0); }, "CDF", "Rmath " R_VERSION_STRING);
644 non_central_chi_squared.run_timed_tests([](const std::vector<double>& v, double x) { return qnchisq(x, v[0], v[1], 1, 0); }, "quantile", "Rmath " R_VERSION_STRING, true);
645 // exponential
646 exponential.run_timed_tests([](const std::vector<double>& v, double x) { return dexp(x, 1 / v[0], 0); }, "PDF", "Rmath " R_VERSION_STRING);
647 exponential.run_timed_tests([](const std::vector<double>& v, double x) { return pexp(x, 1 / v[0], 1, 0); }, "CDF", "Rmath " R_VERSION_STRING);
648 exponential.run_timed_tests([](const std::vector<double>& v, double x) { return qexp(x, 1 / v[0], 1, 0); }, "quantile", "Rmath " R_VERSION_STRING, true);
649 // F
650 fisher.run_timed_tests([](const std::vector<double>& v, double x) { return df(x, v[0], v[1], 0); }, "PDF", "Rmath " R_VERSION_STRING);
651 fisher.run_timed_tests([](const std::vector<double>& v, double x) { return pf(x, v[0], v[1], 1, 0); }, "CDF", "Rmath " R_VERSION_STRING);
652 fisher.run_timed_tests([](const std::vector<double>& v, double x) { return qf(x, v[0], v[1], 1, 0); }, "quantile", "Rmath " R_VERSION_STRING, true);
653 // non central F
654 non_central_f.run_timed_tests([](const std::vector<double>& v, double x) { return dnf(x, v[0], v[1], v[2], 0); }, "PDF", "Rmath " R_VERSION_STRING);
655 non_central_f.run_timed_tests([](const std::vector<double>& v, double x) { return pnf(x, v[0], v[1], v[2], 1, 0); }, "CDF", "Rmath " R_VERSION_STRING);
656 non_central_f.run_timed_tests([](const std::vector<double>& v, double x) { return qnf(x, v[0], v[1], v[2], 1, 0); }, "quantile", "Rmath " R_VERSION_STRING, true);
657 // gamma
658 gamma.run_timed_tests([](const std::vector<double>& v, double x) { return dgamma(x, v[0], v[1], 0); }, "PDF", "Rmath " R_VERSION_STRING);
659 gamma.run_timed_tests([](const std::vector<double>& v, double x) { return pgamma(x, v[0], v[1], 1, 0); }, "CDF", "Rmath " R_VERSION_STRING);
660 gamma.run_timed_tests([](const std::vector<double>& v, double x) { return qgamma(x, v[0], v[1], 1, 0); }, "quantile", "Rmath " R_VERSION_STRING, true);
661 // geometric
662 geometric.run_timed_tests([](const std::vector<double>& v, double x) { return dgeom(x, v[0], 0); }, "PDF", "Rmath " R_VERSION_STRING);
663 geometric.run_timed_tests([](const std::vector<double>& v, double x) { return pgeom(x, v[0], 1, 0); }, "CDF", "Rmath " R_VERSION_STRING);
664 geometric.run_timed_tests([](const std::vector<double>& v, double x) { return qgeom(x, v[0], 1, 0); }, "quantile", "Rmath " R_VERSION_STRING, true);
665 // hypergeometric
666 hypergeometric.run_timed_tests([](const std::vector<double>& v, double x) { return dhyper(x, v[0], v[2] - v[0], v[1], 0); }, "PDF", "Rmath " R_VERSION_STRING);
667 hypergeometric.run_timed_tests([](const std::vector<double>& v, double x) { return phyper(x, v[0], v[2] - v[0], v[1], 1, 0); }, "CDF", "Rmath " R_VERSION_STRING);
668 hypergeometric.run_timed_tests([](const std::vector<double>& v, double x) { return qhyper(x, v[0], v[2] - v[0], v[1], 1, 0); }, "quantile", "Rmath " R_VERSION_STRING, true);
669 // logistic
670 logistic.run_timed_tests([](const std::vector<double>& v, double x) { return dlogis(x, v[0], v[1], 0); }, "PDF", "Rmath " R_VERSION_STRING);
671 logistic.run_timed_tests([](const std::vector<double>& v, double x) { return plogis(x, v[0], v[1], 1, 0); }, "CDF", "Rmath " R_VERSION_STRING);
672 logistic.run_timed_tests([](const std::vector<double>& v, double x) { return qlogis(x, v[0], v[1], 1, 0); }, "quantile", "Rmath " R_VERSION_STRING, true);
673 // lognormal
674 lognormal.run_timed_tests([](const std::vector<double>& v, double x) { return dlnorm(x, v[0], v[1], 0); }, "PDF", "Rmath " R_VERSION_STRING);
675 lognormal.run_timed_tests([](const std::vector<double>& v, double x) { return plnorm(x, v[0], v[1], 1, 0); }, "CDF", "Rmath " R_VERSION_STRING);
676 lognormal.run_timed_tests([](const std::vector<double>& v, double x) { return qlnorm(x, v[0], v[1], 1, 0); }, "quantile", "Rmath " R_VERSION_STRING, true);
677 // negative_binomial
678 negative_binomial.run_timed_tests([](const std::vector<double>& v, double x) { return dnbinom(x, v[0], v[1], 0); }, "PDF", "Rmath " R_VERSION_STRING);
679 negative_binomial.run_timed_tests([](const std::vector<double>& v, double x) { return pnbinom(x, v[0], v[1], 1, 0); }, "CDF", "Rmath " R_VERSION_STRING);
680 negative_binomial.run_timed_tests([](const std::vector<double>& v, double x) { return qnbinom(x, v[0], v[1], 1, 0); }, "quantile", "Rmath " R_VERSION_STRING, true);
681 // normal
682 n.run_timed_tests([](const std::vector<double>& v, double x) { return dnorm(x, v[0], v[1], 0); }, "PDF", "Rmath " R_VERSION_STRING);
683 n.run_timed_tests([](const std::vector<double>& v, double x) { return pnorm(x, v[0], v[1], 1, 0); }, "CDF", "Rmath " R_VERSION_STRING);
684 n.run_timed_tests([](const std::vector<double>& v, double x) { return qnorm(x, v[0], v[1], 1, 0); }, "quantile", "Rmath " R_VERSION_STRING, true);
685 // poisson
686 poisson.run_timed_tests([](const std::vector<double>& v, double x) { return dpois(x, v[0], 0); }, "PDF", "Rmath " R_VERSION_STRING);
687 poisson.run_timed_tests([](const std::vector<double>& v, double x) { return ppois(x, v[0], 1, 0); }, "CDF", "Rmath " R_VERSION_STRING);
688 poisson.run_timed_tests([](const std::vector<double>& v, double x) { return qpois(x, v[0], 1, 0); }, "quantile", "Rmath " R_VERSION_STRING, true);
689 // T
690 students_t.run_timed_tests([](const std::vector<double>& v, double x) { return dt(x, v[0], 0); }, "PDF", "Rmath " R_VERSION_STRING);
691 students_t.run_timed_tests([](const std::vector<double>& v, double x) { return pt(x, v[0], 1, 0); }, "CDF", "Rmath " R_VERSION_STRING);
692 students_t.run_timed_tests([](const std::vector<double>& v, double x) { return qt(x, v[0], 1, 0); }, "quantile", "Rmath " R_VERSION_STRING, true);
693 // non central T
694 non_central_t.run_timed_tests([](const std::vector<double>& v, double x) { return dnt(x, v[0], v[1], 0); }, "PDF", "Rmath " R_VERSION_STRING);
695 non_central_t.run_timed_tests([](const std::vector<double>& v, double x) { return pnt(x, v[0], v[1], 1, 0); }, "CDF", "Rmath " R_VERSION_STRING);
696 non_central_t.run_timed_tests([](const std::vector<double>& v, double x) { return qnt(x, v[0], v[1], 1, 0); }, "quantile", "Rmath " R_VERSION_STRING, true);
697 // weibull
698 weibull.run_timed_tests([](const std::vector<double>& v, double x) { return dweibull(x, v[0], v[1], 0); }, "PDF", "Rmath " R_VERSION_STRING);
699 weibull.run_timed_tests([](const std::vector<double>& v, double x) { return pweibull(x, v[0], v[1], 1, 0); }, "CDF", "Rmath " R_VERSION_STRING);
700 weibull.run_timed_tests([](const std::vector<double>& v, double x) { return qweibull(x, v[0], v[1], 1, 0); }, "quantile", "Rmath " R_VERSION_STRING, true);
701
702 #endif
703
704 #ifdef TEST_DCDFLIB
705 n.run_timed_tests([](const std::vector<double>& v, double x) { return dcdflib_norm_cdf(x, v[0], v[1]); }, "CDF", "DCDFLIB");
706 n.run_timed_tests([](const std::vector<double>& v, double x) { return dcdflib_norm_quantile(x, v[0], v[1]); }, "quantile", "DCDFLIB", true);
707
708 beta.run_timed_tests([](const std::vector<double>& v, double x) { return dcdflib_beta_cdf(x, v[0], v[1]); }, "CDF", "DCDFLIB");
709 beta.run_timed_tests([](const std::vector<double>& v, double x) { return dcdflib_beta_quantile(x, v[0], v[1]); }, "quantile", "DCDFLIB", true);
710
711 binomial.run_timed_tests([](const std::vector<double>& v, double x) { return dcdflib_binomial_cdf(x, v[0], v[1]); }, "CDF", "DCDFLIB");
712 binomial.run_timed_tests([](const std::vector<double>& v, double x) { return dcdflib_binomial_quantile(x, v[0], v[1]); }, "quantile", "DCDFLIB", true);
713
714 chi_squared.run_timed_tests([](const std::vector<double>& v, double x) { return dcdflib_chi_cdf(x, v[0]); }, "CDF", "DCDFLIB");
715 chi_squared.run_timed_tests([](const std::vector<double>& v, double x) { return dcdflib_chi_quantile(x, v[0]); }, "quantile", "DCDFLIB", true);
716
717 non_central_chi_squared.run_timed_tests([](const std::vector<double>& v, double x) { return dcdflib_chi_n_cdf(x, v[0], v[1]); }, "CDF", "DCDFLIB");
718 non_central_chi_squared.run_timed_tests([](const std::vector<double>& v, double x) { return dcdflib_chi_n_quantile(x, v[0], v[1]); }, "quantile", "DCDFLIB", true);
719
720 fisher.run_timed_tests([](const std::vector<double>& v, double x) { return dcdflib_f_cdf(x, v[0], v[1]); }, "CDF", "DCDFLIB");
721 fisher.run_timed_tests([](const std::vector<double>& v, double x) { return dcdflib_f_quantile(x, v[0], v[1]); }, "quantile", "DCDFLIB", true);
722
723 non_central_f.run_timed_tests([](const std::vector<double>& v, double x) { return dcdflib_f_n_cdf(x, v[0], v[1], v[2]); }, "CDF", "DCDFLIB");
724 non_central_f.run_timed_tests([](const std::vector<double>& v, double x) { return dcdflib_f_n_quantile(x, v[0], v[1], v[2]); }, "quantile", "DCDFLIB", true);
725
726 gamma.run_timed_tests([](const std::vector<double>& v, double x) { return dcdflib_gamma_cdf(x, v[0], v[1]); }, "CDF", "DCDFLIB");
727 gamma.run_timed_tests([](const std::vector<double>& v, double x) { return dcdflib_gamma_quantile(x, v[0], v[1]); }, "quantile", "DCDFLIB", true);
728
729 negative_binomial.run_timed_tests([](const std::vector<double>& v, double x) { return dcdflib_nbin_cdf(x, v[0], v[1]); }, "CDF", "DCDFLIB");
730 negative_binomial.run_timed_tests([](const std::vector<double>& v, double x) { return dcdflib_nbin_quantile(x, v[0], v[1]); }, "quantile", "DCDFLIB", true);
731
732 poisson.run_timed_tests([](const std::vector<double>& v, double x) { return dcdflib_poisson_cdf(x, v[0]); }, "CDF", "DCDFLIB");
733 poisson.run_timed_tests([](const std::vector<double>& v, double x) { return dcdflib_poisson_quantile(x, v[0]); }, "quantile", "DCDFLIB", true);
734
735 students_t.run_timed_tests([](const std::vector<double>& v, double x) { return dcdflib_t_cdf(x, v[0]); }, "CDF", "DCDFLIB");
736 students_t.run_timed_tests([](const std::vector<double>& v, double x) { return dcdflib_t_quantile(x, v[0]); }, "quantile", "DCDFLIB", true);
737
738 //non_central_t.run_timed_tests([](const std::vector<double>& v, double x) { return dcdflib_t_n_cdf(x, v[0], v[1]); }, "CDF", "DCDFLIB");
739 //non_central_t.run_timed_tests([](const std::vector<double>& v, double x) { return dcdflib_t_n_quantile(x, v[0], v[1]); }, "quantile", "DCDFLIB", true);
740 #endif
741
742 }
743 catch(const std::exception& e)
744 {
745 std::cout << "Test run aborted due to thrown exception: " << e.what() << std::endl;
746 return 1;
747 }
748 return 0;
749 }
750
751