Multilevel approximation of Gaussian random fields: Fast simulation

We propose and analyze several multilevel algorithms for the fast simulation of possibly nonstationary Gaussian random fields (GRFs) indexed, for example, by the closure of a bounded domain [Formula: see text] or, more generally, by a compact metric space [Formula: see text] such as a compact [Formula: see text]-manifold [Formula: see text]. A colored GRF [Formula: see text], admissible for our algorithms, solves the stochastic fractional-order equation [Formula: see text] for some [Formula: see text], where [Formula: see text] is a linear, local, second-order elliptic self-adjoint differential operator in divergence form and [Formula: see text] is white noise on [Formula: see text]. We thus consider GRFs on [Formula: see text] with covariance operators of the form [Formula: see text]. The proposed algorithms numerically approximate samples of [Formula: see text] on nested sequences [Formula: see text] of regular, simplicial partitions [Formula: see text] of [Formula: see text] and [Formula: see text], respectively. Work and memory to compute one approximate realization of the GRF [Formula: see text] on the triangulation [Formula: see text] of [Formula: see text] with consistency [Formula: see text], for some consistency order [Formula: see text], scale essentially linearly in [Formula: see text], independent of the possibly low regularity of the GRF. The algorithms are based on a sinc quadrature for an integral representation of (the application of) the negative fractional-order elliptic “coloring” operator [Formula: see text] to white noise [Formula: see text]. For the proposed numerical approximation, we prove bounds of the computational cost and the consistency error in various norms.