Block-Quadrature Reduced-Order Modeling

• Block-quadrature reduced-order modeling (BQ-ROM) is defined as an approach that creates efficient, low-dimensional surrogates for high-fidelity models by decomposing the system into linear and quadratic blocks.
• It employs block-structured quadrature rules alongside projection techniques, such as Loewner interpolation and least-squares fitting, to ensure moment-matching and provide rigorous error bounds.
• BQ-ROM methods are practically applied in quadratic control, large-scale PDEs, and hyperreduction for fluid dynamics, leading to significant computational speed-up and accuracy.

Block-quadrature reduced-order modeling (BQ-ROM) refers to a class of techniques that produce efficient, low-dimensional surrogates of high-fidelity dynamical or PDE models by leveraging block-structured quadrature rules in conjunction with subspace projection, system identification, or data-driven fitting procedures. These methods systematically address the computational bottlenecks associated with large-scale, nonlinear, or parameter-dependent models by decomposing the problem into linear and non-linear (typically quadratic) "blocks," each handled with tailored quadrature or projection strategies. BQ-ROMs are notable for their broad applicability, including quadratic control systems, large-scale electromagnetic PDEs, and hyperreduced convection-dominated flows, and for their ability to guarantee interpolation or moment-matching constraints alongside principled error bounds.

1. Foundational Concepts and Mathematical Formulation

The foundational concept in BQ-ROMs is to approximate the high-dimensional system's behavior using a reduced subspace, while capturing both linear dynamics (via block-projected operators or transfer functions) and nonlinear effects such as harmonics or quadratic terms, using block least-squares fits or moment-matching over carefully chosen quadrature nodes.

A canonical example is the quadratic control system: $$\dot{x}(t) = A x(t) + H[x(t) \otimes x(t)] + B u(t),\quad y(t) = C x(t),$$ where $x(t)\in\mathbb{R}^n$, $A\in\mathbb{R}^{n\times n}$, $H$ is a symmetric third-order tensor ($\mathbb{R}^{n\times n^2}$), and $u$ and $y$ are input-output variables (Gosea et al., 2020). For PDE-based models, the block-structured finite-volume or finite-element discretizations lead to equations such as: $(A + \imath\,\omega I) x = b,$ where $A$ arises from self-adjoint mass/scattering operators and $B$ encodes source/receiver structure, with transfer functions written as block quadratic forms (Zimmerling et al., 22 Nov 2025).

Key unifying themes:

• Projection of operators and data on block-Krylov or balancing subspaces.
• Quadrature (Gaussian, Gauss–Radau) used for transfer function or Gramian approximation.
• Block-wise decomposition: a (typically) linear "block" extracted via Loewner, Krylov, or quadrature; and a nonlinear (quadratic or parameter-dependent) block learned via least squares or auxiliary harmonic information.
• Iterative refinement—particularly for nonlinear/quadratic blocks—possibly involving higher-order harmonics.

2. Block-quadrature ROM Workflow in Quadratic Control and Harmonic Fitting

Block-Quadrature ROMs for nonlinear (especially quadratic) control systems proceed in modular steps (Gosea et al., 2020):

1. Linear block extraction: Linear system identification is performed from input-output measurements or transfer function samples $G^{(1)}(j\omega_k)$ using Loewner interpolation. This provides reduced operators $(\tilde{A}, \tilde{B}, \tilde{C})$ that interpolate the first-harmonic data exactly.
2. Quadratic block fitting: The residual quadratic effect appears as higher harmonics (e.g., second, third, ...). The quadratic tensor $H_r$ is learned via convex least-squares fitting to measured or simulated higher-harmonic outputs using the reduced (Zauner or Loewner) surrogate.
3. Iterative block-quadrature coupling: If higher-harmonic (third, etc.) data are included, a fixed-point or alternating minimization is used, alternating between least-squares for the (linear) second-harmonic block and a linearized update for the third-harmonic response. This increases data efficiency and provides robustness to noise.
4. Synthesis: At convergence, the full reduced-order model generically takes the form:

$$\dot{x}_r = \tilde{A} x_r + H_r (x_r \otimes x_r) + \tilde{B} u,\quad y_r = \tilde{C} x_r$$

with $r \ll n$.

This is strictly non-intrusive in that it requires only input-output data, and permits convex regularization or stabilization for the quadratic block.

3. Block-quadrature-based ROMs for Large-scale PDEs and Block Quadrature Rules

When applied to large-scale parameterized PDEs, such as frequency-domain Maxwell's equations, BQ-ROM leverages block-Krylov subspaces and associated block Gaussian quadrature rules to compress the transfer function or solution operators efficiently (Zimmerling et al., 22 Nov 2025). This entails:

• Block Krylov subspace construction: Forming $\mathcal{K}_m(A,B)$, the span of $\{B, AB, ..., A^{m-1}B\}$, where $B$ aggregates multiple sources/receivers.
• Block Lanczos recursion: Generating an orthonormal basis with block-tridiagonal projected matrices $T_m$, enabling dimensionality reduction by projecting the dynamics/state evolution onto this subspace.
• Block Gaussian/block Gauss–Radau quadrature: Projecting the transfer function and computing frequency responses or state solutions in the reduced space; error control and stopping is based on convergence of block-Gauss to block-Gauss–Radau approximations.
• Adjoint calculations: The adjoint-state sensitivity needed for inverse problems is efficiently handled within the block-reduced space.
• Error estimation and efficiency: The residual transfer function difference under Gauss–Radau modification provides a rigorous bound for the reduction error, and computational/flux savings scale superlinearly with the reduction ratio.

