Nested Operator Inference (OpInf)

• Nested Operator Inference is a reduced-order modeling framework that organizes operator learning hierarchically to optimize accuracy and conditioning.
• It employs warm-start strategies by using previously learned dominant modes as initial guesses, focusing regularization on newly introduced dynamics.
• The approach is applied both in standalone hierarchical formulations and embedded within Schwarz methods, enabling efficient coupling of high-fidelity and reduced models.

Searching arXiv for the cited Nested OpInf papers to ground the article in the current literature. Nested Operator Inference (OpInf) denotes a family of reduced-order modeling procedures in which operator inference is organized through an explicit hierarchy rather than executed as a single monolithic regression. In the hierarchical formulation, the method exploits the inherent hierarchy within the reduced space to iteratively construct initial guesses for the OpInf learning problem that prioritize the interactions of the dominant modes, and the initial guess computed for any target reduced dimension corresponds to a ROM with provably smaller or equal snapshot reconstruction error than with standard OpInf (Aretz et al., 15 Aug 2025). A related formulation embeds pre-trained subdomain-local OpInf ROMs inside an overlapping Schwarz alternating method, so that reduced models and high-fidelity models are alternated within each domain-decomposition iteration rather than trained or advanced in isolation (Moore et al., 6 Oct 2025). In these 2025 formulations, “nested” therefore refers either to progression along a POD hierarchy or to embedding OpInf advances inside an outer Schwarz loop.

1. Terminological scope and problem setting

The term has two distinct uses in the cited literature. In "Nested Operator Inference for Adaptive Data-Driven Learning of Reduced-order Models" (Aretz et al., 15 Aug 2025), nestedness is internal to the reduced space: one learns ROM operators block by block, starting with the smallest subspace $V_1$, then $V_2,\dots,V_r$, and uses operators from $V_{s-1}$ as an initial guess for $V_s$. In "Domain Decomposition-Based Coupling of High-Fidelity Finite Element and Reduced Order Operator Inference Models Using the Schwarz Alternating Method" (Moore et al., 6 Oct 2025), nestedness is algorithmic: high-fidelity FE solves and low-dimensional OpInf ROM advances are interleaved within each Schwarz iteration. The later solid-dynamics work "Hybrid coupling with operator inference and the overlapping Schwarz alternating method" extends this hybrid viewpoint to FOM–ROM and ROM–ROM couplings on overlapping subdomains (Tezaur et al., 20 Nov 2025).

This distinction matters because the two formulations address different bottlenecks. The hierarchical formulation is motivated by ill-conditioning in large-$r$ OpInf regression systems and by the hierarchy of POD modes. The Schwarz-embedded formulation is motivated by multiscale coupling, geometry and mesh flexibility, and the localization of expensive high-fidelity solves to difficult regions. A common ambiguity is therefore terminological rather than mathematical: nested OpInf is not a single regression formula, but a broader organizational principle for how OpInf is trained or deployed.

2. Operator Inference background and the source of nesting

The hierarchical paper formulates the full-order model from snapshot data $x(t_1),\dots,x(t_K)\in\mathbb R^n$ of a high-dimensional dynamical system

$$\dot x(t)=\mathbf c+\mathbf A x(t)+\mathbf H[x(t)\otimes x(t)],\qquad x(0)=x_0,$$

with $\mathbf c\in\mathbb R^n$, $\mathbf A\in\mathbb R^{n\times n}$, and $\mathbf H\in\mathbb R^{n\times n^2}$. A POD basis

$V_2,\dots,V_r$0

is chosen so that the projection error $V_2,\dots,V_r$1 is minimal among all $V_2,\dots,V_r$2-dimensional subspaces. The reduced state $V_2,\dots,V_r$3 then satisfies a ROM with operators $V_2,\dots,V_r$4, $V_2,\dots,V_r$5, and $V_2,\dots,V_r$6 obtained by projection or learned from data (Aretz et al., 15 Aug 2025).

Standard OpInf constructs a data matrix

$V_2,\dots,V_r$7

