Vegas algorithm

MCIntegration.Vegas.montecarloMethod
function montecarlo(config::Configuration{Ni,V,P,O,T}, integrand::Function, neval,
    verbose=0, timer=[], debug=false;
    measure::Union{Nothing,Function}=nothing, measurefreq::Int=1) where {Ni,V,P,O,T}

This function implements the Vegas algorithm, a Monte Carlo method specifically designed for multi-dimensional integration. The underlying principle and methodology of the algorithm can be explored further in the Vegas documentation.

Overview

The Vegas algorithm employs an importance sampling scheme. For a one-dimensional integral with the integrand $f(x)$, the algorithm constructs an optimized distribution $\rho(x)$ that approximates the integrand as closely as possible (a strategy known as the Vegas map trick; refer to Dist.Continuous for more details).

The variable $x$ is then sampled using the distribution $\rho(x)$, and the integral is estimated by averaging the estimator $f(x)/\rho(x)$.

Note

  • If there are multiple integrals, all of them are sampled and measured at each Monte Carlo step.
  • This algorithm is particularly efficient for low-dimensional integrations but might be less efficient and robust than the Markov-chain Monte Carlo solvers for high-dimensional integrations.

Arguments

  • integrand: A user-defined function evaluating the integrand. The function should be either integrand(var, config) or integrand(var, weights, config) depending on whether inplace is false or true respectively. Here, var are the random variables and weights is an output array to store the calculated weights. The last parameter passes the MC Configuration struct to the integrand, so that user has access to userdata, etc.

  • measure: An optional user-defined function to accumulate the integrand weights into the observable. The function signature should be measure(var, obs, relative_weights, config). Here, obs is a vector of observable values for each component of the integrand and relative_weights are the weights calculated from the integrand multiplied by the probability of the corresponding variables.

The following are the snippets of the integrand and measure functions:

function integrand(var, config)
    # calculate your integrand values
    # return integrand1, integrand2, ...
end
function measure(var, obs, weights, config)
    # accumulates the weight into the observable
    # For example,
    # obs[1] = weights[1] # integral 1
    # obs[2] = weights[2] # integral 2
    # ...
end

Examples

The following command calls the Vegas solver,

julia> integrate((x, c)->(x[1]^2+x[2]^2); var = Continuous(0.0, 1.0), dof = [[2,],], verbose=-1, solver=:vegas)
Integral 1 = 0.667203631824444 ± 0.0005046485925614018   (reduced chi2 = 1.46)
source