SubApSnap: Solving parameter-dependent linear systems with a snapshot and subsampling

TL;DR

SubApSnap addresses the challenge of solving many related linear systems $A(p)x(p)=b(p)$ by combining snapshot-based subspace projection with row subsampling. Offline, it builds a low-dimensional span $X$ from solutions at a set of snapshot parameters, and online it solves a small, subsampled least-squares problem to obtain $x(p_*)$ without reading the full matrix $A(p_*)$. The method, which generalises DEIM to parameter-dependent $A(p)$, achieves substantial speedups (often orders of magnitude) while maintaining accuracy, demonstrated across parameter-dependent PDEs, transfer-function evaluation in model order reduction, large-scale delay systems, and kernel ridge regression hyperparameter searches. Theoretical bounds connect subsampling quality to the full problem through polar decompositions and subspace embeddings, while practical guidance on snapshot placement and sampling strategies supports robust performance. Overall, SubApSnap offers a scalable, data-efficient framework for rapid parameter-sweep computations in scientific computing and data science contexts.

Abstract

A growing number of problems in computational mathematics can be reduced to the solution of many linear systems that are related, often depending smoothly or slowly on a parameter $p$, that is, $A(p)x(p)=b(p)$. We introduce an efficient algorithm for solving such parameter-dependent linear systems for many values of $p$. The algorithm, which we call SubApSnap (for \emph{Sub}sampled $A(p)$ times \emph{Snap}shot), is based on combining ideas from model order reduction and randomised linear algebra: namely, taking a snapshot matrix, and solving the resulting tall-skinny least-squares problems using a subsampling-based dimension-reduction approach. We show that SubApSnap is a strict generalisation of the popular DEIM algorithm in nonlinear model order reduction. SubApSnap is a sublinear-time algorithm, and once the snapshot and subsampling are determined, it solves $A(p_*)x(p_*)=b(p_*)$ for a new value of $p_*$ at a dramatically improved speed: it does not even need to read the whole matrix $A(p_*)$ to solve the linear system for a new value of $p_*$. We prove under natural assumptions that, given a good subsampling and snapshot, SubApSnap yields solutions with small residual for all parameter values of interest. We illustrate the efficiency and performance of the algorithm with problems arising in PDEs, model reduction, and kernel ridge regression, where SubApSnap achieves speedups of many orders of magnitude over a standard solution; for example over $20,000\times$ for a $10^7\times 10^7$ problem, while providing good accuracy.

Figures (7)

• Figure 1: SubApSnap residual and its bounds \ref{['eqn:boundclosest']} (shown as bound 1) and \ref{['eq:bound1']} (bound 2), for the example \ref{['ex:tridiag']}. The plots are 0 at the snapshot points $p_i$ by construction.
• Figure 2: Left: Comparison of subsampling methods for a model order reduction problem; details are given in Section \ref{['sec:mor']}. Right: same plot only showing leverage score sampling, together with its estimate and heuristic confidence interval.
• Figure 3: Solution and error for 2-d heat equation with SubApSnap for two values of $p$.
• Figure 4: Transfer function $H(p)$ (above) and errors of their approximations with SubApSnap with four subsampling strategies (Random, LU with partial pivoting, ARP, and leverage-score sampling), and ApSnap (below).
• Figure 5: Model order reduction problem. Left: runtime for different values of $p$. Note that snapshots are taken only at $r=30$ points $p_i$, whereas other methods are used for $500$ values of $p$. Right: singular values of the snapshot matrix $X$ (before QR or SVD preprocessing).

Theorems & Definitions (3)

• Lemma 5.1
• Theorem 5.2
• Corollary 5.3