Multidimensional integration using Vegas Monte Carlo (GSL) More...
#include <mcarlo_vegas.h>
The output options are a little different than the original GSL routine. The default setting of mcarlo::verbose is 0, which turns off all output. A verbose value of 1 prints summary information about the weighted average and final result, while a value of 2 also displays the grid coordinates. A value of 3 prints information from the rebinning procedure for each iteration.
Some original documentation from GSL:
The input coordinates are x[j], with upper and lower limits xu[j] and xl[j]. The integration length in the j-th direction is delx[j]. Each coordinate x[j] is rescaled to a variable y[j] in the range 0 to 1. The range is divided into bins with boundaries xi[i][j], where i=0 corresponds to y=0 and i=bins to y=1. The grid is refined (ie, bins are adjusted) using d[i][j] which is some variation on the squared sum. A third parameter used in defining the real coordinate using random numbers is called z. It ranges from 0 to bins. Its integer part gives the lower index of the bin into which a call is to be placed, and the remainder gives the location inside the bin. When stratified sampling is used the bins are grouped into boxes, and the algorithm allocates an equal number of function calls to each box. The variable alpha controls how "stiff" the rebinning algorithm is. alpha = 0 means never change the grid. Alpha is typically set between 1 and 2.
Based on Lepage78 . The current version of the algorithm was described in the Cornell preprint CLNS-80/447 of March,
Definition at line 115 of file mcarlo_vegas.h.
Public Types | |
typedef boost::numeric::ublas::vector< double > | ubvector |
typedef boost::numeric::ublas::vector< size_t > | ubvector_size_t |
typedef boost::numeric::ublas::vector< int > | ubvector_int |
Integration mode (default is mode_importance) | |
int | mode |
double | result |
Result from last iteration. | |
double | sigma |
Uncertainty from last iteration. | |
double | alpha |
The stiffness of the rebinning algorithm (default 1.5) More... | |
unsigned int | iterations |
Set the number of iterations (default 5) | |
double | chisq |
The chi-squared per degree of freedom for the weighted estimate of the integral. More... | |
std::ostream * | outs |
The output stream to send output information (default std::cout ) | |
static const int | mode_importance =1 |
static const int | mode_importance_only =0 |
static const int | mode_stratified =-1 |
static const size_t | bins_max =50 |
Maximum number of bins. | |
size_t | dim |
Number of dimensions. | |
unsigned int | bins |
Number of bins. | |
unsigned int | boxes |
Number of boxes. | |
ubvector | xi |
Boundaries for each bin. | |
ubvector | xin |
Storage for grid refinement. | |
ubvector | delx |
The iteration length in each direction. | |
ubvector | weight |
The weight for each bin. | |
double | vol |
The volume of the current bin. | |
ubvector_int | bin |
The bins for each direction. | |
ubvector_int | box |
The boxes for each direction. | |
ubvector | d |
Distribution. | |
Scratch variables preserved between calls to | |
double | jac |
double | wtd_int_sum |
double | sum_wgts |
double | chi_sum |
unsigned int | it_start |
The starting iteration number. | |
unsigned int | it_num |
The total number of iterations. | |
unsigned int | samples |
Number of samples for computing chi squared. | |
unsigned int | calls_per_box |
Number of function calls per box. | |
vec_t | x |
Point for function evaluation. | |
virtual void | init_box_coord (ubvector_int &boxt) |
Initialize box coordinates. | |
int | change_box_coord (ubvector_int &boxt) |
Change box coordinates. More... | |
virtual void | init_grid (const vec_t &xl, const vec_t &xu, size_t ldim) |
Initialize grid. | |
virtual void | reset_grid_values () |
Reset grid values. | |
void | accumulate_distribution (ubvector_int &lbin, double y) |
Add the most recently generated result to the distribution. More... | |
void | random_point (vec_t &lx, ubvector_int &lbin, double &bin_vol, const ubvector_int &lbox, const vec_t &xl, const vec_t &xu) |
Generate a random position in a given box. More... | |
virtual void | resize_grid (unsigned int lbins) |
Resize the grid. | |
virtual void | refine_grid () |
Refine the grid. | |
virtual void | print_lim (const vec_t &xl, const vec_t &xu, unsigned long ldim) |
Print limits of integration. | |
virtual void | print_head (unsigned long num_dim, unsigned long calls, unsigned int lit_num, unsigned int lbins, unsigned int lboxes) |
Print header. | |
virtual void | print_res (unsigned int itr, double res, double err, double cum_res, double cum_err, double chi_sq) |
Print results. | |
virtual void | print_dist (unsigned long ldim) |
Print distribution. | |
virtual void | print_grid (unsigned long ldim) |
Print grid. | |
mcarlo_vegas () | |
virtual int | allocate (size_t ldim) |
Allocate memory. | |
virtual int | vegas_minteg_err (int stage, func_t &func, size_t ndim, const vec_t &xl, const vec_t &xu, double &res, double &err) |
Integrate function func from x=a to x=b. More... | |
virtual | ~mcarlo_vegas () |
virtual int | minteg_err (func_t &func, size_t ndim, const vec_t &a, const vec_t &b, double &res, double &err) |
Integrate function func from x=a to x=b. | |
virtual double | minteg (func_t &func, size_t ndim, const vec_t &a, const vec_t &b) |
Integrate function func over the hypercube from ![]() ![]() ![]() | |
virtual const char * | type () |
Return string denoting type ("mcarlo_vegas") | |
Additional Inherited Members | |
![]() | |
virtual double | minteg (multi_funct &func, size_t ndim, const boost::numeric::ublas::vector< double > &a, const boost::numeric::ublas::vector< double > &b) |
Integrate function func over the hypercube from ![]() ![]() ![]() | |
virtual int | minteg_err (multi_funct &func, size_t ndim, const boost::numeric::ublas::vector< double > &a, const boost::numeric::ublas::vector< double > &b, double &res, double &err)=0 |
Integrate function func over the hypercube from ![]() ![]() ![]() | |
double | get_error () |
Return the error in the result from the last call to minteg() or minteg_err() More... | |
const char * | type () |
Return string denoting type ("inte_multi") | |
![]() | |
unsigned long | n_points |
Number of integration points (default 1000) | |
rng_gsl_uniform_real | rng_dist |
The random number distribution. | |
rng_gsl | rng |
The random number generator. | |
![]() | |
bool | err_nonconv |
If true, call the error handler if the routine does not "converge". | |
int | verbose |
Verbosity. | |
double | tol_rel |
The maximum "uncertainty" in the value of the integral (default ![]() | |
![]() | |
double | interror |
The uncertainty for the last integration computation. | |
|
inlineprotected |
This is among the member functions that is not virtual because it is part of the innermost loop.
Definition at line 282 of file mcarlo_vegas.h.
|
inlineprotected |
Steps through the box coordinates, e.g.
{0,0}, {0,1}, {0,2}, {0,3}, {1,0}, {1,1}, {1,2}, ...
This is among the member functions that is not virtual because it is part of the innermost loop.
Definition at line 230 of file mcarlo_vegas.h.
|
inlineprotected |
Use the random number generator to return a random position x in a given box. The value of bin gives the bin location of the random position (there may be several bins within a given box)
This is among the member functions that is not virtual because it is part of the innermost loop.
Definition at line 301 of file mcarlo_vegas.h.
|
inlinevirtual |
Original documentation from GSL:
Normally, stage = 0
which begins with a new uniform grid and empty weighted average. Calling vegas with stage = 1
retains the grid from the previous run but discards the weighted average, so that one can "tune" the grid using a relatively small number of points and then do a large run with stage = 1
on the optimized grid. Setting stage = 2
keeps the grid and the weighted average from the previous run, but may increase (or decrease) the number of histogram bins in the grid depending on the number of calls available. Choosing stage = 3
enters at the main loop, so that nothing is changed, and is equivalent to performing additional iterations in a previous call.
Should stage be passed by reference?
There was an update between gsl-1.12 and 1.15 which has not been implemented here yet.
Definition at line 608 of file mcarlo_vegas.h.
double o2scl::mcarlo_vegas< func_t, vec_t, rng_t >::alpha |
This usual range is between 1 and 2.
Definition at line 141 of file mcarlo_vegas.h.
double o2scl::mcarlo_vegas< func_t, vec_t, rng_t >::chisq |
After an integration, this should be close to 1. If it is not, then this indicates that the values of the integral from different iterations are inconsistent, and the error may be underestimated. Further iterations of the algorithm may enable one to obtain more reliable results.
Definition at line 155 of file mcarlo_vegas.h.
Documentation generated with Doxygen. Provided under the
GNU Free Documentation License (see License Information).