Scalable matrix-free adaptive product-convolution approximation for locally translation-invariant operators

TL;DR

This paper introduces a scalable, adaptive, matrix-free approximation method for locally translation-invariant operators, enabling efficient application, inversion, and preconditioning in high-dimensional PDE-related problems.

Contribution

The paper presents a novel product-convolution approximation scheme that is adaptive, scalable, and effective for high-rank operators arising in PDEs and inverse problems, with boundary artifact elimination.

Findings

Outperforms existing approximation methods in numerical tests.

Enables efficient application and inversion of operators using FFT and H-matrix techniques.

Provides automated construction of preconditioners for complex PDE operators.

Abstract

We present an adaptive grid matrix-free operator approximation scheme based on a "product-convolution" interpolation of convolution operators. This scheme is appropriate for operators that are locally translation-invariant, even if these operators are high-rank or full-rank. Such operators arise in Schur complement methods for solving partial differential equations (PDEs), as Hessians in PDE-constrained optimization and inverse problems, as integral operators, as covariance operators, and as Dirichlet-to-Neumann maps. Constructing the approximation requires computing the impulse responses of the operator to point sources centered on nodes in an adaptively refined grid of sample points. A randomized a-posteriori error estimator drives the adaptivity. Once constructed, the approximation can be efficiently applied to vectors using the fast Fourier transform. The approximation can be efficiently converted to hierarchical matrix ($H$-matrix) format, then inverted or factorized using scalable $H$-matrix arithmetic. The quality of the approximation degrades gracefully as fewer sample points are used, allowing cheap lower quality approximations to be used as preconditioners. This yields an automated method to construct preconditioners for locally translation-invariant Schur complements. We directly address issues related to boundaries and prove that our scheme eliminates boundary artifacts. We test the scheme on a spatially varying blurring kernel, on the non-local component of an interface Schur complement for the Poisson operator, and on the data misfit Hessian for an advection dominated advection-diffusion inverse problem. Numerical results show that the scheme outperforms existing methods.