Operator Inference

This page gives a short explanation of the mathematical details behind Operator Inference. For a full treatment, see [1]. Note that some notation has been altered for coding convenience and clarity.

Contents

• Problem Statement
• Projection-based Model Reduction
• Operator Inference via Least Squares
• Extensions and Variations

Problem Statement

Consider the (possibly nonlinear) system of n ordinary differential equations with state variable x, input (control) variable u, and independent variable t:

where

This system is called the full-order model (FOM). If n is large, as it often is in high-consequence engineering applications, it is computationally expensive to numerically solve the FOM. This package provides tools for constructing a reduced-order model (ROM) that is up to quadratic in the state x with optional linear control inputs u. The procedure is data-driven, non-intrusive, and relatively inexpensive. In the most general case, the code can construct and solve a reduced-order system with the polynomial form

where now

This reduced low-dimensional system approximates the original high-dimensional system, but it is much easier (faster) to solve because of its low dimension r << n.

Projection-based Model Reduction

Model reduction via projection occurs in three steps:

1. Data Collection: Gather snapshot data, i.e., solutions to the full-order model (the FOM) at various times / parameters.
2. Compression: Compute a low-rank basis (which defines a low-dimensional linear subspace) that captures most of the behavior of the snapshots.
3. Projection: Use the low-rank basis to construct a low-dimensional ODE (the ROM) that approximates the FOM.

This package focuses mostly on step 3 and provides a few light tools for step 2.

Let X be the n x k matrix whose k columns are each solutions to the FOM of length n (step 1), and let V_r be an orthonormal n x r matrix representation for an r-dimensional subspace (step 2). A common choice for V_r is the POD basis of rank r, the matrix whose columns are the first r singular vectors of X. We call X the snapshot matrix and V_r the basis matrix.

The classical intrusive approach to the projection step is to make the Ansatz

Inserting this into the FOM and multiplying both sides by the transpose of V_r (Galerkin projection) yields

This new system is r-dimensional in the sense that

If the FOM operator f is known and has a nice structure, this reduced system can be solved cheaply by precomputing any involved matrices and then applying a time-stepping scheme. For example, if f is linear in x and there is no input u, then

However, this approach breaks down if the FOM operator f is unknown, uncertain, or highly nonlinear.

Operator Inference via Least Squares

Instead of directly computing the reduced operators, the Operator Inference framework takes a data-driven approach: assuming a specific structure of the ROM (linear, quadratic, etc.), solve for the involved operators that best fit the data. For example, suppose that we seek a ROM of the form

We start with k snapshots x_j and inputs u_j. That is, x_j is an approximate solution to the FOM at time t_j with input u_j = u(t_j). We compute the basis matrix V_r from the snapshots (e.g., by taking the SVD of the matrix whose columns are the x_j) and project the snapshots onto the r-dimensional subspace defined by the basis via

We also require time derivative information for the snapshots. These may be provided by the FOM solver or estimated, for example with finite differences of the projected snapshots. With projected snapshots, inputs, and time derivative information in hand, we then solve the least-squares problem

Note that this minimum-residual problem is not (yet) in a typical linear least-squares form, as the unknown quantities are the matrices, not the vectors. Recalling that the vector 2-norm is related to the matrix Frobenius norm, i.e.,

we can rewrite the residual objective function in the more typical matrix form:

where

and where 1_k is a k-vector of 1's and d(r,m) = 1 + r + r² + m. For our purposes, the ⊗ operator between matrices denotes a column-wise Kronecker product, sometimes called the Khatri-Rao product.

It can be shown [1] that, under some idealized assumptions, the operators inferred by solving this data-driven minimization problem converge to the operators computed by explicit projection. The key idea, however, is that the inferred operators can be cheaply computed without knowing the full-order model. This is very convenient in "glass box" situations where the FOM is given by a legacy code for complex simulations but the target dynamics are known.

Regularization

The minimization problem given above decouples into r independent ordinary least-squares problems, one for each of the columns of O^T. Though each of the independent sub-problems are typically well-posed, they are also susceptible to noise from model misspecification, the truncation of the basis, numerical estimation of time derivatives, etc., and can therefore suffer from overfitting. To combat this, the problems can be regularized with a Tikhonov penalization. In this case, the complete minimization problem is given by

where o_i is row i of O and r_i is row i of R. Selecting an appropriate regularizer 𝚪 is somewhat problem dependent; see [10] for a principled approach.

Implementation Note: The Kronecker Product

The Kronecker product ⊗ introduces some redundancies. For example, x ⊗ x contains both x₁x₂ and x₂x₁. To avoid these redundancies, we introduce a "compact" Kronecker product which only computes the unique terms of the usual Kronecker product:

The dimension r(r+1)/2 arises because we choose 2 of r entries without replacement, i.e., this is a multiset coefficient:

We similarly define a cubic compact product recursively with the quadratic compact product.

Extensions and Variations

The Discrete Setting

The framework described above can also be used to construct reduced-order models for approximating discrete dynamical systems, as may arise from discretizing PDEs in both space and time. For instance, we can learn a discrete ROM of the form

The procedure is the same as described in the previous section with the exception that the snapshot matrix X has columns x_j for j = 0,1,...,k-1, while the right-hand side matrix R has columns x_j for j = 1,2,...,k. Thus, the (not yet regularized) least-squares problem to be solved has the form