where $V_2,\dots,V_r$8, and solves a regularized least-squares problem for $V_2,\dots,V_r$9. Under mild conditions this recovers the Galerkin-projected operators. The central difficulty is that, as $V_{s-1}$0 grows, the smallest singular value of $V_{s-1}$1 decays, making OpInf unstable or requiring heavy regularization. The nested formulation is introduced precisely to exploit the fact that the POD projection error decays rapidly for small $V_{s-1}$2 and more slowly for larger $V_{s-1}$3, so that dominant-mode interactions can be learned first on better-conditioned subproblems (Aretz et al., 15 Aug 2025).

3. Hierarchical mathematical formulation

For $V_{s-1}$4, nested OpInf writes the operators at level $V_{s-1}$5 as zero-padded expansions of the operators at level $V_{s-1}$6 plus increments:

$V_{s-1}$7

and similarly for $V_{s-1}$8. The zero-padded quantities define the natural initial guess $V_{s-1}$9. Instead of regularizing directly toward zero, the nested objective regularizes the deviation from this inherited model:

$V_s$0

Equivalently, one may optimize the increment $V_s$1 toward matching $V_s$2 (Aretz et al., 15 Aug 2025).

This reformulation changes the role of regularization. In standard OpInf, all blocks are regularized uniformly toward the zero operator. In nested OpInf, previously learned dominant interactions are retained as the default state, while regularization is concentrated on the newly introduced directions. The paper explicitly characterizes this as a structured form of regularization. A plausible implication is that the method separates stabilization of sensitive high-index modes from preservation of already learned low-index dynamics.

4. Algorithmic realization, warm-starts, and guarantees

Algorithmically, nested OpInf accepts a target dimension $V_s$3, the POD basis $V_s$4, projected snapshots $V_s$5, derivative matrix $V_s$6, a weight grid $V_s$7, and optionally a starting dimension $V_s$8 with pre-existing initial guesses. For each $V_s$9, the algorithm restricts $r$0 and $r$1 to the first $r$2 columns, forms $r$3, constructs the base initial guess by zero-padding the accepted $r$4 model, evaluates the base ROM error $r$5, solves the nested least-squares problem for each weight triple, integrates the resulting ROM to compute the error $r$6 and conditioning $r$7, and selects the index

$r$8

that is, the best stability among those achieving near-optimal error (Aretz et al., 15 Aug 2025).

The same framework supports warm-starts. One can choose $r$9 and feed a previously learned ROM at $x(t_1),\dots,x(t_K)\in\mathbb R^n$0 as the initial guess for iteration $x(t_1),\dots,x(t_K)\in\mathbb R^n$1. The paper also states that the nested framework naturally supports on-the-fly expansion of $x(t_1),\dots,x(t_K)\in\mathbb R^n$2, model-form updates such as polynomial order, and greedy snapshot selection.

Theoretical properties are given in Proposition 3.1 and Corollary 3.2. Proposition 3.1 relates the reconstruction error of the padded $x(t_1),\dots,x(t_K)\in\mathbb R^n$3-dimensional model to that of the $x(t_1),\dots,x(t_K)\in\mathbb R^n$4-dimensional model through an additive term $x(t_1),\dots,x(t_K)\in\mathbb R^n$5, and the paper states that the error does not blow up uncontrollably and near $x(t_1),\dots,x(t_K)\in\mathbb R^n$6 it even improves. Corollary 3.2 gives a “do-no-harm guarantee”: if one were to pad zeros all the way to dimension $x(t_1),\dots,x(t_K)\in\mathbb R^n$7, the ROM error would be no worse than the trivial zero operator ROM, so the nested initial guess at $x(t_1),\dots,x(t_K)\in\mathbb R^n$8 is at least as good as Tikhonov toward zero, with strict improvement except in a worst-case scenario. A worst-case snapshot-error bound is also reported. On the computational side, each of the $x(t_1),\dots,x(t_K)\in\mathbb R^n$9 least-squares solves at level $$\dot x(t)=\mathbf c+\mathbf A x(t)+\mathbf H[x(t)\otimes x(t)],\qquad x(0)=x_0,$$0 scales as $$\dot x(t)=\mathbf c+\mathbf A x(t)+\mathbf H[x(t)\otimes x(t)],\qquad x(0)=x_0,$$1, each ROM integration as $$\dot x(t)=\mathbf c+\mathbf A x(t)+\mathbf H[x(t)\otimes x(t)],\qquad x(0)=x_0,$$2, and the total cost summed over $$\dot x(t)=\mathbf c+\mathbf A x(t)+\mathbf H[x(t)\otimes x(t)],\qquad x(0)=x_0,$$3 is $$\dot x(t)=\mathbf c+\mathbf A x(t)+\mathbf H[x(t)\otimes x(t)],\qquad x(0)=x_0,$$4. All offline steps use only $$\dot x(t)=\mathbf c+\mathbf A x(t)+\mathbf H[x(t)\otimes x(t)],\qquad x(0)=x_0,$$5, so the offline cost is independent of the full-order dimension $$\dot x(t)=\mathbf c+\mathbf A x(t)+\mathbf H[x(t)\otimes x(t)],\qquad x(0)=x_0,$$6 (Aretz et al., 15 Aug 2025).