This approach enables direct application to electromagnetic inverse problems, such as real-time borehole resistivity inversion and geosteering, and is robust to strong anisotropy, high contrast coefficients, and nonconformal grids, as corroborated by large-scale 3D simulations (Zimmerling et al., 22 Nov 2025).

4. Block-quadrature and Balanced Truncation in Second-order Systems

Block-quadrature methods generalize rational quadrature-based balanced truncation to preserve physical second-order structure in mechanical or electrical networks (Reiter et al., 11 Jun 2025). The workflow is:

• Quadrature-based Gramian approximation: Replace the Lyapunov-contour Gramian integrals with weighted sums using quadrature nodes ($\sigma_i$) and weights ($\omega_i$).
• Low-rank factor construction: Assemble square-root factors of controllability and observability Gramians via the quadrature solves.
• Balancing transformation and projection: Compute the SVD of the cross Gramian, apply balancing, and perform Petrov–Galerkin projection to obtain a reduced-order first-order surrogate.
• Second-order structure recovery: Partition balancing transforms appropriately to recover surrogate mass, damping, stiffness matrices, with proportional damping parameters $(\alpha, \beta)$ determined by least-squares fitting.
• Block hyperreduction: Reduction and optimal weighting are performed in block structured fashion to preserve block-wise system structure and guarantee moment-matching or error control.

This approach enables efficient, structure-preserving surrogates, with documented high accuracy (orders-of-magnitude singular value decay, sub-millisecond online solves), and is applicable to very large $n$ (Reiter et al., 11 Jun 2025).

5. Block-quadrature-driven Hyperreduction with Empirical Quadrature

For classes of nonlinear or convection-dominated PDEs, block-quadrature methodology can be combined with empirical quadrature and manifold learning to overcome the slow decay of Kolmogorov $n$-widths (Mirhoseini et al., 2023):

• Nonlinear trial manifold: State $u$ and domain mapping $x$ are parameterized by low-dimensional affine subspaces, yielding a nonlinear approximation manifold.
• Empirical quadrature hyperreduction: Residual and mesh-distortion objectives are further approximated by sparse weighted sums, with weights $\rho_e$ determined by $\ell_1$-minimization LP to enforce residual accuracy over training samples.
• Greedy offline basis and quadrature training: Offline construction iteratively augments trial subspaces and quadrature weights on representative parameter samples.
• Online phase: The ROM solves mesh-independent, low-dimensional optimization problems with computation localized to the smallest subset of elements with nonzero weights.
• Overcoming slow $n$-width decay: This approach decouples shock/transport features from the linear subspace, allowing the reduced model to retain accuracy using minimal modes and a small fraction of mesh elements (Mirhoseini et al., 2023).

Empirical results demonstrate two–three orders-of-magnitude speed-up and robust accuracy on shock-dominated fluid problems.

6. Comparative Summary and Properties

Aspect	Quadratic Control BQ-ROM (Gosea et al., 2020)	PDE Block Quadrature (Zimmerling et al., 22 Nov 2025)	Balanced2nd Block Quadrature (Reiter et al., 11 Jun 2025)	EQ Hyperreduction (Mirhoseini et al., 2023)
Linear block extraction	Loewner interpolation	Block Lanczos, Krylov projection	Quadrature balancing	N/A
Nonlinear/Quadratic block fitting	Least squares on higher harmonics	Not focal (applies to linear Maxwell)	Not used	Nonlinear manifold, residual minimization
Quadrature role	Harmonic moment sampling	Matrix Gaussian/Gauss-Radau integration	Rational quadrature for Gramian integrals	Empirical quadrature for element selection
Error/stopping criteria	Harmonic LS residual	Gauss–Radau transfer function difference	Singular value decay, Gramian error	Residual norm, greedy convergence
Data required	Transfer function and higher harmonics	PDE operator and source terms	System matrices or data-driven surrogates	Solution snapshots, Jacobian evaluations
Adjoint/inversion	Not explicit	Krylov-based adjoint evaluation	N/A	N/A

A plausible implication is that the block-quadrature paradigm flexibly adapts to both intrusive (PDE-discretization) and non-intrusive (data-driven, control) settings, maintaining strong interpolation and structure-preservation guarantees.

7. Applications, Robustness, and Computational Characteristics

Block-quadrature ROMs are deployed in multiple domains:

• Quadratic control model discovery: Non-intrusive system identification from oscillatory input-output measurements, yielding reduced quadratic control systems for embedded control (Gosea et al., 2020).
• Electromagnetic and geophysical inversion: Fast, accurate simulation and adjoint sensitivity computation for large-scale, anisotropic, nonconformal 3D Maxwell PDEs, facilitating real-time inversion and field operations (Zimmerling et al., 22 Nov 2025).
• Mechanical and structural modeling: Data-driven surrogates for structural dynamics with proportional damping constraints, significantly accelerating design iteration and optimization (Reiter et al., 11 Jun 2025).
• Hyperreduced model evaluation for fluid and transport phenomena: Effective ROMs for problems with slow $n$-width decay due to transport or shocks, substantially lowering computational cost and requirements (Mirhoseini et al., 2023).

Numerical results confirm superlinear scaling in run time and memory with increasing reduction ratios, monotonic error control, robustness to noise and parameter variation, and preservation of reciprocity and other physical invariants. The integration with effective-medium averaging permits reliable performance on heterogeneous and anisotropic media (Zimmerling et al., 22 Nov 2025).

Empirical quadrature as a block-quadrature extension achieves mesh-independent online complexity, with offline cost dominated by a modest number of high-dimensional solves and convex optimization problems. This suggests strong suitability for deployment in real-time, embedded, or resource-constrained environments.