Compressing Large-Scale Wave Propagation Models via Phase-Preconditioned Rational Krylov Subspaces

Abstract

Rational Krylov subspace (RKS) techniques are well-established and powerful tools for projection-based model reduction of time-invariant dynamic systems. For hyperbolic wavefield problems, such techniques perform well in configurations where only a few modes contribute to the field. RKS methods, however, are fundamentally limited by the Nyquist-Shannon sampling rate, making them unsuitable for the approximation of wavefields in configuration characterized by large travel times and propagation distances, since wavefield responses in such configurations are highly oscillatory in the frequency-domain. To overcome this limitation, we propose to precondition the RKSs by factoring out the rapidly varying frequency-domain field oscillations. The remaining amplitude functions are generally slowly varying functions of source position and spatial coordinate and allow for a significant compression of the approximation subspace. Our one-dimensional analysis together with numerical experiments for large scale 2D acoustic models show superior approximation properties of preconditioned RKS compared with the standard RKS model-order reduction. The preconditioned RKS results in a reduction of the frequency sampling well below the Nyquist-Shannon rate, a weak dependence of the RKS size on the number of inputs and outputs for multiple-input/multiple-output (MIMO) problems, and, most importantly, in a significant coarsening of the finite-difference grid used to generate the RKS. A prototype implementation indicates that the preconditioned RKS algorithm is competitive in the modern high performance computing environment.AMS subject classifications. 65M60, 65M80, 93B11, 37M05, 78A05 1. Introduction. Numerical modeling of wave propagation is fundamental to many applications in design optimization and wavefield imaging. In the oil and gas industry, for instance, the solution of the Maxwell equations is required to invert electromagnetic measurements, while in seismic imaging the solution to the elastodynamic wave equation is needed to ultimately image the subsurface of the Earth.Finite difference discretization of the governing wave equations leads to large-scale linear systems, whose solution is computationally intense. Imaging and optimization often use multiple frequencies, sources, and receivers, which leads to systems that need to be evaluated for multiple right-hand sides, time-steps or frequencies, depending on whether the problem is solved in the time-or frequency-domain. Therefore, these so-called multiple-input/multiple-output (MIMO) systems have a high demand on memory and computational power, causing long runtimes. To be more specific, let us consider a surface seismic imaging problem in a k dimensional space (1 ≤ k ≤ 3), with maximal propagation distance of N wavelengths. This would require the solution of a discretized system with O(N k ) state variables, O(N k−1 ) sources and receivers, and O(N ) frequencies or time steps [22]. Model-order reduction aims to reduce the complexity and computational burden of large-scale pro...

A conjugate-gradient-type rational Krylov subspace method for ill-posed problems

“…As a possible application, rational Krylov subspaces have been successfully used to accelerate computations related to seismic imaging (e.g. [4,19,31,32]), which is known to be an ill-posed inverse problem.…”