Vegas algorithm
MCIntegration.Vegas.montecarlo — Methodfunction 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 eitherintegrand(var, config)orintegrand(var, weights, config)depending on whetherinplaceisfalseortruerespectively. Here,varare the random variables andweightsis an output array to store the calculated weights. The last parameter passes the MCConfigurationstruct 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 bemeasure(var, obs, relative_weights, config). Here,obsis a vector of observable values for each component of the integrand andrelative_weightsare 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, ...
endfunction 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
# ...
endExamples
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)