123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195 |
- /* test_real_distribution.ipp
- *
- * Copyright Steven Watanabe 2011
- * Distributed under the Boost Software License, Version 1.0. (See
- * accompanying file LICENSE_1_0.txt or copy at
- * http://www.boost.org/LICENSE_1_0.txt)
- *
- * $Id$
- *
- */
- #ifndef BOOST_MATH_DISTRIBUTION_INIT
- #ifdef BOOST_RANDOM_ARG2_TYPE
- #define BOOST_MATH_DISTRIBUTION_INIT (BOOST_RANDOM_ARG1_NAME, BOOST_RANDOM_ARG2_NAME)
- #else
- #define BOOST_MATH_DISTRIBUTION_INIT (BOOST_RANDOM_ARG1_NAME)
- #endif
- #endif
- #ifndef BOOST_RANDOM_DISTRIBUTION_INIT
- #ifdef BOOST_RANDOM_ARG2_TYPE
- #define BOOST_RANDOM_DISTRIBUTION_INIT (BOOST_RANDOM_ARG1_NAME, BOOST_RANDOM_ARG2_NAME)
- #else
- #define BOOST_RANDOM_DISTRIBUTION_INIT (BOOST_RANDOM_ARG1_NAME)
- #endif
- #endif
- #ifndef BOOST_RANDOM_P_CUTOFF
- #define BOOST_RANDOM_P_CUTOFF 0.99
- #endif
- #include <boost/random/mersenne_twister.hpp>
- #include <boost/lexical_cast.hpp>
- #include <boost/exception/diagnostic_information.hpp>
- #include <boost/preprocessor/stringize.hpp>
- #include <boost/range/numeric.hpp>
- #include <boost/numeric/conversion/cast.hpp>
- #include <iostream>
- #include <vector>
- #include "statistic_tests.hpp"
- #include "chi_squared_test.hpp"
- bool do_test(BOOST_RANDOM_ARG1_TYPE BOOST_RANDOM_ARG1_NAME,
- #ifdef BOOST_RANDOM_ARG2_TYPE
- BOOST_RANDOM_ARG2_TYPE BOOST_RANDOM_ARG2_NAME,
- #endif
- long long max, boost::mt19937& gen) {
- std::cout << "running " BOOST_PP_STRINGIZE(BOOST_RANDOM_DISTRIBUTION_NAME) "("
- << BOOST_RANDOM_ARG1_NAME;
- #ifdef BOOST_RANDOM_ARG2_NAME
- std::cout << ", " << BOOST_RANDOM_ARG2_NAME;
- #endif
- std::cout << ")" << " " << max << " times: " << std::flush;
- BOOST_MATH_DISTRIBUTION expected BOOST_MATH_DISTRIBUTION_INIT;
-
- BOOST_RANDOM_DISTRIBUTION dist BOOST_RANDOM_DISTRIBUTION_INIT;
- #ifdef BOOST_RANDOM_DISTRIBUTION_MAX
- BOOST_RANDOM_DISTRIBUTION::result_type max_value = BOOST_RANDOM_DISTRIBUTION_MAX;
- std::vector<double> expected_pdf(max_value+1);
- {
- for(int i = 0; i <= max_value; ++i) {
- expected_pdf[i] = pdf(expected, i);
- }
- expected_pdf.back() += 1 - boost::accumulate(expected_pdf, 0.0);
- }
-
- std::vector<long long> results(max_value + 1);
- for(long long i = 0; i < max; ++i) {
- ++results[(std::min)(dist(gen), max_value)];
- }
- long long sum = boost::accumulate(results, 0ll);
- if(sum != max) {
- std::cout << "*** Failed: incorrect total: " << sum << " ***" << std::endl;
- return false;
- }
- double prob = chi_squared_test(results, expected_pdf, max);
- #else
- kolmogorov_experiment test(boost::numeric_cast<int>(max));
- boost::variate_generator<boost::mt19937&, BOOST_RANDOM_DISTRIBUTION > vgen(gen, dist);
- double prob = test.probability(test.run(vgen, expected));
- #endif
- bool result = prob < BOOST_RANDOM_P_CUTOFF;
- const char* err = result? "" : "*";
- std::cout << std::setprecision(17) << prob << err << std::endl;
- std::cout << std::setprecision(6);
- return result;
- }
- template<class Dist1
- #ifdef BOOST_RANDOM_ARG2_NAME
- , class Dist2
- #endif
- >
- bool do_tests(int repeat, Dist1 d1,
- #ifdef BOOST_RANDOM_ARG2_NAME
- Dist2 d2,
- #endif
- long long trials) {
- boost::mt19937 gen;
- int errors = 0;
- for(int i = 0; i < repeat; ++i) {
- if(!do_test(d1(gen),
- #ifdef BOOST_RANDOM_ARG2_NAME
- d2(gen),
- #endif
- trials, gen)) {
- ++errors;
- }
- }
- if(errors != 0) {
- std::cout << "*** " << errors << " errors detected ***" << std::endl;
- }
- return errors == 0;
- }
- int usage() {
- std::cerr << "Usage: test_" BOOST_PP_STRINGIZE(BOOST_RANDOM_DISTRIBUTION_NAME)
- " -r <repeat>"
- " -" BOOST_PP_STRINGIZE(BOOST_RANDOM_ARG1_NAME)
- " <max " BOOST_PP_STRINGIZE(BOOST_RANDOM_ARG1_NAME) ">"
- #ifdef BOOST_RANDOM_ARG2_NAME
- " -" BOOST_PP_STRINGIZE(BOOST_RANDOM_ARG2_NAME)
- " <max " BOOST_PP_STRINGIZE(BOOST_RANDOM_ARG2_NAME) ">"
- #endif
- " -t <trials>" << std::endl;
- return 2;
- }
- template<class T>
- bool handle_option(int& argc, char**& argv, const char* opt, T& value) {
- if(std::strcmp(argv[0], opt) == 0 && argc > 1) {
- --argc;
- ++argv;
- value = boost::lexical_cast<T>(argv[0]);
- return true;
- } else {
- return false;
- }
- }
- int main(int argc, char** argv) {
- int repeat = 1;
- BOOST_RANDOM_ARG1_TYPE max_arg1 = BOOST_RANDOM_ARG1_DEFAULT;
- #ifdef BOOST_RANDOM_ARG2_TYPE
- BOOST_RANDOM_ARG2_TYPE max_arg2 = BOOST_RANDOM_ARG2_DEFAULT;
- #endif
- long long trials = 100000;
- if(argc > 0) {
- --argc;
- ++argv;
- }
- while(argc > 0) {
- if(argv[0][0] != '-') return usage();
- else if(!handle_option(argc, argv, "-r", repeat)
- && !handle_option(argc, argv, "-" BOOST_PP_STRINGIZE(BOOST_RANDOM_ARG1_NAME), max_arg1)
- #ifdef BOOST_RANDOM_ARG2_TYPE
- && !handle_option(argc, argv, "-" BOOST_PP_STRINGIZE(BOOST_RANDOM_ARG2_NAME), max_arg2)
- #endif
- && !handle_option(argc, argv, "-t", trials)) {
- return usage();
- }
- --argc;
- ++argv;
- }
- try {
- if(do_tests(repeat,
- BOOST_RANDOM_ARG1_DISTRIBUTION(max_arg1),
- #ifdef BOOST_RANDOM_ARG2_TYPE
- BOOST_RANDOM_ARG2_DISTRIBUTION(max_arg2),
- #endif
- trials)) {
- return 0;
- } else {
- return EXIT_FAILURE;
- }
- } catch(...) {
- std::cerr << boost::current_exception_diagnostic_information() << std::endl;
- return EXIT_FAILURE;
- }
- }
|