This is the mail archive of the gsl-discuss@sources.redhat.com mailing list for the GSL project.


Index Nav: [Date Index] [Subject Index] [Author Index] [Thread Index]
Message Nav: [Date Prev] [Date Next] [Thread Prev] [Thread Next]
Other format: [Raw text]

MISER integration problem


Hello,

I am trying to use the gsl Monte Carlo integration routines...
I hope this is the right place to ask my question!


Unfortunately, depending on the value of certain parameter I get the
following error when using the MISER routines:

gsl: miser.c:115: ERROR: insufficient calls for subvolume
Aborted


As you can see below, I effectively copied the sample program from the
documentation. The error arrises when I change the 3rd line in main() 

from :  double xu[2] = { 1000, 1000};
to   :  double xu[2] = { 100, 100};

I hope someone can help me and that the error is not due to some stupid
mistake  of mine.

Thanks a million,

Gianguido Cianci

PS:

Attached is a copy of the code I used:

/****************************************************/
#include <stdlib.h>
#include <iostream>
#include <math.h>
#include <gsl/gsl_math.h>
#include <gsl/gsl_monte.h>
#include <gsl/gsl_monte_plain.h>
#include <gsl/gsl_monte_miser.h>
#include <gsl/gsl_monte_vegas.h>



//double exact = 1.3932039296856768591842462603255;
double exact;

double
test (double *coord, size_t dim, void *params)
{
  double a(1.0), b(1.0);

  double omega(exp(-1.0*(coord[0]*coord[0]/a + coord[1]*coord[1]/b)));
  exact = M_PI/4.0*sqrt(a*b);
  return omega;

    
}

void
display_results2 (char *title, double result, double error)
{
  cout.precision(12);
  cout.setf(ios::scientific, ios::floatfield);
  cout << title << " ================" << endl 
       << "result = " << result << endl
       << "sigma = " << error << endl
       << "exact = " << exact << endl
       << "error = " << result-exact << " = " 
       <<fabs (result - exact) / error << " sigma" 
       << endl;

}




int
main (void)
{
  double res, err;

  double xl[2] = { 0, 0 };
  double xu[2] = { 100, 100};

  const gsl_rng_type *T;
  gsl_rng *r;

  gsl_monte_function G = { &test, 2, 0 };

  size_t calls = 500000;

  gsl_rng_env_setup ();

  T = gsl_rng_default;
  r = gsl_rng_alloc (T);

  {
    gsl_monte_plain_state *s = gsl_monte_plain_alloc (2);
    gsl_monte_plain_integrate (&G, xl, xu, 2, calls, r, s,
                               &res, &err);
    gsl_monte_plain_free (s);

    display_results2 ("\nplain", res, err);
  }

  {
    gsl_monte_miser_state *s = gsl_monte_miser_alloc (2);
    gsl_monte_miser_integrate (&G, xl, xu, 2, calls, r, s,
			       &res, &err);
    gsl_monte_miser_free (s);

    display_results2 ("\nmiser", res, err);
  }

  {
    gsl_monte_vegas_state *s = gsl_monte_vegas_alloc (2);

    gsl_monte_vegas_integrate (&G, xl, xu, 2, 10000, r, s,
                               &res, &err);
    display_results2 ("\nvegas warm-up", res, err);

    printf ("converging...\n");

    do
      {
        gsl_monte_vegas_integrate (&G, xl, xu, 2, calls/5, r, s,
                                   &res, &err);
        printf ("result = % .6f sigma = % .6f "
                "chisq/dof = %.1f\n", res, err, s->chisq);
      }
    while (fabs (s->chisq - 1.0) > 0.5);

    display_results2 ("\nvegas final", res, err);

    gsl_monte_vegas_free (s);
  }
  return 0;
}


/**************************************************************/

=======================================================================
Gianguido C. Cianci                      <gianguido.cianci@physics.org> 
Physics Department
Emory University
1510 Clifton Road NE                    Microscopy Lab:  (404) 712 8669
Atlanta, GA 30322                       Chemistry Lab:   (404) 712 8670
USA                                     Fax:             (404) 727 0873

=======================================================================



Index Nav: [Date Index] [Subject Index] [Author Index] [Thread Index]
Message Nav: [Date Prev] [Date Next] [Thread Prev] [Thread Next]