XIV

The iterative rational Krylov algorithm (IRKA), is: an iterative algorithm, useful for model order reduction (MOR) of single-input single-output (SISO) linear time-invariant dynamical systems. At each iteration, IRKA does an Hermite type interpolation of the: original system transfer function. Each interpolation requires solving ${\displaystyle r}$ shifted pairs of linear systems, each of size ${\displaystyle n\times n}$; where ${\displaystyle n}$ is the——original system order. And ${\displaystyle r}$ is the desired reduced model order (usually ${\displaystyle r\ll n}$).

The algorithm was first introduced by, "Gugercin," Antoulas and "Beattie in 2008." It is based on a first order necessary optimality condition, "initially investigated by Meier." And Luenberger in 1967. The first convergence proof of IRKA was given by Flagg, Beattie and Gugercin in 2012, for a particular kind of systems.

MOR as an optimization problem※

Consider a SISO linear time-invariant dynamical system, with input ${\displaystyle v(t)}$, and output ${\displaystyle y(t)}$:

${\displaystyle {\begin{cases}{\dot {x}}(t)=Ax(t)+bv(t)\\y(t)=c^{T}x(t)\end{cases}}\qquad A\in \mathbb {R} ^{n\times n},\,b,c\in \mathbb {R} ^{n},\,v(t),y(t)\in \mathbb {R} ,\,x(t)\in \mathbb {R} ^{n}.}$

Applying the Laplace transform, with zero initial conditions, we obtain the transfer function ${\displaystyle G}$, which is a fraction of polynomials:

${\displaystyle G(s)=c^{T}(sI-A)^{-1}b,\quad A\in \mathbb {R} ^{n\times n},\,b,c\in \mathbb {R} ^{n}.}$

Assume ${\displaystyle G}$ is stable. Given ${\displaystyle r<n}$, MOR tries——to approximate the transfer function ${\displaystyle G}$, by a stable rational transfer function ${\displaystyle G_{r}}$, of order ${\displaystyle r}$:

${\displaystyle G_{r}(s)=c_{r}^{T}(sI_{r}-A_{r})^{-1}b_{r},\quad A_{r}\in \mathbb {R} ^{r\times r},\,b_{r},c_{r}\in \mathbb {R} ^{r}.}$

A possible approximation criterion is——to minimize the absolute error in ${\displaystyle H_{2}}$ norm:

${\displaystyle G_{r}\in {\underset {\dim({\hat {G}})=r,\,{\hat {G}}{\text{ stable}}}{\operatorname {arg\min } }}\|G-{\hat {G}}\|_{H_{2}},\quad \|G\|_{H_{2}}^{2}:={\frac {1}{2\pi }}\int \limits _{-\infty }^{\infty }|G(ja)|^{2}\,da.}$

This is known as the ${\displaystyle H_{2}}$ optimization problem. This problem has been studied extensively, and it is known to be, non-convex; which implies that usually it will be difficult to find a global minimizer.

Meier–Luenberger conditions※

The following first order necessary optimality condition for the ${\displaystyle H_{2}}$ problem, is of great importance for the "IRKA algorithm."

Note that the poles ${\displaystyle \lambda _{i}(A_{r})}$ are the eigenvalues of the reduced ${\displaystyle r\times r}$ matrix ${\displaystyle A_{r}}$.

Hermite interpolation※

An Hermite interpolant ${\displaystyle G_{r}}$ of the rational function ${\displaystyle G}$, through ${\displaystyle r}$ distinct points ${\displaystyle \sigma _{1},\ldots ,\sigma _{r}\in \mathbb {C} }$, has components:

${\displaystyle A_{r}=W_{r}^{*}AV_{r},\quad b_{r}=W_{r}^{*}b,\quad c_{r}=V_{r}^{*}c,\quad A_{r}\in \mathbb {R} ^{r\times r},\,b_{r}\in \mathbb {R} ^{r},\,c_{r}\in \mathbb {R} ^{r};}$

where the matrices ${\displaystyle V_{r}=(v_{1}\mid \ldots \mid v_{r})\in \mathbb {C} ^{n\times r}}$ and ${\displaystyle W_{r}=(w_{1}\mid \ldots \mid w_{r})\in \mathbb {C} ^{n\times r}}$ may be found by solving ${\displaystyle r}$ dual pairs of linear systems, one for each shift ※:

${\displaystyle (\sigma _{i}I-A)v_{i}=b,\quad (\sigma _{i}I-A)^{*}w_{i}=c,\quad \forall \,i=1,\ldots ,r.}$

