Class mcarlo_vegas (o2scl)¶
-
template<class func_t = multi_funct, class vec_t = boost::numeric::ublas::vector<double>, class rng_t = rng<>>
class mcarlo_vegas : public o2scl::mcarlo<multi_funct, boost::numeric::ublas::vector<double>, rng<>>¶ Multidimensional integration using Vegas Monte Carlo (GSL)
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.
Todo
CLass mcarlo_vegas: Mode = importance only doesn’t give the same answer as GSL yet.
- Idea for Future:
Prettify the verbose output
- Idea for Future:
Allow the user to get information about the how the sampling was done, possibly by converting the bins and boxes into a structure or class.
- Idea for Future:
Allow the user to change the maximum number of bins.
Based on [Lepage78]
The current version of the algorithm was described in the Cornell preprint CLNS-80/447 of March,
The GSL code follows most closely the C version by D. R. Yennie, coded in 1984.
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)
This usual range is between 1 and 2.
-
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.
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.
-
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.
-
double vol¶
The volume of the current bin.
-
ubvector_int bin¶
The bins for each direction.
-
ubvector_int box¶
The boxes for each direction.
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.
-
inline virtual void init_box_coord(ubvector_int &boxt)¶
Initialize box coordinates.
-
inline int change_box_coord(ubvector_int &boxt)¶
Change box coordinates.
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.
-
inline virtual void reset_grid_values()¶
Reset grid values.
-
inline void accumulate_distribution(ubvector_int &lbin, double y)¶
Add the most recently generated result to the distribution.
This is among the member functions that is not virtual because it is part of the innermost loop.
-
inline 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.
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.
-
inline virtual void resize_grid(unsigned int lbins)¶
Resize the grid.
-
inline virtual void refine_grid()¶
Refine the grid.
-
inline virtual void print_lim(const vec_t &xl, const vec_t &xu, unsigned long ldim)¶
Print limits of integration.
-
inline virtual void print_head(unsigned long num_dim, unsigned long calls, unsigned int lit_num, unsigned int lbins, unsigned int lboxes)¶
Print header.
-
inline virtual void print_res(unsigned int itr, double res, double err, double cum_res, double cum_err, double chi_sq)¶
Print results.
-
inline virtual void print_dist(unsigned long ldim)¶
Print distribution.
-
inline virtual void print_grid(unsigned long ldim)¶
Print grid.
-
inline mcarlo_vegas()¶
-
inline virtual int allocate(size_t ldim)¶
Allocate memory.
-
inline 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.Original documentation from GSL:
Normally,
stage = 0
which begins with a new uniform grid and empty weighted average. Calling vegas withstage = 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 withstage = 1
on the optimized grid. Settingstage = 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. Choosingstage = 3
enters at the main loop, so that nothing is changed, and is equivalent to performing additional iterations in a previous call.Todo
Function mcarlo_vegas::vegas_minteg_err():
Should stage be passed by reference?
There was an update between gsl-1.12 and 1.15 which has not been implemented here yet.
-
inline virtual ~mcarlo_vegas()¶
-
inline 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.The result of the integral is stored in
res
and the error estimate inerr
.
-
inline virtual double minteg(func_t &func, size_t ndim, const vec_t &a, const vec_t &b)¶
Integrate function
func
over the hypercube from \( x_i=a_i \) to \( x_i=b_i \) for \( 0<i< \) ndim-1.
-
inline virtual const char *type()¶
Return string denoting type (“mcarlo_vegas”)