Class mcarlo_miser (o2scl)¶

O2scl : Class List

template<class func_t = multi_funct, class vec_t = boost::numeric::ublas::vector<double>, class rng_t = rng<>> class mcarlo_miser : public o2scl::mcarlo<multi_funct, boost::numeric::ublas::vector<double>, rng<>>¶

Multidimensional integration using the MISER Monte Carlo algorithm (GSL)

This class uses recursive stratified sampling to estimate the value of an integral over a hypercubic region.

By default the minimum number of calls to estimate the variance is 16 times the number of dimensions. This ratio is stored in calls_per_dim. By default the minimum number of calls per bisection is 32 times calls_per_dim times the number of dimensions. This ratio is stored in bisection_ratio. These ratios are employed by minteg_err().

Alternatively, the user can directly set these minimums by set_min_calls() and set_min_calls_per_bisection() and use miser_minteg_err() which ignores calls_per_dim and bisection_ratio.

If mcarlo::verbose is greater than zero, then the status of the integration is output at every level of bisection less than n_levels_out. If it is greater than 1, then the boundaries of the current region are also output. Finally, if it is greater than 2, a keypress is required after each output.

Based on [Press90].

Arrays which contain a value for each dimension

ubvector xmid¶: The current midpoint.

ubvector sigma_l¶: The left variance.

ubvector sigma_r¶: The right variance.

ubvector fmax_l¶: The maximum function value in the left half.

ubvector fmax_r¶: The maximum function value in the right half.

ubvector fmin_l¶: The minimum function value in the left half.

ubvector fmin_r¶: The minimum function value in the right half.

ubvector fsum_l¶: The sum in the left half.

ubvector fsum_r¶: The sum in the right half.

ubvector fsum2_l¶: The sum of the squares in the left half.

ubvector fsum2_r¶: The sum of the squares in the right half.

ubvector_size_t hits_l¶: The number of evaluation points in the left half.

ubvector_size_t hits_r¶: The number of evaluation points in the right half.

vec_t x¶: The most recent integration point.

inline mcarlo_miser()¶

inline virtual ~mcarlo_miser()¶

inline virtual int allocate(size_t ldim)¶: Allocate memory.

inline virtual int miser_minteg_err(func_t &func, size_t ndim, const vec_t &xl, const vec_t &xu, size_t calls, size_t level, double &res, double &err)¶: Integrate function func over the hypercube from \( x_i=\mathrm{xl}_i \) to \( x_i=\mathrm{xu}_i \) for \( 0<i< \) ndim-1.

Note

The values of min_calls and min_calls_per_bisection should be set before calling this function. The default values if not set are 100 and 3000 respectively, which correspond to the GSL default setting for a 6 dimensional problem.

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.

This function is just a wrapper to miser_minteg_err() which allocates the memory if necessary, sets min_calls and min_calls_per_bisection, calls miser_minteg_err(), and then frees the previously allocated memory.

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.

This function is just a wrapper to minteg_err() which allocates the memory, sets min_calls and min_calls_per_bisection, calls miser_minteg_err(), and then frees the previously allocated memory.

inline virtual const char *type()¶: Return string denoting type (“mcarlo_miser”)

inline virtual int estimate_corrmc(func_t &func, size_t ndim, const vec_t &xl, const vec_t &xu, size_t calls, double &res, double &err, const ubvector &lxmid, ubvector &lsigma_l, ubvector &lsigma_r)¶

Estimate the variance.

Idea for Future:: Remove the reference to GSL_POSINF and replace with a function parameter.

Public Types

typedef boost::numeric::ublas::vector<double> ubvector¶

typedef boost::numeric::ublas::vector<size_t> ubvector_size_t¶

Public Functions

inline int set_min_calls(size_t &mc)¶

Minimum number of calls to estimate the variance.

This is set by minteg() and minteg_err() to be calls_per_dim times the number of dimensions in the problem. The default value of calls_per_dim is 16 (which is the GSL default).

From GSL documentation:

This parameter specifies the minimum number of function calls
required for each estimate of the variance. If the number of
function calls allocated to the estimate using ESTIMATE_FRAC falls
below MIN_CALLS then MIN_CALLS are used instead.  This ensures
that each estimate maintains a reasonable level of accuracy.

inline int set_min_calls_per_bisection(size_t &mcb)¶

Minimum number of calls required to proceed with bisection.

This is set by minteg() and minteg_err() to be calls_per_dim times bisection_ratio times the number of dimensions in the problem. The default values give 512 times the number of dimensions (also the GSL default).

From GSL documentation:

This parameter specifies the minimum number of function calls
required to proceed with a bisection step.  When a recursive step
has fewer calls available than MIN_CALLS_PER_BISECTION it performs
a plain Monte Carlo estimate of the current sub-region and
terminates its branch of the recursion.

Public Members

size_t calls_per_dim¶: Number of calls per dimension (default 16)

size_t bisection_ratio¶: Factor to set min_calls_per_bisection (default 32)

double dither¶

Introduce random variation into bisection (default 0.0)

From GSL documentation:

This parameter introduces a random fractional variation of size
DITHER into each bisection, which can be used to break the
symmetry of integrands which are concentrated near the exact
center of the hypercubic integration region.  The default value of
dither is zero, so no variation is introduced. If needed, a
typical value of DITHER is 0.1.

double estimate_frac¶

Specify fraction of function calls for estimating variance (default 0.1)

From GSL documentation:

This parameter specifies the fraction of the currently available
number of function calls which are allocated to estimating the
variance at each recursive step. The default value is 0.1.

double alpha¶

How estimated variances for two sub-regions are combined (default 2.0)

The error handler will be called if this is less than zero.

From GSL documentation:

This parameter controls how the estimated variances for the two
sub-regions of a bisection are combined when allocating points.
With recursive sampling the overall variance should scale better
than 1/N, since the values from the sub-regions will be obtained
using a procedure which explicitly minimizes their variance.  To
accommodate this behavior the MISER algorithm allows the total
variance to depend on a scaling parameter \alpha,

\Var(f) = {\sigma_a \over N_a^\alpha} + {\sigma_b \over N_b^\alpha}.

The authors of the original paper describing MISER recommend the
value \alpha = 2 as a good choice, obtained from numerical
experiments, and this is used as the default value in this
implementation.

size_t n_levels_out¶: The number of recursive levels to output when verbose is greater than zero (default 5)

Protected Attributes

size_t min_calls¶: Minimum number of calls to estimate the variance.

size_t min_calls_per_bisection¶: Minimum number of calls required to proceed with bisection.

size_t dim¶: The number of dimensions of memory allocated.