#include "performance_measure.hpp"

#include <boost/math/distributions/beta.hpp>
#include <boost/math/distributions/binomial.hpp>
#include <boost/math/distributions/cauchy.hpp>
#include <boost/math/distributions/chi_squared.hpp>
#include <boost/math/distributions/exponential.hpp>
#include <boost/math/distributions/fisher_f.hpp>
#include <boost/math/distributions/gamma.hpp>
#include <boost/math/distributions/lognormal.hpp>
#include <boost/math/distributions/negative_binomial.hpp>
#include <boost/math/distributions/normal.hpp>
#include <boost/math/distributions/poisson.hpp>
#include <boost/math/distributions/students_t.hpp>
#include <boost/math/distributions/weibull.hpp>

double probabilities[] = {
   1e-5,
   1e-4,
   1e-3,
   1e-2,
   0.05,
   0.1,
   0.2,
   0.3,
   0.4,
   0.5,
   0.6,
   0.7,
   0.8,
   0.9,
   0.95,
   1-1e-5,
   1-1e-4,
   1-1e-3,
   1-1e-2
};

int int_values[] = {
   1,
   2,
   3,
   5,
   10,
   20,
   50,
   100,
   1000,
   10000,
   100000
};

double real_values[] = {
   1e-5,
   1e-4,
   1e-2,
   1e-1,
   1,
   10,
   100,
   1000,
   10000,
   100000
};

#define BOOST_MATH_DISTRIBUTION2_TEST(name, param1_table, param2_table, random_variable_table, probability_table) \
   BOOST_MATH_PERFORMANCE_TEST(BOOST_JOIN(dist, name), "dist-" BOOST_STRINGIZE(name) "-cdf")\
   {\
   double result = 0;\
   unsigned a_size = sizeof(param1_table)/sizeof(param1_table[0]);\
   unsigned b_size = sizeof(param2_table)/sizeof(param2_table[0]);\
   unsigned c_size = sizeof(random_variable_table)/sizeof(random_variable_table[0]);\
   \
   for(unsigned i = 0; i < a_size; ++i)\
   {\
      for(unsigned j = 0; j < b_size; ++j)\
      {\
         for(unsigned k = 0; k < c_size; ++k)\
         {\
            result += cdf(boost::math:: BOOST_JOIN(name, _distribution) <>(param1_table[i], param2_table[j]), random_variable_table[k]);\
         }\
      }\
   }\
   \
   consume_result(result);\
   set_call_count(a_size * b_size * c_size);\
   }\
   BOOST_MATH_PERFORMANCE_TEST(BOOST_JOIN(dist_pdf_, name), "dist-" BOOST_STRINGIZE(name) "-pdf")\
   {\
   double result = 0;\
   unsigned a_size = sizeof(param1_table)/sizeof(param1_table[0]);\
   unsigned b_size = sizeof(param2_table)/sizeof(param2_table[0]);\
   unsigned c_size = sizeof(random_variable_table)/sizeof(random_variable_table[0]);\
   \
   for(unsigned i = 0; i < a_size; ++i)\
   {\
      for(unsigned j = 0; j < b_size; ++j)\
      {\
         for(unsigned k = 0; k < c_size; ++k)\
         {\
            result += pdf(boost::math:: BOOST_JOIN(name, _distribution) <>(param1_table[i], param2_table[j]), random_variable_table[k]);\
         }\
      }\
   }\
   \
   consume_result(result);\
   set_call_count(a_size * b_size * c_size);\
   }\
   BOOST_MATH_PERFORMANCE_TEST(BOOST_JOIN(dist_quant, name), "dist-" BOOST_STRINGIZE(name) "-quantile")\
   {\
      double result = 0;\
      unsigned a_size = sizeof(param1_table)/sizeof(param1_table[0]);\
      unsigned b_size = sizeof(param2_table)/sizeof(param2_table[0]);\
      unsigned c_size = sizeof(probability_table)/sizeof(probability_table[0]);\
      \
      for(unsigned i = 0; i < a_size; ++i)\
      {\
         for(unsigned j = 0; j < b_size; ++j)\
         {\
            for(unsigned k = 0; k < c_size; ++k)\
            {\
               result += quantile(boost::math:: BOOST_JOIN(name, _distribution) <>(param1_table[i], param2_table[j]), probability_table[k]);\
            }\
         }\
      }\
      \
      consume_result(result);\
      set_call_count(a_size * b_size * c_size);\
   }