5. Numerical behavior and reported performance

The hierarchical nested OpInf paper reports two main numerical studies: a cubic heat conduction problem and a large-scale parameterized model of the Greenland ice sheet. In the cubic heat conduction case, the PDE is $$\dot x(t)=\mathbf c+\mathbf A x(t)+\mathbf H[x(t)\otimes x(t)],\qquad x(0)=x_0,$$7 on $$\dot x(t)=\mathbf c+\mathbf A x(t)+\mathbf H[x(t)\otimes x(t)],\qquad x(0)=x_0,$$8 with Dirichlet boundary conditions, discretized by FE with $$\dot x(t)=\mathbf c+\mathbf A x(t)+\mathbf H[x(t)\otimes x(t)],\qquad x(0)=x_0,$$9 and Crank–Nicolson time stepping $\mathbf c\in\mathbb R^n$0. Training uses three $\mathbf c\in\mathbb R^n$1 values $\mathbf c\in\mathbb R^n$2, $\mathbf c\in\mathbb R^n$3 snapshots up to $\mathbf c\in\mathbb R^n$4, and a POD basis of size $\mathbf c\in\mathbb R^n$5 that captures $\mathbf c\in\mathbb R^n$6 energy. Nested OpInf uses $\mathbf c\in\mathbb R^n$7, iterative updates $\mathbf c\in\mathbb R^n$8, and an offline cost of $\mathbf c\in\mathbb R^n$9 s, while the standard OpInf grid search over $\mathbf A\in\mathbb R^{n\times n}$0 weight combinations costs $\mathbf A\in\mathbb R^{n\times n}$1 s (Aretz et al., 15 Aug 2025).

In that benchmark, the reported online speed-up is $\mathbf A\in\mathbb R^{n\times n}$2 ($\mathbf A\in\mathbb R^{n\times n}$3 ms versus $\mathbf A\in\mathbb R^{n\times n}$4 s). On the training interval $\mathbf A\in\mathbb R^{n\times n}$5, standard OpInf has mean and max reconstruction errors of $\mathbf A\in\mathbb R^{n\times n}$6 and $\mathbf A\in\mathbb R^{n\times n}$7, whereas nested OpInf has $\mathbf A\in\mathbb R^{n\times n}$8 and $\mathbf A\in\mathbb R^{n\times n}$9, described as approximately $\mathbf H\in\mathbb R^{n\times n^2}$0 smaller mean error. Extrapolated to $\mathbf H\in\mathbb R^{n\times n^2}$1, standard OpInf yields mean $\mathbf H\in\mathbb R^{n\times n^2}$2 and max $\mathbf H\in\mathbb R^{n\times n^2}$3, while nested OpInf yields mean $\mathbf H\in\mathbb R^{n\times n^2}$4 and max $\mathbf H\in\mathbb R^{n\times n^2}$5. The effectivity $\mathbf H\in\mathbb R^{n\times n^2}$6 is reported as up to $\mathbf H\in\mathbb R^{n\times n^2}$7, mean $\mathbf H\in\mathbb R^{n\times n^2}$8 for standard OpInf, and up to $\mathbf H\in\mathbb R^{n\times n^2}$9, mean $V_2,\dots,V_r$00 for nested OpInf. The paper further states that, for $V_2,\dots,V_r$01, standard OpInf stagnates at $V_2,\dots,V_r$02, while nested OpInf continues to improve (Aretz et al., 15 Aug 2025).

