Implementation and analysis of an adaptive multilevel Monte Carlo algorithm

We present an adaptive multilevel Monte Carlo (MLMC) method for weak approximations of solutions to Itô stochastic differential equations (SDE). The work [Oper. Res. 56 (2008), 607-617] proposed and analyzed an MLMC method based on a hierarchy of uniform time discretizations and control variates to reduce the computational effort required by a single level Euler-Maruyama Monte Carlo method from ( TOL -3 )${{{\mathcal {O}}({\mathrm {TOL}}^{-3})}}$ to ( TOL -2 log( TOL -1 ) 2 )${{{\mathcal {O}}({\mathrm {TOL}}^{-2}\log ({\mathrm {TOL}}^{-1})^{2})}}$ for a mean square error of ( TOL 2 )${{{\mathcal {O}}({\mathrm {TOL}}^2)}}$ . Later, the work [Lect. Notes Comput. Sci. Eng. 82, Springer-Verlag, Berlin (2012), 217-234] presented an MLMC method using a hierarchy of adaptively refined, non-uniform time discretizations, and, as such, it may be considered a generalization of the uniform time discretization MLMC method. This work improves the adaptive MLMC algorithms presented in [Lect. Notes Comput. Sci. Eng. 82, Springer-Verlag, Berlin (2012), 217-234] and it also provides mathematical analysis of the improved algorithms. In particular, we show that under some assumptions our adaptive MLMC algorithms are asymptotically accurate and essentially have the correct complexity but with improved control of the complexity constant factor in the asymptotic analysis. Numerical tests include one case with singular drift and one with stopped diffusion, where the complexity of a uniform single level method is ( TOL -4 )${{{\mathcal {O}}({\mathrm {TOL}}^{-4})}}$ . For both these cases the results confirm the theory, exhibiting savings in the computational cost for achieving the accuracy ( TOL )${{{\mathcal {O}}({\mathrm {TOL}})}}$ from ( TOL -3 )${{{\mathcal {O}}({\mathrm {TOL}}^{-3})}}$ for the adaptive single level algorithm to essentially ( TOL -2 log( TOL -1 ) 2 )${{{\mathcal {O}}({\mathrm {TOL}}^{-2}\log ({\mathrm {TOL}}^{-1})^2)}}$ for the adaptive MLMC algorith

Rare Event Simulation and Splitting for Discontinuous Random Variables

Multilevel Splitting methods, also called Sequential Monte-Carlo or \emph{Subset Simulation}, are widely used methods for estimating extreme probabilities of the form $P[S(\mathbf{U}) > q]$ where $S$ is a deterministic real-valued function and $\mathbf{U}$ can be a random finite- or infinite-dimensional vector. Very often, $X := S(\mathbf{U})$ is supposed to be a continuous random variable and a lot of theoretical results on the statistical behaviour of the estimator are now derived with this hypothesis. However, as soon as some threshold effect appears in $S$ and/or $\mathbf{U}$ is discrete or mixed discrete/continuous this assumption does not hold any more and the estimator is not consistent. In this paper, we study the impact of discontinuities in the \emph{cdf} of $X$ and present three unbiased \emph{corrected} estimators to handle them. These estimators do not require to know in advance if $X$ is actually discontinuous or not and become all equal if $X$ is continuous. Especially, one of them has the same statistical properties in any case. Efficiency is shown on a 2-D diffusive process as well as on the \emph{Boolean SATisfiability problem} (SAT).Comment: 16 pages (12 + Appendix 4 pages), 6 figure

Optimization of mesh hierarchies in Multilevel Monte Carlo samplers

We perform a general optimization of the parameters in the Multilevel Monte Carlo (MLMC) discretization hierarchy based on uniform discretization methods with general approximation orders and computational costs. We optimize hierarchies with geometric and non-geometric sequences of mesh sizes and show that geometric hierarchies, when optimized, are nearly optimal and have the same asymptotic computational complexity as non-geometric optimal hierarchies. We discuss how enforcing constraints on parameters of MLMC hierarchies affects the optimality of these hierarchies. These constraints include an upper and a lower bound on the mesh size or enforcing that the number of samples and the number of discretization elements are integers. We also discuss the optimal tolerance splitting between the bias and the statistical error contributions and its asymptotic behavior. To provide numerical grounds for our theoretical results, we apply these optimized hierarchies together with the Continuation MLMC Algorithm. The first example considers a three-dimensional elliptic partial differential equation with random inputs. Its space discretization is based on continuous piecewise trilinear finite elements and the corresponding linear system is solved by either a direct or an iterative solver. The second example considers a one-dimensional It\^o stochastic differential equation discretized by a Milstein scheme

Multilevel Sparse Grid Methods for Elliptic Partial Differential Equations with Random Coefficients

Stochastic sampling methods are arguably the most direct and least intrusive means of incorporating parametric uncertainty into numerical simulations of partial differential equations with random inputs. However, to achieve an overall error that is within a desired tolerance, a large number of sample simulations may be required (to control the sampling error), each of which may need to be run at high levels of spatial fidelity (to control the spatial error). Multilevel sampling methods aim to achieve the same accuracy as traditional sampling methods, but at a reduced computational cost, through the use of a hierarchy of spatial discretization models. Multilevel algorithms coordinate the number of samples needed at each discretization level by minimizing the computational cost, subject to a given error tolerance. They can be applied to a variety of sampling schemes, exploit nesting when available, can be implemented in parallel and can be used to inform adaptive spatial refinement strategies. We extend the multilevel sampling algorithm to sparse grid stochastic collocation methods, discuss its numerical implementation and demonstrate its efficiency both theoretically and by means of numerical examples

Large deviations principle for the Adaptive Multilevel Splitting Algorithm in an idealized setting

The Adaptive Multilevel Splitting (AMS) algorithm is a powerful and versatile method for the simulation of rare events. It is based on an interacting (via a mutation-selection procedure) system of replicas, and depends on two integer parameters: n $\in$ N * the size of the system and the number k $\in$ {1, . . . , n -- 1} of the replicas that are eliminated and resampled at each iteration. In an idealized setting, we analyze the performance of this algorithm in terms of a Large Deviations Principle when n goes to infinity, for the estimation of the (small) probability P(X \textgreater{} a) where a is a given threshold and X is real-valued random variable. The proof uses the technique introduced in [BLR15]: in order to study the log-Laplace transform, we rely on an auxiliary functional equation. Such Large Deviations Principle results are potentially useful to study the algorithm beyond the idealized setting, in particular to compute rare transitions probabilities for complex high-dimensional stochastic processes

Computation of Electromagnetic Fields Scattered From Objects With Uncertain Shapes Using Multilevel Monte Carlo Method

Computational tools for characterizing electromagnetic scattering from objects with uncertain shapes are needed in various applications ranging from remote sensing at microwave frequencies to Raman spectroscopy at optical frequencies. Often, such computational tools use the Monte Carlo (MC) method to sample a parametric space describing geometric uncertainties. For each sample, which corresponds to a realization of the geometry, a deterministic electromagnetic solver computes the scattered fields. However, for an accurate statistical characterization the number of MC samples has to be large. In this work, to address this challenge, the continuation multilevel Monte Carlo (CMLMC) method is used together with a surface integral equation solver. The CMLMC method optimally balances statistical errors due to sampling of the parametric space, and numerical errors due to the discretization of the geometry using a hierarchy of discretizations, from coarse to fine. The number of realizations of finer discretizations can be kept low, with most samples computed on coarser discretizations to minimize computational cost. Consequently, the total execution time is significantly reduced, in comparison to the standard MC scheme.Comment: 25 pages, 10 Figure