test.c

来自「math library from gnu」· C语言 代码 · 共 1,309 行 · 第 1/5 页

  return gsl_ran_logarithmic_pdf (k, 0.9);
}

double
test_logarithmic_cdf_P (unsigned int k)
{
  return gsl_cdf_logarithmic_P (k, 0.9);
}

double
test_logarithmic_cdf_Q (unsigned int k)
{
  return gsl_cdf_logarithmic_Q (k, 0.9);
}
#endif

void
test_discrete_cdf_P (double (*pdf)(unsigned int), 
                     double (*cdf_P)(unsigned int), 
                     struct range (*range)(void), 
                     const char * desc)
{
  double sum;
  double tol = TEST_TOL2;
  int i, min, max;

  struct range r = range();

  min = r.min;
  max = r.max;
  sum = 0.0;

  for (i = min; i <= max; i++)
    {
      double pi = pdf(i);
      double Pi = cdf_P(i);
      sum += pi;
      gsl_test_rel (Pi, sum, tol, desc, i);
    }
}

void
test_discrete_cdf_Q (double (*pdf)(unsigned int), 
                     double (*cdf_Q)(unsigned int), 
                     struct range (*range)(void), 
                     const char * desc)
{
  double sum;
  double tol = TEST_TOL2;
  int i, min, max;

  struct range r = range();

  min = r.min;
  max = r.max;
  sum = cdf_Q(max);

  for (i = max; i >= min; i--)
    {
      double pi = pdf(i);
      double Qi = cdf_Q(i);
      gsl_test_rel (Qi, sum, tol, desc, i);
      sum += pi;
    }
}

void
test_discrete_cdf_PQ (double (*cdf_P)(unsigned int), 
                      double (*cdf_Q)(unsigned int), 
                      struct range (*range)(void), 
                      const char * desc)
{
  double sum;
  double tol = GSL_DBL_EPSILON;
  int i, min, max;

  struct range r = range();

  min = r.min;
  max = r.max;

  for (i = min; i <= max; i++)
    {
      double Pi = cdf_P(i);
      double Qi = cdf_Q(i);
      sum = Pi + Qi;
      gsl_test_rel (sum, 1.0, tol, desc, i);
      {
        int s1 = (Pi<0 || Pi>1);     
        int s2 = (Qi<0 || Qi>1);     
        gsl_test(s1, "Pi in range [0,1] (%.18e)", Pi);
        gsl_test(s2, "Qi in range [0,1] (%.18e)", Qi);
      }
    }

}