For the Greenland ice sheet model, the state dimension is $V_2,\dots,V_r$03 for ice thickness, plus $V_2,\dots,V_r$04 for velocity, with adaptive time stepping from 2015 to 2050. Training uses 17 parameter values and approximately 7000 snapshots of thickness change from 2020. The thickness POD basis with $V_2,\dots,V_r$05 captures $V_2,\dots,V_r$06 of snapshot energy, while the velocity basis $V_2,\dots,V_r$07 captures $V_2,\dots,V_r$08. The reduced model has the affine-parameterized form

$V_2,\dots,V_r$09

with $V_2,\dots,V_r$10 and $V_2,\dots,V_r$11. The paper reports learning $V_2,\dots,V_r$12 constant vectors $V_2,\dots,V_r$13 and $V_2,\dots,V_r$14 matrices $V_2,\dots,V_r$15, for approximately $V_2,\dots,V_r$16 degrees of freedom, from a data matrix $V_2,\dots,V_r$17 with $V_2,\dots,V_r$18 (Aretz et al., 15 Aug 2025).

The corresponding performance numbers are an offline cost of $V_2,\dots,V_r$19 min, online speed-up from approximately $V_2,\dots,V_r$20 in the worst case to approximately $V_2,\dots,V_r$21 on average, and time-averaged relative errors of $V_2,\dots,V_r$22 on training parameters and $V_2,\dots,V_r$23 on test parameters, with average $V_2,\dots,V_r$24. The paper states that effectivity remains below $V_2,\dots,V_r$25 for all times and parameters and that, despite model-form approximations such as stationary velocity and annual SMB, nested OpInf yields ROM accuracy close to POD projection limits. A plausible implication is that the hierarchical strategy is especially useful when the reduced dimension is large enough for monolithic OpInf to become ill-conditioned, but the dominant POD subspace remains strongly informative.

6. Schwarz-embedded nested OpInf in domain decomposition

A second line of work places OpInf inside an overlapping Schwarz alternating method rather than inside a hierarchy of reduced dimensions. For the unsteady convection–diffusion–reaction problem on $V_2,\dots,V_r$26,

$V_2,\dots,V_r$27

the FE semi-discretization on a subdomain $V_2,\dots,V_r$28 is

$V_2,\dots,V_r$29

while the reduced OpInf model is

$V_2,\dots,V_r$30

The reduced operators are inferred from snapshot data $V_2,\dots,V_r$31 by solving a Tikhonov-regularized least-squares problem. The online algorithm then nests two loops: a time-stepping loop and a Schwarz iteration loop. At each Schwarz iteration, each subdomain is solved sequentially as FE or ROM depending on its designation, the solution trace on the interface is extracted, and that trace becomes Dirichlet input for the neighboring subdomain at the next subsolve. No re-training occurs online (Moore et al., 6 Oct 2025).

The reported numerical example uses $V_2,\dots,V_r$32, $V_2,\dots,V_r$33, $V_2,\dots,V_r$34, $V_2,\dots,V_r$35, $V_2,\dots,V_r$36, $V_2,\dots,V_r$37, homogeneous $V_2,\dots,V_r$38, four overlapping square subdomains, and a uniform FE mesh with $V_2,\dots,V_r$39 in each subdomain. The upper-right block $V_2,\dots,V_r$40, which contains the boundary layer, is kept as FE, while $V_2,\dots,V_r$41 are ROM. ROM training uses the time window $V_2,\dots,V_r$42, snapshots every $V_2,\dots,V_r$43, $V_2,\dots,V_r$44, POD basis dimension $V_2,\dots,V_r$45, and unregularized inference $V_2,\dots,V_r$46. For comparison, the monolithic OpInf ROM uses $V_2,\dots,V_r$47 and $V_2,\dots,V_r$48. At $V_2,\dots,V_r$49, the reported CPU times and time-averaged relative pointwise errors versus a monolithic FE run are $V_2,\dots,V_r$50 s and $V_2,\dots,V_r$51 for all-FE Schwarz, $V_2,\dots,V_r$52 s and $V_2,\dots,V_r$53 for OpInf-FE Schwarz, and $V_2,\dots,V_r$54 s and $V_2,\dots,V_r$55 for monolithic OpInf. The paper states that the hybrid OpInf-FE method is approximately $V_2,\dots,V_r$56 faster than all-FE Schwarz, retains high accuracy without any regularization, and that the monolithic OpInf ROM is unstable or inaccurate even with large basis and tuning of $V_2,\dots,V_r$57 (Moore et al., 6 Oct 2025).