#define BOOST_MATH_DISTRIBUTION1_TEST(name, param1_table, random_variable_table, probability_table) \
   BOOST_MATH_PERFORMANCE_TEST(BOOST_JOIN(dist, name), "dist-" BOOST_STRINGIZE(name) "-cdf")\
   {\
   double result = 0;\
   unsigned a_size = sizeof(param1_table)/sizeof(param1_table[0]);\
   unsigned c_size = sizeof(random_variable_table)/sizeof(random_variable_table[0]);\
   \
   for(unsigned i = 0; i < a_size; ++i)\
   {\
         for(unsigned k = 0; k < c_size; ++k)\
         {\
            result += cdf(boost::math:: BOOST_JOIN(name, _distribution) <>(param1_table[i]), random_variable_table[k]);\
         }\
   }\
   \
   consume_result(result);\
   set_call_count(a_size * c_size);\
   }\
   BOOST_MATH_PERFORMANCE_TEST(BOOST_JOIN(dist_pdf_, name), "dist-" BOOST_STRINGIZE(name) "-pdf")\
   {\
   double result = 0;\
   unsigned a_size = sizeof(param1_table)/sizeof(param1_table[0]);\
   unsigned c_size = sizeof(random_variable_table)/sizeof(random_variable_table[0]);\
   \
   for(unsigned i = 0; i < a_size; ++i)\
   {\
         for(unsigned k = 0; k < c_size; ++k)\
         {\
            result += pdf(boost::math:: BOOST_JOIN(name, _distribution) <>(param1_table[i]), random_variable_table[k]);\
         }\
   }\
   \
   consume_result(result);\
   set_call_count(a_size * c_size);\
   }\
   BOOST_MATH_PERFORMANCE_TEST(BOOST_JOIN(dist_quant, name), "dist-" BOOST_STRINGIZE(name) "-quantile")\
   {\
      double result = 0;\
      unsigned a_size = sizeof(param1_table)/sizeof(param1_table[0]);\
      unsigned c_size = sizeof(probability_table)/sizeof(probability_table[0]);\
      \
      for(unsigned i = 0; i < a_size; ++i)\
      {\
            for(unsigned k = 0; k < c_size; ++k)\
            {\
               result += quantile(boost::math:: BOOST_JOIN(name, _distribution) <>(param1_table[i]), probability_table[k]);\
            }\
      }\
      \
      consume_result(result);\
      set_call_count(a_size * c_size);\
   }

BOOST_MATH_DISTRIBUTION2_TEST(beta, probabilities, probabilities, probabilities, probabilities)
BOOST_MATH_DISTRIBUTION2_TEST(binomial, int_values, probabilities, int_values, probabilities)
BOOST_MATH_DISTRIBUTION2_TEST(cauchy, int_values, real_values, int_values, probabilities)
BOOST_MATH_DISTRIBUTION1_TEST(chi_squared, int_values, real_values, probabilities)
BOOST_MATH_DISTRIBUTION1_TEST(exponential, real_values, real_values, probabilities)
BOOST_MATH_DISTRIBUTION2_TEST(fisher_f, int_values, int_values, real_values, probabilities)
BOOST_MATH_DISTRIBUTION2_TEST(gamma, real_values, real_values, real_values, probabilities)
BOOST_MATH_DISTRIBUTION2_TEST(lognormal, real_values, real_values, real_values, probabilities)
BOOST_MATH_DISTRIBUTION2_TEST(negative_binomial, int_values, probabilities, int_values, probabilities)
BOOST_MATH_DISTRIBUTION2_TEST(normal, real_values, real_values, real_values, probabilities)
BOOST_MATH_DISTRIBUTION1_TEST(poisson, real_values, real_values, probabilities)
BOOST_MATH_DISTRIBUTION1_TEST(students_t, int_values, real_values, probabilities)
BOOST_MATH_DISTRIBUTION2_TEST(weibull, real_values, real_values, real_values, probabilities)