#define TEST_DISCRETE(name) do {  \
  test_discrete_cdf_P(&test_ ## name ## _pdf, &test_ ## name ## _cdf_P, &test_ ## name ## _range, "test gsl_cdf_" #name "_P (k=%d)") ; \
  test_discrete_cdf_Q(&test_ ## name ## _pdf, &test_ ## name ## _cdf_Q, &test_ ## name ## _range, "test gsl_cdf_" #name "_Q (k=%d)") ; \
  test_discrete_cdf_PQ(&test_ ## name ## _cdf_P, &test_ ## name ## _cdf_Q, &test_ ## name ## _range, "test gsl_cdf_" #name "_P+Q (k=%d)") ; \
} while (0);

int
main (void)
{
  gsl_ieee_env_setup ();

  TEST_DISCRETE(binomial);
  TEST_DISCRETE(poisson);
  TEST_DISCRETE(geometric);
  TEST_DISCRETE(negative_binomial);
  TEST_DISCRETE(pascal);
  TEST_DISCRETE(hypergeometric);
#ifdef HYPERGEOMETRIC2
  TEST_DISCRETE(hypergeometric2);
#endif
#ifdef LOGARITHMIC
  TEST_DISCRETE(logarithmic);
#endif

  test_discrete_cdf_PQ(&test_hypergeometric2_cdf_P, 
                       &test_hypergeometric2_cdf_Q, 
                       &test_hypergeometric2a_range, 
                       "test gsl_cdf_hypergeometric_P+Q (k=%d)") ; 

  test_discrete_cdf_PQ(&test_hypergeometric2_cdf_P, 
                       &test_hypergeometric2_cdf_Q, 
                       &test_hypergeometric2b_range, 
                       "test gsl_cdf_hypergeometric_P+Q (k=%d)") ; 


  /* exit (gsl_test_summary ()); */

  /* Tests for gaussian cumulative distribution function 
     Function values computed with PARI, 28 digits precision */
  
  test_ugaussian ();
  test_exponential ();
  test_exppow ();
  test_tdist (); 
  test_fdist (); 
  test_gamma ();
  test_chisq (); 
  test_beta (); 

  test_ugaussianinv ();
  test_exponentialinv ();
  test_gammainv (); 
  test_chisqinv (); 
  test_tdistinv (); 
  test_betainv ();
  test_finv ();

  test_auto_beta ();
  test_auto_fdist ();
  test_auto_cauchy ();
  test_auto_gaussian ();
  test_auto_laplace ();
  test_auto_rayleigh ();
  test_auto_flat ();
  test_auto_lognormal ();
  test_auto_gamma ();
  test_auto_chisq ();
  test_auto_tdist ();
  test_auto_gumbel1 ();
  test_auto_gumbel2 ();
  test_auto_weibull ();
  test_auto_pareto ();
  test_auto_logistic ();
  test_auto_gammalarge ();
  
  exit (gsl_test_summary ());
}

void test_ugaussian (void)
{
  TEST (gsl_cdf_ugaussian_P, (0.0), 0.5, TEST_TOL0);
  TEST (gsl_cdf_ugaussian_P, (1e-32), 0.5, TEST_TOL0);
  TEST (gsl_cdf_ugaussian_P, (1e-16), 0.5000000000000000398942280401, TEST_TOL0);
  TEST (gsl_cdf_ugaussian_P, (1e-8), 0.5000000039894228040143267129, TEST_TOL0);
  TEST (gsl_cdf_ugaussian_P, (0.5), 0.6914624612740131036377046105, TEST_TOL0);
  TEST (gsl_cdf_ugaussian_P, (0.7), 0.7580363477769269852506495717, TEST_TOL0);
  TEST (gsl_cdf_ugaussian_P, (5.0), 0.9999997133484281208060883262, TEST_TOL0);
  TEST (gsl_cdf_ugaussian_P, (10.0), 0.9999999999999999999999923801, TEST_TOL0);
  TEST (gsl_cdf_ugaussian_P, (30.0), 1.000000000000000000000000000, TEST_TOL0);
  TEST (gsl_cdf_ugaussian_P, (40.0), 1.000000000000000000000000000, TEST_TOL0);
  TEST (gsl_cdf_ugaussian_P, (1e10), 1.000000000000000000000000000, TEST_TOL0);

  TEST (gsl_cdf_ugaussian_P, (-1e-32), 0.5, TEST_TOL0);
  TEST (gsl_cdf_ugaussian_P, (-1e-16), 0.4999999999999999601057719598, TEST_TOL0);
  TEST (gsl_cdf_ugaussian_P, (-1e-8), 0.4999999960105771959856732870, TEST_TOL0);
  TEST (gsl_cdf_ugaussian_P, (-0.5), 0.3085375387259868963622953894, TEST_TOL0);
  TEST (gsl_cdf_ugaussian_P, (-0.7), 0.2419636522230730147493504282, TEST_TOL0);
  TEST (gsl_cdf_ugaussian_P, (-5.0), 0.0000002866515718791939116737523329, TEST_TOL1);
  TEST (gsl_cdf_ugaussian_P, (-10.0), 7.619853024160526065973343257e-24, TEST_TOL3);
  TEST (gsl_cdf_ugaussian_P, (-30.0), 4.906713927148187059533809288e-198, TEST_TOL3);
  TEST (gsl_cdf_ugaussian_P, (-1e10), 0.0, 0.0);

  TEST (gsl_cdf_ugaussian_Q, (0.0), 0.5, TEST_TOL0);
  TEST (gsl_cdf_ugaussian_Q, (1e-32), 0.5, TEST_TOL0);
  TEST (gsl_cdf_ugaussian_Q, (1e-16), 0.4999999999999999601057719598, TEST_TOL0);
  TEST (gsl_cdf_ugaussian_Q, (1e-8), 0.4999999960105771959856732870, TEST_TOL0);
  TEST (gsl_cdf_ugaussian_Q, (0.5), 0.3085375387259868963622953894, TEST_TOL0);
  TEST (gsl_cdf_ugaussian_Q, (0.7), 0.2419636522230730147493504282, TEST_TOL0);
  TEST (gsl_cdf_ugaussian_Q, (5.0), 0.0000002866515718791939116737523329, TEST_TOL3);
  TEST (gsl_cdf_ugaussian_Q, (10.0), 7.619853024160526065973343257e-24, TEST_TOL3);
  TEST (gsl_cdf_ugaussian_Q, (30.0), 4.906713927148187059533809288e-198, TEST_TOL3);
  TEST (gsl_cdf_ugaussian_Q, (1e10), 0.0, 0.0);

  TEST (gsl_cdf_ugaussian_Q, (-1e-32), 0.5, TEST_TOL0);
  TEST (gsl_cdf_ugaussian_Q, (-1e-16), 0.5000000000000000398942280401, TEST_TOL0);
  TEST (gsl_cdf_ugaussian_Q, (-1e-8), 0.5000000039894228040143267129, TEST_TOL0);
  TEST (gsl_cdf_ugaussian_Q, (-0.5), 0.6914624612740131036377046105, TEST_TOL0);
  TEST (gsl_cdf_ugaussian_Q, (-0.7), 0.7580363477769269852506495717, TEST_TOL0);
  TEST (gsl_cdf_ugaussian_Q, (-5.0), 0.9999997133484281208060883262, TEST_TOL0);
  TEST (gsl_cdf_ugaussian_Q, (-10.0), 0.9999999999999999999999923801, TEST_TOL0);
  TEST (gsl_cdf_ugaussian_Q, (-30.0), 1.000000000000000000000000000, TEST_TOL0);
  TEST (gsl_cdf_ugaussian_Q, (-40.0), 1.000000000000000000000000000, TEST_TOL0);
  TEST (gsl_cdf_ugaussian_Q, (-1e10), 1.000000000000000000000000000, TEST_TOL0);
}
  /* Test values from Abramowitz & Stegun, Handbook of Mathematical
     Functions, Table 26.1.  Error term is given by dx = dP / Z(x) */

void test_ugaussianinv (void) {
  TEST (gsl_cdf_ugaussian_Pinv, (0.9999997133), 5.0, 1e-4);
  TEST (gsl_cdf_ugaussian_Pinv, (0.9999683288), 4.0, 1e-6);
  TEST (gsl_cdf_ugaussian_Pinv, (0.9986501020), 3.0, 1e-8);
  TEST (gsl_cdf_ugaussian_Pinv, (0.977249868051821), 2.0, 1e-14);
  TEST (gsl_cdf_ugaussian_Pinv, (0.841344746068543), 1.0, TEST_TOL3);
  TEST (gsl_cdf_ugaussian_Pinv, (0.691462461274013), 0.5, TEST_TOL2);
  TEST (gsl_cdf_ugaussian_Pinv, (0.655421741610324), 0.4, TEST_TOL1);
  TEST (gsl_cdf_ugaussian_Pinv, (0.617911422188953), 0.3, TEST_TOL1);
  TEST (gsl_cdf_ugaussian_Pinv, (0.579259709439103), 0.2, TEST_TOL1);
  TEST (gsl_cdf_ugaussian_Pinv, (0.539827837277029), 0.1, TEST_TOL1);
  TEST (gsl_cdf_ugaussian_Pinv, (0.5), 0.0, TEST_TOL0);
  TEST (gsl_cdf_ugaussian_Pinv, (4.60172162722971e-1), -0.1, TEST_TOL1);
  TEST (gsl_cdf_ugaussian_Pinv, (4.20740290560897e-1), -0.2, TEST_TOL1);
  TEST (gsl_cdf_ugaussian_Pinv, (3.82088577811047e-1), -0.3, TEST_TOL1);
  TEST (gsl_cdf_ugaussian_Pinv, (3.44578258389676e-1), -0.4, TEST_TOL1);
  TEST (gsl_cdf_ugaussian_Pinv, (3.08537538725987e-1), -0.5, TEST_TOL2);
  TEST (gsl_cdf_ugaussian_Pinv, (1.58655253931457e-1), -1.0, TEST_TOL3);
  TEST (gsl_cdf_ugaussian_Pinv, (2.2750131948179e-2), -2.0, 1e-14);
  TEST (gsl_cdf_ugaussian_Pinv, (1.349898e-3), -3.0, 1e-8);
  TEST (gsl_cdf_ugaussian_Pinv, (3.16712e-5), -4.0, 1e-6);
  TEST (gsl_cdf_ugaussian_Pinv, (2.86648e-7), -5.0, 1e-4);

  TEST (gsl_cdf_ugaussian_Pinv, (7.61985302416052e-24), -10.0, 1e-4);

  TEST (gsl_cdf_ugaussian_Qinv, (7.61985302416052e-24), 10.0, 1e-4);

  TEST (gsl_cdf_ugaussian_Qinv, (2.86648e-7), 5.0, 1e-4);
  TEST (gsl_cdf_ugaussian_Qinv, (3.16712e-5), 4.0, 1e-6);
  TEST (gsl_cdf_ugaussian_Qinv, (1.349898e-3), 3.0, 1e-8);
  TEST (gsl_cdf_ugaussian_Qinv, (2.2750131948179e-2), 2.0, 1e-14);
  TEST (gsl_cdf_ugaussian_Qinv, (1.58655253931457e-1), 1.0, TEST_TOL3);
  TEST (gsl_cdf_ugaussian_Qinv, (3.08537538725987e-1), 0.5, TEST_TOL2);
  TEST (gsl_cdf_ugaussian_Qinv, (3.44578258389676e-1), 0.4, TEST_TOL1);
  TEST (gsl_cdf_ugaussian_Qinv, (3.82088577811047e-1), 0.3, TEST_TOL1);
  TEST (gsl_cdf_ugaussian_Qinv, (4.20740290560897e-1), 0.2, TEST_TOL1);
  TEST (gsl_cdf_ugaussian_Qinv, (4.60172162722971e-1), 0.1, TEST_TOL1);
  TEST (gsl_cdf_ugaussian_Qinv, (0.5), 0.0, TEST_TOL0);