[Click here for a PDF of this post with nicer formatting]


Peeter’s lecture notes from class. These may be incoherent and rough.

Model order reduction (cont).

An approximation of the following system is sought

\BG \Bx(t) + C \dot{\Bx}(t) = B \Bu(t)
\By(t) = \BL^\T \Bx(t).

The strategy is to attempt to find a \( N \times q \) projector \( \BV \) of the form

\BV =
\Bv_1 & \Bv_2 & \cdots & \Bv_q

so that the solution of the constrained q-variable state vector \( \Bx_q \) is sought after letting

\Bx(t) = \BV \Bx_q(t).

Moment matching

&= \lr{ \BG + s \BC }^{-1} \BB \\
&= \BM_0 + \BM_1 s + \BM_2 s^2 + \cdots + \BM_{q-1} s^{q-1} + M_q s^q + \cdots

The reduced model is created such that

\BM_0 + \BM_1 s + \BM_2 s^2 + \cdots + \BM_{q-1} s^{q-1} + \tilde{\BM}_q s^q.

Form an \( N \times q \) projection matrix

\BV_q \equiv
\BM_0 & \BM_1 & \cdots & \BM_{q-1}

With the substitution of fig. 1, becomes


This is a system of \( N \) equations, in \( q \) unknowns. A set of moments from the frequency domain have been used to project the time domain system. This relies on the following unproved theorem (references to come)


If \( \text{span}\{ \Bv_q \} = \text{span} \{ \BM_0, \BM_1, \cdots, \BM_{q-1} \} \), then the reduced model will match the first \( q \) moments.

Left multiplication by \( \BV_q^\T \) yields fig. 2.


This is now a system of \( q \) equations in \( q \) unknowns.


\BG_q = \BV_q^\T \BG \BV_q
\BC_q = \BV_q^\T \BC \BV_q
\BB_q = \BV_q^\T \BB
\BL_q^\T = \BL^\T \BV_q

the system is reduced to

\BG_q \Bx_q(t) + \BC_q \dot{\Bx}_q(t) = \BB_q \Bu(t).
\By(t) \approx \BL_q^\T \Bx_q(t)

Moments calculation

\lr{ \BG + s \BC }^{-1} \BB = \BM_0 + \BM_1 s + \BM_2 s^2 + \cdots

\BB &=
\lr{ \BG + s \BC }
\BM_0 +
\lr{ \BG + s \BC }
\BM_1 s +
\lr{ \BG + s \BC }
\BM_2 s^2 + \cdots \\
\BG \BM_0
+ s \lr{ C \BM_0 + \BG \BM_1 }
+ s^2 \lr{ C \BM_1 + \BG \BM_2 }
+ \cdots

Since \( \BB \) is a zeroth order matrix, setting all the coefficients of \( s \) equal to zero provides a method to solve for the moments

\BB &= \BG \BM_0 \\
-\BC \BM_0 &= \BG \BM_1 \\
-\BC \BM_1 &= \BG \BM_2 \\

The moment \( \BM_0 \) can be found with LU of \( \BG \), plus the forward and backward substitutions. Proceeding recursively, using the already computed LU factorization, each subsequent moment calculation requires only one pair of forward and backward substitutions.

Numerically, each moment has the exact value

\BM_q = \lr{- \BG^{-1} \BC }^q \BM_0.

As \( q \rightarrow \infty \), this goes to some limit, say \( \Bw \). The value \( \Bw \) is related to the largest eigenvalue of \( -\BG^{-1} \BC \). Incidentally, this can be used to find the largest eigenvalue of \( -\BG^{-1} \BC \).

The largest eigenvalue of this matrix will dominate these factors, and can cause some numerical trouble. For this reason it is desirable to avoid such explicit moment determination, instead using implicit methods.

The key is to utilize the theorem above, and look instead for an alternate basis \( \{ \Bv_q \} \) that also spans \( \{ \BM_0, \BM_1, \cdots, \BM_q \} \).

Space generate by the moments


\BM_q = \BA^q \BR,

where in this case

\BA &= – \BG^{-1} \BC \\
\BR &= \BM_0 = -\BG \BB

The span of interest is

\text{span} \{ \BR, \BA \BR, \BA^2 \BR, \cdots, \BA^{q-1} \BR \}.

Such a sequence is called a Krylov subspace. One method to compute such a basis, the Arnoldi process, relies on Gram-Schmidt orthonormalization methods:


Some numerical examples and plots on the class slides.