The 3D solid-dynamics extension formulates each subdomain-local FOM as

$V_2,\dots,V_r$58

augmented by Dirichlet data on $V_2,\dots,V_r$59, and postulates an up-to-cubic reduced model

$V_2,\dots,V_r$60

Snapshots are gathered from a FOM–FOM Schwarz-coupled run; subdomain and boundary POD bases are computed by thin SVD; and the reduced operators are learned through regularized least squares using finite-difference approximations of acceleration. The overlapping Schwarz algorithm then alternates subdomain solves on controller time steps, with spatial and temporal projection operators transferring interface data across nonmatching meshes and time grids (Tezaur et al., 20 Nov 2025).

For that solid-dynamics setting, the paper reports four test problems—clamped beam, bolted joint, torsion bar, and tension specimen—and defines the time-integrated relative error

$V_2,\dots,V_r$61

The reported observations are that FOM–FOM couplings typically require approximately $V_2,\dots,V_r$62 Schwarz iterations, ROM–ROM couplings often converge in $V_2,\dots,V_r$63 iterations, FOM–ROM couplings achieve $V_2,\dots,V_r$64 speedups, and ROM–ROM couplings attain $V_2,\dots,V_r$65 speedups at less than or equal to $V_2,\dots,V_r$66 global error. The paper attributes the reduction in Schwarz iterations to smoother subdomain solutions under OpInf ROM coupling. This suggests that, in the Schwarz setting, nestedness is not only a deployment strategy for local surrogates but also a mechanism that can empirically alter interface convergence behavior (Tezaur et al., 20 Nov 2025).

7. Limitations, generalizations, and open directions

The hierarchical nested OpInf paper states several limitations directly. It requires a nested basis hierarchy such as POD; heavily nonlinear systems may still require additional stabilization such as reprojection; and the method is most beneficial when the target reduced dimension is large enough to cause ill-conditioning in a monolithic OpInf solve, when a natural hierarchy of modes is available, and when one needs on-the-fly expansion or warm-start functionality. Reported extensions include block-wise regularization, adaptive weight grids, multi-fidelity snapshot selection, application to non-polynomial operators via EIM/DEIM, and the use of IterativeUpdates to correct early-iteration divergence by matching ROM roll-outs (Aretz et al., 15 Aug 2025).

The Schwarz-based studies report a different set of limitations. The FE–OpInf convection–diffusion–reaction study requires an offline FE training phase on each ROM subdomain, is presently linear, and notes that extension to non-linear PDEs will require inferring nonlinear operators. It also remarks that the Schwarz method is sequential in subdomains, and identifies parallel overlap Schwarz or optimized transmission conditions such as Robin and Neumann as possible accelerations. Listed future developments include non-overlapping Schwarz via Mortar methods, parallel asynchronous Schwarz, data-adaptive basis enrichment or online operator update for strong parameter variations, and extension to multi-physics problems such as Navier–Stokes and thermo-elasticity (Moore et al., 6 Oct 2025).

The solid-dynamics coupling paper likewise points to non-overlapping and optimized Schwarz, additive Schwarz for parallel scalability, online switching between OpInf ROM and high-fidelity FOM driven by error indicators, and alternative subdomain surrogates including symplectic or Hamiltonian-preserving neural surrogates, PGD, POD–Galerkin, and physics-informed networks. Taken together, these directions indicate that nested OpInf has become a label for two converging research programs: one centered on hierarchical regularization and warm-start learning in reduced coordinates, and the other centered on hybrid multi-fidelity deployment within domain decomposition. A plausible implication is that future work may combine both senses of nesting simultaneously, using hierarchical operator learning inside each subdomain ROM while embedding those ROMs inside optimized Schwarz couplings (Tezaur et al., 20 Nov 2025).