Multilevel Monte Carlo methods for stochastic convection-diffusion eigenvalue problems

TL;DR

This work develops multilevel Monte Carlo (MLMC) methods to estimate the expected smallest eigenvalue $\mathbb{E}[\lambda]$ of a stochastic convection–diffusion operator, using a variational finite element framework and two discretizations (Galerkin FEM and SUPG). It provides a rigorous error/complexity analysis, incorporating Rayleigh quotient iteration and implicitly restarted Arnoldi, and introduces two extensions: homotopy MLMC and multilevel quasi-Monte Carlo (MLQMC) to further enhance efficiency. Numerical results on three test cases demonstrate the effectiveness of SUPG in stabilizing convection-dominated problems, verify the predicted MLMC rates, and show that MLQMC offers the best computational efficiency among the tested strategies. Overall, the paper advances uncertainty quantification for stochastic eigenvalue problems by delivering practical MLMC strategies with solid theoretical guarantees and compelling numerical evidence.

Abstract

We develop new multilevel Monte Carlo (MLMC) methods to estimate the expectation of the smallest eigenvalue of a stochastic convection-diffusion operator with random coefficients. The MLMC method is based on a sequence of finite element (FE) discretizations of the eigenvalue problem on a hierarchy of increasingly finer meshes. For the discretized, algebraic eigenproblems we use both the Rayleigh quotient (RQ) iteration and implicitly restarted Arnoldi (IRA), providing an analysis of the cost in each case. By studying the variance on each level and adapting classical FE error bounds to the stochastic setting, we are able to bound the total error of our MLMC estimator and provide a complexity analysis. As expected, the complexity bound for our MLMC estimator is superior to plain Monte Carlo. To improve the efficiency of the MLMC further, we exploit the hierarchy of meshes and use coarser approximations as starting values for the eigensolvers on finer ones. To improve the stability of the MLMC method for convection-dominated problems, we employ two additional strategies. First, we consider the streamline upwind Petrov--Galerkin formulation of the discrete eigenvalue problem, which allows us to start the MLMC method on coarser meshes than is possible with standard FEs. Second, we apply a homotopy method to add stability to the eigensolver for each sample. Finally, we present a multilevel quasi-Monte Carlo method that replaces Monte Carlo with a quasi-Monte Carlo (QMC) rule on each level. Due to the faster convergence of QMC, this improves the overall complexity. We provide detailed numerical results comparing our different strategies to demonstrate the practical feasibility of the MLMC method in different use cases. The results support our complexity analysis and further demonstrate the superiority over plain Monte Carlo in all cases.

Figures (11)

• Figure 1: The first 20 computed eigenvalues of the SUPG (left) and FEM (right) discretizations of the convection-diffusion problem for $\kappa(\mathbf{x})=1$ and $\mathbf{a}=[50,0]^T$ using mesh sizes $h=2^{-3}, 2^{-4}, 2^{-5}$.
• Figure 2: MLMC method using tgRQI for Test Case I with $\mathbf{a}=[20;0]^T$ and Galerkin FEM: (a) Mean (blue) and variance (red) of the eigenvalue $\lambda_\ell$ (dashed) and of $\lambda_\ell-\lambda_{\ell-1}$ (solid); (b) computational times for one multilevel difference (blue) and average number of Rayleigh quotient iterations (red) on each level. Where shown, the error bars represent $\pm$ one standard deviation.
• Figure 3: MLMC method using homotopy and tgRQI for Test Case I with $\mathbf{a}=[20;0]^T$ and Galerkin FEM: (a) Mean (blue) and variance (red) of the eigenvalue $\lambda_\ell$ (dashed) and of $\lambda_\ell-\lambda_{\ell-1}$ (solid); (b) computational times for one multilevel difference (blue) and average number of RQIs (red) on each level. Where shown, the error bars represent $\pm$ one standard deviation.
• Figure 4: MLMC method using IRAr for Test Case I with $\mathbf{a}=[20;0]^T$ and Galerkin FEM, both without (a) and with (b) homotopy: average computational cost (blue) and average number of matrix-vector products (red) per sample of $\lambda_\ell-\lambda_{\ell-1}$. The error bars represent $\pm$ one standard deviation.
• Figure 5: CPU time vs. root mean square error of all estimators in Test Case I.

Theorems & Definitions (15)

• Proposition 1
• proof
• Theorem 1
• proof
• Theorem 2
• proof
• Theorem 3
• Corollary 1: Order of convergence
• proof
• Remark 1