Discrete Empirical Interpolation Method

DEIM [1]

A nonlinear function $f(t) \in \mathbb{R}^n$ can be approximated in a finite space spanned by a set of basis functions, $V^f = [u_1, \dots, u_m] \in \mathbb{R}^{n \times m}$ ( $n \gg m$ ). The approximation of $f(t)$ is of the form

Copy LaTeX

We can use the orthogonality to determine $\alpha^f(t)$ . Assume that the residual $f(t) - V^f \alpha^f(t)$ are orthogonal to $V^f$ .
$Copy LaTeX$
We have $\alpha^f(t) = ((V^f)^T V^f)^{-1} (V^f)^T f(t)$ and
$Copy LaTeX$
$f_*(t)$ is the orthogonal projection of $f(t)$ onto $\text{span}(V^f)$ .
We can also use the interpolation property to determine $\alpha^f(t)$ . Assume that some components of $f(t)$ and $V^f \alpha^f(t)$ are strictly equal.
$Copy LaTeX$
where
$Copy LaTeX$
Suppose that $(P^f)^T V^f \in \mathbb{R}^{m \times m}$ is non-singular, we have $\alpha^f(t) = ((P^f)^T V^f)^{-1} (P^f)^T f(t)$ and
$Copy LaTeX$
$\hat{f}(t)$ is an oblique projection of $f(t)$ onto $\text{span}(V^f)$ .

While the orthogonal projection $f_*(t)$ provides the optimal approximation with respect to norms,

Copy LaTeX

the oblique projection $\hat{f}(t)$ constructed using interpolation properties can achieve approximation and reduce computational cost: only $m$ components of $f(t) \in \mathbb{R}^n$ need to be calculated. ( $m \ll n$ )

Copy LaTeX

To obtain the approximation $\hat{f}(t)$ , we must specify the interpolation indices $\{\varphi_i\}_{i=1}^m$ . The DEIM algorithm selects an index corresponding to one of these discrete spatial points at each iteration to limit growth of an error bound. The main steps of DEIM are described in the DEIM algorithm.

DEIM algorithm

Input: linearly independent vectors $u_1,\ldots,u_m\in\mathbb R^n$ .

Choose $\varphi_1\in\operatorname*{arg\,max}_i |(u_1)_i|$ , and set $V_1=[u_1]$ , $P_1=[e_{\varphi_1}]$ .
For :
- Solve $(P_{\ell-1}^T V_{\ell-1})\alpha=P_{\ell-1}^T u_\ell$ for $\alpha\in\mathbb R^{\ell-1}$ .
- Compute $r_\ell=u_\ell-V_{\ell-1}\alpha$ .
- Choose $\varphi_\ell\in\operatorname*{arg\,max}_i |(r_\ell)_i|$ .
- Set $V_\ell=[V_{\ell-1},u_\ell]$ and $P_\ell=[P_{\ell-1},e_{\varphi_\ell}]$ .
Return the indices $\varphi_1,\ldots,\varphi_m$ .

We need to explain some details of the DEIM algorithm to ensure that the DEIM algorithm is well-defined. Since $\{u_\ell\}_{\ell=1}^m \subset \mathbb{R}^n$ are linearly independent, we have $u_\ell \not\in \text{span}(V_{\ell-1})$ ,

Copy LaTeX

Using interpolation properties, it can be proven that the components of $r_\ell$ corresponding to the selected indices are strictly zero,

Copy LaTeX

Thus, $\varphi_\ell$ exists and $\varphi_\ell \notin \{\varphi_1, \dots, \varphi_{\ell-1}\}$ . We have $P_{1}^T V_{1} = u_{1,\varphi_1} \neq 0$ because $\|u_1\|_\infty = |u_{1,\varphi_1}| \neq 0$ ; If $P_{\ell-1}^T V_{\ell-1}$ is non-singular, we have

Copy LaTeX

where $e_{\varphi_\ell}^T r_\ell = r_{\ell,\varphi_\ell} \neq 0$ because $\|r_\ell\|_{\infty} = |r_{\ell,\varphi_\ell}| \neq 0$ . Thus $P_{\ell}^T V_{\ell}$ is non-singular for $\ell = 1, \dots, m$ .

Lemma

Let $f \in \mathbb{R}^n$ be an arbitrary vector. Let $\{u_\ell\}_{\ell=1}^m \subset \mathbb{R}^n$ be a given orthonormal set of vectors. The DEIM approximation of order $m \le n$ for $f$ in the space spanned by $\{u_\ell\}_{\ell=1}^m$ is

Copy LaTeX

where $V = [u_1, \dots, u_m]$ and $P = [e_{\varphi_1}, \dots, e_{\varphi_m}]$ with $\{\varphi_1, \dots, \varphi_m\}$ being the output from the DEIM algorithm with the input basis $\{u_\ell\}_{\ell=1}^m$ . An error bound for $\hat{f}$ is given by

Copy LaTeX

where $f_* = V V^T f$ is the orthogonal projection of $f$ . Vector norms are Euclidean norms and matrix norms are the induced spectral norms.

Proof

Define the projector $\mathbb{P} = V (P^T V)^{-1} P^T$ . We have

Copy LaTeX

The bound for the error $\|f - \hat{f}\|$ is given by