#ifdef TEST_R

extern "C" {
#include "Rmath.h"
}

#define BOOST_MATH_R_DISTRIBUTION2_TEST(name, param1_table, param2_table, random_variable_table, probability_table) \
   BOOST_MATH_PERFORMANCE_TEST(BOOST_JOIN(dist, name), "dist-" #name "-R-cdf")\
   {\
   double result = 0;\
   unsigned a_size = sizeof(param1_table)/sizeof(param1_table[0]);\
   unsigned b_size = sizeof(param2_table)/sizeof(param2_table[0]);\
   unsigned c_size = sizeof(random_variable_table)/sizeof(random_variable_table[0]);\
   \
   for(unsigned i = 0; i < a_size; ++i)\
   {\
      for(unsigned j = 0; j < b_size; ++j)\
      {\
         for(unsigned k = 0; k < c_size; ++k)\
         {\
         result += Rf_p##name (random_variable_table[k], param1_table[i], param2_table[j], 1, 0);\
         }\
      }\
   }\
   \
   consume_result(result);\
   set_call_count(a_size * b_size * c_size);\
   }\
   BOOST_MATH_PERFORMANCE_TEST(BOOST_JOIN(dist_pdf_, name), "dist-" #name "-R-pdf")\
   {\
   double result = 0;\
   unsigned a_size = sizeof(param1_table)/sizeof(param1_table[0]);\
   unsigned b_size = sizeof(param2_table)/sizeof(param2_table[0]);\
   unsigned c_size = sizeof(random_variable_table)/sizeof(random_variable_table[0]);\
   \
   for(unsigned i = 0; i < a_size; ++i)\
   {\
      for(unsigned j = 0; j < b_size; ++j)\
      {\
         for(unsigned k = 0; k < c_size; ++k)\
         {\
            result += Rf_d##name (random_variable_table[k], param1_table[i], param2_table[j], 0);\
         }\
      }\
   }\
   \
   consume_result(result);\
   set_call_count(a_size * b_size * c_size);\
   }\
   BOOST_MATH_PERFORMANCE_TEST(BOOST_JOIN(dist_quant, name), "dist-" #name "-R-quantile")\
   {\
      double result = 0;\
      unsigned a_size = sizeof(param1_table)/sizeof(param1_table[0]);\
      unsigned b_size = sizeof(param2_table)/sizeof(param2_table[0]);\
      unsigned c_size = sizeof(probability_table)/sizeof(probability_table[0]);\
      \
      for(unsigned i = 0; i < a_size; ++i)\
      {\
         for(unsigned j = 0; j < b_size; ++j)\
         {\
            for(unsigned k = 0; k < c_size; ++k)\
            {\
               result += Rf_q##name (probability_table[k], param1_table[i], param2_table[j], 1, 0);\
            }\
         }\
      }\
      \
      consume_result(result);\
      set_call_count(a_size * b_size * c_size);\
   }

#define BOOST_MATH_R_DISTRIBUTION1_TEST(name, param1_table, random_variable_table, probability_table) \
   BOOST_MATH_PERFORMANCE_TEST(BOOST_JOIN(dist, name), "dist-" #name "-R-cdf")\
   {\
   double result = 0;\
   unsigned a_size = sizeof(param1_table)/sizeof(param1_table[0]);\
   unsigned c_size = sizeof(random_variable_table)/sizeof(random_variable_table[0]);\
   \
   for(unsigned i = 0; i < a_size; ++i)\
   {\
         for(unsigned k = 0; k < c_size; ++k)\
         {\
         result += Rf_p##name (random_variable_table[k], param1_table[i], 1, 0);\
         }\
   }\
   \
   consume_result(result);\
   set_call_count(a_size * c_size);\
   }\
   BOOST_MATH_PERFORMANCE_TEST(BOOST_JOIN(dist_pdf_, name), "dist-" #name "-R-pdf")\
   {\
   double result = 0;\
   unsigned a_size = sizeof(param1_table)/sizeof(param1_table[0]);\
   unsigned c_size = sizeof(random_variable_table)/sizeof(random_variable_table[0]);\
   \
   for(unsigned i = 0; i < a_size; ++i)\
   {\
         for(unsigned k = 0; k < c_size; ++k)\
         {\
            result += Rf_d##name (random_variable_table[k], param1_table[i], 0);\
         }\
   }\
   \
   consume_result(result);\
   set_call_count(a_size * c_size);\
   }\
   BOOST_MATH_PERFORMANCE_TEST(BOOST_JOIN(dist_quant, name), "dist-" #name "-R-quantile")\
   {\
      double result = 0;\
      unsigned a_size = sizeof(param1_table)/sizeof(param1_table[0]);\
      unsigned c_size = sizeof(probability_table)/sizeof(probability_table[0]);\
      \
      for(unsigned i = 0; i < a_size; ++i)\
      {\
            for(unsigned k = 0; k < c_size; ++k)\
            {\
               result += Rf_q##name (probability_table[k], param1_table[i], 1, 0);\
            }\
      }\
      \
      consume_result(result);\
      set_call_count(a_size * c_size);\
   }
//
// The R lib actually gives these non-obvious names that we can't deduce
// in our own macro code, so create unline forwarders and let R's own
// macros have control over dispatch to the right algorithm:
//
inline double Rf_dnorm(double x, double mu, double sigma, int give_log)
{
   return dnorm(x, mu, sigma, give_log);
}

inline double Rf_pnorm(double x, double mu, double sigma, int lower_tail, int give_log)
{
   return pnorm(x, mu, sigma, lower_tail, give_log);
}

inline double Rf_qnorm(double p, double mu, double sigma, int lower_tail, int log_p)
{
   return qnorm(p, mu, sigma, lower_tail, log_p);
}

BOOST_MATH_R_DISTRIBUTION2_TEST(beta, probabilities, probabilities, probabilities, probabilities)
BOOST_MATH_R_DISTRIBUTION2_TEST(binom, int_values, probabilities, int_values, probabilities)
BOOST_MATH_R_DISTRIBUTION2_TEST(cauchy, int_values, real_values, int_values, probabilities)
BOOST_MATH_R_DISTRIBUTION1_TEST(chisq, int_values, real_values, probabilities)
BOOST_MATH_R_DISTRIBUTION1_TEST(exp, real_values, real_values, probabilities)
BOOST_MATH_R_DISTRIBUTION2_TEST(f, int_values, int_values, real_values, probabilities)
BOOST_MATH_R_DISTRIBUTION2_TEST(gamma, real_values, real_values, real_values, probabilities)
BOOST_MATH_R_DISTRIBUTION2_TEST(lnorm, real_values, real_values, real_values, probabilities)
BOOST_MATH_R_DISTRIBUTION2_TEST(nbinom, int_values, probabilities, int_values, probabilities)
BOOST_MATH_R_DISTRIBUTION2_TEST(norm, real_values, real_values, real_values, probabilities)
BOOST_MATH_R_DISTRIBUTION1_TEST(pois, real_values, real_values, probabilities)
BOOST_MATH_R_DISTRIBUTION1_TEST(t, int_values, real_values, probabilities)
BOOST_MATH_R_DISTRIBUTION2_TEST(weibull, real_values, real_values, real_values, probabilities)

#endif