MDFEM: Multivariate decomposition finite element method for elliptic PDEs with lognormal diffusion coefficients using higher-order QMC and FEM

Abstract

We introduce the \emph{multivariate decomposition finite element method} (MDFEM) for elliptic PDEs with lognormal diffusion coefficients, that is, when the diffusion coefficient has the form $a=\exp(Z)$ where $Z$ is a Gaussian random field defined by an infinite series expansion $Z(\bsy) = \sum_{j \ge 1} y_j \, \phi_j$ with $y_j \sim \calN(0,1)$ and a given sequence of functions $\{\phi_j\}_{j \ge 1}$. We use the MDFEM to approximate the expected value of a linear functional of the solution of the PDE which is an infinite-dimensional integral over the parameter space. The proposed algorithm uses the \emph{multivariate decomposition method} (MDM) to compute the infinite-dimensional integral by a decomposition into finite-dimensional integrals, which we resolve using \emph{quasi-Monte Carlo} (QMC) methods, and for which we use the \emph{finite element method} (FEM) to solve different instances of the PDE. We develop higher-order quasi-Monte Carlo rules for integration over the finite-di\-men\-si\-onal Euclidean space with respect to the Gaussian distribution by use of a truncation strategy. By linear transformations of interlaced polynomial lattice rules from the unit cube to a multivariate box of the Euclidean space we achieve higher-order convergence rates for functions belonging to a class of \emph{anchored Gaussian Sobolev spaces} while taking into account the truncation error. These cubature rules are then used in the MDFEM algorithm. Under appropriate conditions, the MDFEM achieves higher-order convergence rates in term of error versus cost, i.e., to achieve an accuracy of $O(\epsilon)$ the computational cost is $O(\epsilon^{-1/\lambda-\dd/\lambda}) = O(\epsilon^{-(p^* + \dd/\tau)/(1-p^*)})$ where $\epsilon^{-1/\lambda}$ and $\epsilon^{-\dd/\lambda}$ are respectively the cost of the quasi-Monte Carlo cubature and the finite element approximations, with $\dd = d \, (1+\ddelta)$ for some $\ddelta \ge 0$ and $d$ the physical dimension, and $0 < p^* \le (2 + \dd/\tau)^{-1}$ is a parameter representing the sparsity of $\{\phi_j\}_{j \ge 1}$.

Scaled lattice rules for integration on ℝ^{𝕕} achieving higher-order convergence with error analysis in terms of orthogonal projections onto periodic spaces

“…), obtained by the Bernoulli polynomial method on the box [a, b], cf. (21), is equivalent to the orthogonal projection on H(K…”

“…Proof. (i) Using the definition of F [a,b] , see (21), the bound on the maximum magnitude of the Bernoulli polynomials (6) and the decay condition (22) we obtain…”

“…Our major concern however will be how to analyse the error of integrating the non-periodic integrand on the box [a, b] by a lattice rule. Different approaches exist in the literature to apply quasi-Monte Carlo methods designed for integration over the unit cube [0, 1] d to integration over R d , e.g., [5,22,20,17,14,21,16]. In [20] analytic functions with exponential decay in the physical and in the Fourier space are studied and exponential convergence is obtained using a regular grid with appropriate scaling in each direction.…”

“…Therein, higher-order convergence is achieved using linearly transformed higher-order digital nets on a box. Similarly, in [21] higher-order convergence is achieved with interlaced polynomial lattice rules where a different type of function spaces, namely anchored Sobolev spaces, are considered. In [22] randomly shifted lattice rules are studied in combination with an inverse mapping of the considered cumulative density function, and first order convergence is achieved.…”

Goal-Oriented Adaptive Finite Element Multilevel Monte Carlo with Convergence Rates

“…N (0, 1) random variables. For the well-posedness and solution regularity for problem (1) with a lognormal coefficient field (2), we refer to [5,8,17,36,68]. There, well posedness requires control on the pathwise coercivity cf.…”