Public Types | Public Member Functions | Public Attributes | Protected Member Functions | Protected Attributes | List of all members
o2scl::mcarlo_miser< func_t, vec_t, rng_t > Class Template Reference

Multidimensional integration using the MISER Monte Carlo algorithm (GSL) More...

#include <mcarlo_miser.h>

Inheritance diagram for o2scl::mcarlo_miser< func_t, vec_t, rng_t >:
o2scl::mcarlo< func_t, vec_t, rng_t > o2scl::inte_multi< func_t, vec_t >

Public Types

typedef boost::numeric::ublas::vector< double > ubvector
 
typedef boost::numeric::ublas::vector< size_t > ubvector_size_t
 

Public Member Functions

int set_min_calls (size_t &mc)
 Minimum number of calls to estimate the variance. More...
 
int set_min_calls_per_bisection (size_t &mcb)
 Minimum number of calls required to proceed with bisection. More...
 
virtual int allocate (size_t ldim)
 Allocate memory.
 
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. More...
 
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. More...
 
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. More...
 
virtual const char * type ()
 Return string denoting type ("mcarlo_miser")
 
- Public Member Functions inherited from o2scl::inte_multi< func_t, vec_t >
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")
 

Public Attributes

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) More...
 
double estimate_frac
 Specify fraction of function calls for estimating variance (default 0.1) More...
 
double alpha
 How estimated variances for two sub-regions are combined (default 2.0) More...
 
size_t n_levels_out
 The number of recursive levels to output when verbose is greater than zero (default 5)
 
- Public Attributes inherited from o2scl::mcarlo< func_t, vec_t, rng_t >
unsigned long n_points
 Number of integration points (default 1000)
 
rng_gsl_uniform_real rng_dist
 The random number distribution.
 
rng_t rng
 The random number generator.
 
- Public Attributes inherited from o2scl::inte_multi< func_t, vec_t >
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 $ 10^{-8} $).
 

Protected Member Functions

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. More...
 

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.
 
vec_t x
 The most recent integration point.
 
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.
 
- Protected Attributes inherited from o2scl::inte_multi< func_t, vec_t >
double interror
 The uncertainty for the last integration computation.
 

Detailed Description

template<class func_t = multi_funct11, class vec_t = boost::numeric::ublas::vector<double>, class rng_t = rng_gsl>
class o2scl::mcarlo_miser< func_t, vec_t, rng_t >

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 .

Definition at line 94 of file mcarlo_miser.h.

Member Function Documentation

◆ estimate_corrmc()

template<class func_t = multi_funct11, class vec_t = boost::numeric::ublas::vector<double>, class rng_t = rng_gsl>
virtual int o2scl::mcarlo_miser< func_t, vec_t, rng_t >::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 
)
inlineprotectedvirtual
Idea for Future:
Remove the reference to GSL_POSINF and replace with a function parameter.

Definition at line 254 of file mcarlo_miser.h.

◆ minteg()

template<class func_t = multi_funct11, class vec_t = boost::numeric::ublas::vector<double>, class rng_t = rng_gsl>
virtual double o2scl::mcarlo_miser< func_t, vec_t, rng_t >::minteg ( func_t &  func,
size_t  ndim,
const vec_t &  a,
const vec_t &  b 
)
inlinevirtual

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.

Reimplemented from o2scl::inte_multi< func_t, vec_t >.

Definition at line 697 of file mcarlo_miser.h.

◆ minteg_err()

template<class func_t = multi_funct11, class vec_t = boost::numeric::ublas::vector<double>, class rng_t = rng_gsl>
virtual int o2scl::mcarlo_miser< func_t, vec_t, rng_t >::minteg_err ( func_t &  func,
size_t  ndim,
const vec_t &  a,
const vec_t &  b,
double &  res,
double &  err 
)
inlinevirtual

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.

Implements o2scl::inte_multi< func_t, vec_t >.

Definition at line 674 of file mcarlo_miser.h.

◆ miser_minteg_err()

template<class func_t = multi_funct11, class vec_t = boost::numeric::ublas::vector<double>, class rng_t = rng_gsl>
virtual int o2scl::mcarlo_miser< func_t, vec_t, rng_t >::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 
)
inlinevirtual
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.

Definition at line 409 of file mcarlo_miser.h.

◆ set_min_calls()

template<class func_t = multi_funct11, class vec_t = boost::numeric::ublas::vector<double>, class rng_t = rng_gsl>
int o2scl::mcarlo_miser< func_t, vec_t, rng_t >::set_min_calls ( size_t &  mc)
inline

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. 

Definition at line 175 of file mcarlo_miser.h.

◆ set_min_calls_per_bisection()

template<class func_t = multi_funct11, class vec_t = boost::numeric::ublas::vector<double>, class rng_t = rng_gsl>
int o2scl::mcarlo_miser< func_t, vec_t, rng_t >::set_min_calls_per_bisection ( size_t &  mcb)
inline

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. 

Definition at line 196 of file mcarlo_miser.h.

Member Data Documentation

◆ alpha

template<class func_t = multi_funct11, class vec_t = boost::numeric::ublas::vector<double>, class rng_t = rng_gsl>
double o2scl::mcarlo_miser< func_t, vec_t, rng_t >::alpha

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.

Definition at line 158 of file mcarlo_miser.h.

◆ dither

template<class func_t = multi_funct11, class vec_t = boost::numeric::ublas::vector<double>, class rng_t = rng_gsl>
double o2scl::mcarlo_miser< func_t, vec_t, rng_t >::dither

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.

Definition at line 121 of file mcarlo_miser.h.

◆ estimate_frac

template<class func_t = multi_funct11, class vec_t = boost::numeric::ublas::vector<double>, class rng_t = rng_gsl>
double o2scl::mcarlo_miser< func_t, vec_t, rng_t >::estimate_frac

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.

Definition at line 133 of file mcarlo_miser.h.


The documentation for this class was generated from the following file:

Documentation generated with Doxygen. Provided under the GNU Free Documentation License (see License Information).