Copy LaTeX

For $m=n$ , the interpolation is exact and the error bound is immediate. For $1\le m<n$ , note that

Copy LaTeX

The first equality uses $\|\mathbb{I}-\mathbb{P}\| = \|\mathbb{P}\|$ for a projector $\mathbb{P}$ that is neither zero nor the identity, in the spectral norm.

Remark

$\|f - f_*\|$ relates to the approximation quality of the space spanned by $V$ . $\|(P^T V)^{-1}\|$ quantifies the effect of the sampling points and can be bound by

Copy LaTeX

However, the bound is not useful as an a priori estimate, since it is very pessimistic and grows far more rapidly than the actual observed values of $\| (P^T V)^{-1} \|$ .

To enhance robustness, more interpolation points than the basis functions, i.e. $\widetilde{m} > m$ , can be used.

The first $m$ measurement points are found using the classic DEIM procedure.
The next $\widetilde{m}-m$ points are generated from a uniform random distribution where care is taken to not select repeating indices.

This procedure strikes a good balance between efficiency and accuracy and is simple to implement. In this case, the approximation is given by

Copy LaTeX

where $((P^f)^T V^f)^{\dagger}$ is the pseudo-inverse.

DEIM with oversampling

Input: the same basis and an integer $m<\widetilde m\le n$ .

Run the DEIM algorithm to select the first $m$ indices.
Uniformly sample $\widetilde m-m$ further indices without replacement from $\{1,\ldots,n\}\setminus\{\varphi_1,\ldots,\varphi_m\}$ .
Form $P=[e_{\varphi_1},\ldots,e_{\varphi_{\widetilde m}}]$ and use $(P^T V)^{\dagger}$ in the approximation.

POD + DEIM

Consider a system of nonlinear ODEs of the form

Copy LaTeX

with appropriate initial conditions. Here $t \in [0,T]$ denotes time, $\boldsymbol{y}(t) = [y_1(t),\dots,y_n(t)]^T \in \mathbb{R}^n$ , $\boldsymbol{A} \in \mathbb{R}^{n \times n}$ is a constant matrix, and $\boldsymbol{F}$ is a nonlinear function evaluated at $\boldsymbol{y}(t)$ component-wise.

Copy LaTeX

The dimension $n$ reflects the number of spatial grid points used in the FD discretization.

Given a set of snapshots $\{\boldsymbol{y}_1, \dots , \boldsymbol{y}_{n_s} \} \subset \mathbb{R}^n$ , let $\boldsymbol{Y} = [\boldsymbol{y}_1, \dots, \boldsymbol{y}_{n_s}] \in \mathbb{R}^{n \times n_s}$ , $\mathcal{Y} = \text{span}(Y) \subset \mathbb{R}^n$ . A POD basis of dimension $k < r = \operatorname{rank}(Y)$ is a set of orthonormal vectors $\{\phi_i\}_{i=1}^k \subset \mathbb{R}^n$ whose linear span best approximates $\mathcal{Y}$ . Suppose the SVD of $Y$ is

Copy LaTeX

Then the POD basis is $\{\boldsymbol{v}_i\}_{i=1}^k$ . Let $\boldsymbol{V}_k = [\boldsymbol{v}_1,\dots,\boldsymbol{v}_k] \in \mathbb{R}^{n \times k}$ . Then, by replacing $\boldsymbol{y}(t)$ in Eq. $\eqref{eq:ode}$ by $\boldsymbol{V}_k \boldsymbol{\widetilde{y}}(t)$ , ( $\boldsymbol{\widetilde{y}}(t) \in \mathbb{R}^k$ ), and projecting the system onto $\boldsymbol{V}_k$ , the reduced system becomes

Copy LaTeX

where $\widetilde{\boldsymbol{A}} = \boldsymbol{V}_k^T \boldsymbol{A} \boldsymbol{V}_k \in \mathbb{R}^{k \times k}$ , and the nonlinear term

Copy LaTeX

has a computational complexity that depends on $n$ .

These nonlinear snapshots are the sets $\{\boldsymbol{F}(\boldsymbol{y}_1), \dots , \boldsymbol{F}(\boldsymbol{y}_{n_s})\} \subset \mathbb{R}^n$ . Let $\boldsymbol{U}_m = [\boldsymbol{u}_1,\dots,\boldsymbol{u}_m] \in \mathbb{R}^{n \times m}$ where $\{\boldsymbol{u}_j\}$ is the POD basis of nonlinear snapshots. Define $\boldsymbol{P} = [\boldsymbol{e}_{\varphi_1}, \dots, \boldsymbol{e}_{\varphi_m}]$ with $\{\varphi_1, \dots, \varphi_m\}$ being the output from the DEIM algorithm with the input basis $\{\boldsymbol{u}_\ell\}_{\ell=1}^m$ , then the nonlinear function approximated by DEIM can be written as

Copy LaTeX

The nonlinear term can thus be approximated by

Copy LaTeX

The complexity for computing this approximation of the nonlinear term is independent of dimension n of the full-order system.

References

Chaturantabut S, Sorensen DC. Nonlinear model reduction via discrete empirical interpolation. SIAM Journal on Scientific Computing. 2010;32(5):2737–64. Google Scholar Back