Eigensystem Realization Algorithm (ERA)
As shown in the previous section, a structural system’s dynamic behavior
can be represented by the four coefficients
(\(\mathbf{A},\mathbf{B},\mathbf{C},\mathbf{D}\)) of its discrete
LTI state-space representation. The Ho-Kalman
Algorithm, or
Eigensystem Realization Algorithm,
produces a reduced order model for these four coefficients,
(\(\mathbf{\tilde{A}},\mathbf{\tilde{B}},\mathbf{\tilde{C}},\mathbf{\tilde{D}}\)),
based on an impulse input and its corresponding response output. Then,
modal properties can be extracted from \(\mathbf{\tilde{A}}\) and
\(\mathbf{\tilde{C}}\).
With the discrete LTI model, a unit impulse input with zero initial
conditions produces an output of constants
(\(\mathbf{D,CB,CAB,...,CA}^{k-1}\mathbf{B}\)). These constants are
called Markov parameters because they must be unique for a given
system – there is only one possible output for a unit impulse input.
\[\begin{split}\begin{aligned}
\mathbf{x}_{k+1} &= \mathbf{Ax}_{k} + \mathbf{Bu}_{k} \\
\mathbf{y}_{k} &= \mathbf{Cx}_{k} + \mathbf{Du}_{k} \\
\end{aligned}\end{split}\]
\[\begin{split}\begin{aligned}
\mathbf{u}_{0},\mathbf{u}_{1},\mathbf{u}_{2},...,\mathbf{u}_{k} &= \mathbf{I,0,0,...,0} \\
\mathbf{x}_{0},\mathbf{x}_{1},\mathbf{x}_{2},...,\mathbf{x}_{k} &= \mathbf{0,B,AB,...,A}^{k-1}\mathbf{B} \\
\mathbf{y}_{0},\mathbf{y}_{1},\mathbf{y}_{2},...,\mathbf{y}_{k} &= \mathbf{D,CB,CAB,...,CA}^{k-1}\mathbf{B} \\
\end{aligned}\end{split}\]
Knowing that the impulse response output data directly give the Markov
parameters, the data can then be stacked into the generalized blockwise
Hankel matrix \(\mathbf{H}\):
\[\begin{split}\mathbf{H}
=
\begin{bmatrix}
\mathbf{y}_{1} & \mathbf{y}_{2} & \cdots & \mathbf{y}_{m_{c}} \\
\mathbf{y}_{2} & \mathbf{y}_{3} & \cdots & \mathbf{y}_{m_{c}+1} \\
\vdots & \vdots & \ddots & \vdots \\
\mathbf{y}_{m_{o}} & \mathbf{y}_{m_{o}+1} & \cdots & \mathbf{y}_{m_{o}+m_{c}-1} \\
\end{bmatrix}
=
\begin{bmatrix}
\mathbf{CB} & \mathbf{CAB} & \cdots & \mathbf{CA}^{m_{c}-1}\mathbf{B} \\
\mathbf{CAB} & \mathbf{CA}^{2}\mathbf{B} & \cdots & \mathbf{CA}^{m_{c}}\mathbf{B} \\
\vdots & \vdots & \ddots & \vdots \\
\mathbf{CA}^{m_{o}-1}\mathbf{B} & \mathbf{CA}^{m_{o}}\mathbf{B} & \cdots & \mathbf{CA}^{m_{c}+m_{o}-2}\mathbf{B} \\
\end{bmatrix}
=
\mathbf{\mathcal{OC}}\end{split}\]
where \(\mathbf{\mathcal{O}}\) and \(\mathbf{\mathcal{C}}\) are
the observability and controllability matrices of the system:
\[\begin{split}\mathbf{\mathcal{O}} = \begin{bmatrix}
\mathbf{C} \\ \mathbf{CA} \\ \mathbf{CA}^{2} \\ \vdots \\ \mathbf{CA}^{m_{o}-1}
\end{bmatrix}, \hspace{1cm}
\mathbf{\mathcal{C}} = \begin{bmatrix} \mathbf{B} & \mathbf{AB} & \mathbf{A}^{2}\mathbf{B} & \cdots & \mathbf{A}^{m_{c}-1}\mathbf{B} \end{bmatrix}\end{split}\]
The shifted Hankel matrix, \(\mathbf{H'}\) (one time step ahead of
\(\mathbf{H}\)), is shown below:
\[\begin{split}\mathbf{H'}
=
\begin{bmatrix}
\mathbf{y}_{2} & \mathbf{y}_{3} & \cdots & \mathbf{y}_{m_{c}+1} \\
\mathbf{y}_{3} & \mathbf{y}_{4} & \cdots & \mathbf{y}_{m_{c}+2} \\
\vdots & \vdots & \ddots & \vdots \\
\mathbf{y}_{m_{o}+1} & \mathbf{y}_{m_{o}+2} & \cdots & \mathbf{y}_{m_{o}+m_{c}} \\
\end{bmatrix}
=
\begin{bmatrix}
\mathbf{CAB} & \mathbf{CA}^{2}\mathbf{B} & \cdots & \mathbf{CA}^{m_{c}}\mathbf{B} \\
\mathbf{CA}^{2}\mathbf{B} & \mathbf{CA}^{3}\mathbf{B} & \cdots & \mathbf{CA}^{m_{c}+1}\mathbf{B} \\
\vdots & \vdots & \ddots & \vdots \\
\mathbf{CA}^{m_{o}}\mathbf{B} & \mathbf{CA}^{m_{o}+1}\mathbf{B} & \cdots & \mathbf{CA}^{m_{c}+m_{o}-1}\mathbf{B}\\
\end{bmatrix}
=
\mathbf{\mathcal{O}A\mathcal{C}}\end{split}\]
By taking the dominant terms of the singular value decomposition
of \(\mathbf{H}\), and transforming the relationship between
\(\mathbf{H} = \mathbf{\mathcal{OC}}\) and
\(\mathbf{H'} = \mathbf{\mathcal{O}A\mathcal{C}}\), a reduced-order
model is constructed as follows:
\[\begin{split}\mathbf{H} = \mathbf{U}\Sigma\mathbf{V}^{H} =
\begin{bmatrix} \mathbf{\tilde{U}} & \mathbf{U}_{t} \end{bmatrix}
\begin{bmatrix} \tilde{\Sigma} & \mathbf{0} \\ \mathbf{0} & \Sigma_{t} \end{bmatrix}
\begin{bmatrix} \mathbf{\tilde{V}}^{H} \\ \mathbf{V}_{t}^{H} \end{bmatrix}
\approx \mathbf{\tilde{U}}\tilde{\Sigma}\mathbf{\tilde{V}}^{H}\end{split}\]
where the superscript \(H\) denotes conjugate transpose and the
subscript \(t\) indicates elements to be truncated such that only
the first \(r\) dominant singular values in \(\tilde{\Sigma}\)
are retained,
\[\begin{split}\begin{aligned}
\mathbf{\tilde{A}} &= \tilde{\Sigma}^{-1/2}\mathbf{\tilde{U}}^{H}\mathbf{H'\tilde{V}}\tilde{\Sigma}^{-1/2} \\
\mathbf{\tilde{B}} &= \tilde{\Sigma}^{1/2}\mathbf{\tilde{V}}^{H}
\begin{bmatrix} \mathbf{I}_{q} \\ \mathbf{0} \end{bmatrix} \\
\mathbf{\tilde{C}} &= \begin{bmatrix} \mathbf{I}_{p} & \mathbf{0} \end{bmatrix} \mathbf{\tilde{U}}\tilde{\Sigma}^{1/2} \\
\mathbf{\tilde{D}} &= \mathbf{y}_{0}
\end{aligned}\end{split}\]
where \(p\) indicates the number of outputs and \(q\) the number
of inputs, and
\[\begin{split}\begin{aligned}
\mathbf{\tilde{x}}_{k+1} &= \mathbf{\tilde{A}\tilde{x}}_{k} + \mathbf{\tilde{B}u}_{k} \\
\mathbf{y}_{k} &= \mathbf{\tilde{C}\tilde{x}}_{k} + \mathbf{\tilde{D}u}_{k}. \\
\end{aligned}\end{split}\]