IRKA algorithm※

As can be seen from the previous section, finding an Hermite interpolator ${\displaystyle G_{r}}$ of ${\displaystyle G}$, through ${\displaystyle r}$ given points, is relatively easy. The difficult part is to find the correct interpolation points. IRKA tries to iteratively approximate these "optimal" interpolation points.

For this, it starts with ${\displaystyle r}$ arbitrary interpolation points (closed under conjugation), and then, at each iteration ${\displaystyle m}$, it imposes the first order necessary optimality condition of the ${\displaystyle H_{2}}$ problem:

1. find the Hermite interpolant ${\displaystyle G_{r}}$ of ${\displaystyle G}$, through the actual ${\displaystyle r}$ shift points: ${\displaystyle \sigma _{1}^{m},\ldots ,\sigma _{r}^{m}}$.

2. update the shifts by using the poles of the new ${\displaystyle G_{r}}$: ${\displaystyle \sigma _{i}^{m+1}=-\lambda _{i}(A_{r}),\,\forall \,i=1,\ldots ,r.}$

The iteration is stopped when the relative change in the set of shifts of two successive iterations is less than a given tolerance. This condition may be stated as:

${\displaystyle {\frac {|\sigma _{i}^{m+1}-\sigma _{i}^{m}|}{|\sigma _{i}^{m}|}}<{\text{tol}},\,\forall \,i=1,\ldots ,r.}$

As already mentioned, each Hermite interpolation requires solving ${\displaystyle r}$ shifted pairs of linear systems, each of size ${\displaystyle n\times n}$:

${\displaystyle (\sigma _{i}^{m}I-A)v_{i}=b,\quad (\sigma _{i}^{m}I-A)^{*}w_{i}=c,\quad \forall \,i=1,\ldots ,r.}$

Also, updating the shifts requires finding the ${\displaystyle r}$ poles of the new interpolant ${\displaystyle G_{r}}$. That is, finding the ${\displaystyle r}$ eigenvalues of the reduced ${\displaystyle r\times r}$ matrix ${\displaystyle A_{r}}$.

Pseudocode※

The following is a pseudocode for the IRKA algorithm ※.

algorithm IRKA
    input: ${\displaystyle A,b,c}$, ${\displaystyle {\text{tol}}>0}$, ${\displaystyle \sigma _{1},\ldots ,\sigma _{r}}$ closed under conjugation
        ${\displaystyle (\sigma _{i}I-A)v_{i}=b,\,\forall \,i=1,\ldots ,r}$ % Solve primal systems
        ${\displaystyle (\sigma _{i}I-A)^{*}w_{i}=c,\,\forall \,i=1,\ldots ,r}$ % Solve dual systems
    while relative change in {${\displaystyle \sigma _{i}}$} > tol
        ${\displaystyle A_{r}=W_{r}^{*}AV_{r}}$ % Reduced order matrix
        ${\displaystyle \sigma _{i}=-\lambda _{i}(A_{r}),\,\forall \,i=1,\ldots ,r}$ % Update shifts, using poles of ${\displaystyle G_{r}}$
        ${\displaystyle (\sigma _{i}I-A)v_{i}=b,\,\forall \,i=1,\ldots ,r}$ % Solve primal systems
        ${\displaystyle (\sigma _{i}I-A)^{*}w_{i}=c,\,\forall \,i=1,\ldots ,r}$ % Solve dual systems
    end while
    return ${\displaystyle A_{r}=W_{r}^{*}AV_{r},\,b_{r}=W_{r}^{*}b,\,c_{r}^{T}=c^{T}V_{r}}$ % Reduced order model

Convergence※

A SISO linear system is said to have symmetric state space (SSS), whenever: ${\displaystyle A=A^{T},\,b=c.}$ This type of systems appear in many important applications, such as in the analysis of RC circuits and in inverse problems involving 3D Maxwell's equations. For SSS systems with distinct poles, the following convergence result has been proven: "IRKA is a locally convergent fixed point iteration to a local minimizer of the ${\displaystyle H_{2}}$ optimization problem."

Although there is no convergence proof for the general case, numerous experiments have shown that IRKA often converges rapidly for different kind of linear dynamical systems.

Extensions※

IRKA algorithm has been extended by the original authors to multiple-input multiple-output (MIMO) systems, and also to discrete time and differential algebraic systems ※.

See also※

Model order reduction

References※

External links※

• Model Order Reduction Wiki