PDF
\(\newcommand{\footnotename}{footnote}\) \(\def \LWRfootnote {1}\) \(\newcommand {\footnote }[2][\LWRfootnote ]{{}^{\mathrm {#1}}}\) \(\newcommand {\footnotemark }[1][\LWRfootnote ]{{}^{\mathrm {#1}}}\) \(\let \LWRorighspace \hspace \) \(\renewcommand {\hspace }{\ifstar \LWRorighspace \LWRorighspace }\) \(\newcommand {\mathnormal }[1]{{#1}}\) \(\newcommand \ensuremath [1]{#1}\) \(\newcommand {\LWRframebox }[2][]{\fbox {#2}} \newcommand {\framebox }[1][]{\LWRframebox } \) \(\newcommand {\setlength }[2]{}\) \(\newcommand {\addtolength }[2]{}\) \(\newcommand {\setcounter }[2]{}\) \(\newcommand {\addtocounter }[2]{}\) \(\newcommand {\arabic }[1]{}\) \(\newcommand {\number }[1]{}\) \(\newcommand {\noalign }[1]{\text {#1}\notag \\}\) \(\newcommand {\cline }[1]{}\) \(\newcommand {\directlua }[1]{\text {(directlua)}}\) \(\newcommand {\luatexdirectlua }[1]{\text {(directlua)}}\) \(\newcommand {\protect }{}\) \(\def \LWRabsorbnumber #1 {}\) \(\def \LWRabsorbquotenumber "#1 {}\) \(\newcommand {\LWRabsorboption }[1][]{}\) \(\newcommand {\LWRabsorbtwooptions }[1][]{\LWRabsorboption }\) \(\def \mathchar {\ifnextchar "\LWRabsorbquotenumber \LWRabsorbnumber }\) \(\def \mathcode #1={\mathchar }\) \(\let \delcode \mathcode \) \(\let \delimiter \mathchar \) \(\def \oe {\unicode {x0153}}\) \(\def \OE {\unicode {x0152}}\) \(\def \ae {\unicode {x00E6}}\) \(\def \AE {\unicode {x00C6}}\) \(\def \aa {\unicode {x00E5}}\) \(\def \AA {\unicode {x00C5}}\) \(\def \o {\unicode {x00F8}}\) \(\def \O {\unicode {x00D8}}\) \(\def \l {\unicode {x0142}}\) \(\def \L {\unicode {x0141}}\) \(\def \ss {\unicode {x00DF}}\) \(\def \SS {\unicode {x1E9E}}\) \(\def \dag {\unicode {x2020}}\) \(\def \ddag {\unicode {x2021}}\) \(\def \P {\unicode {x00B6}}\) \(\def \copyright {\unicode {x00A9}}\) \(\def \pounds {\unicode {x00A3}}\) \(\let \LWRref \ref \) \(\renewcommand {\ref }{\ifstar \LWRref \LWRref }\) \( \newcommand {\multicolumn }[3]{#3}\) \(\require {textcomp}\) \(\newcommand {\intertext }[1]{\text {#1}\notag \\}\) \(\let \Hat \hat \) \(\let \Check \check \) \(\let \Tilde \tilde \) \(\let \Acute \acute \) \(\let \Grave \grave \) \(\let \Dot \dot \) \(\let \Ddot \ddot \) \(\let \Breve \breve \) \(\let \Bar \bar \) \(\let \Vec \vec \) \(\newcommand {\LWRsubmultirow }[2][]{#2}\) \(\newcommand {\LWRmultirow }[2][]{\LWRsubmultirow }\) \(\newcommand {\multirow }[2][]{\LWRmultirow }\) \(\newcommand {\mrowcell }{}\) \(\newcommand {\mcolrowcell }{}\) \(\newcommand {\STneed }[1]{}\) \(\def \CH {{\mathcal {H}}}\) \(\def \CL {{\mathcal {L}}}\) \(\newcommand {\bbbred }{\ensuremath {\mathbb {b}}}\) \(\newcommand {\bbcred }{\ensuremath {\mathbb {c}}}\) \(\newcommand {\bbmred }{\ensuremath {\mathbb {m}}}\) \(\newcommand {\trans }{\ensuremath {\mkern -1.5mu\mathsf {T}}}\) \(\def \R {{\mathbb R}}\) \(\def \C {{\mathbb C}}\) \(\def \Cpm {\C ^{p\times m}}\) \(\def \Rn {\R ^n}\) \(\def \Rm {\R ^m}\) \(\def \Rp {\R ^p}\) \(\def \Rr {\R ^r}\) \(\def \Rnn {\R ^{n\times n}}\) \(\def \Rnm {\R ^{n\times m}}\) \(\def \Rpn {\R ^{p\times n}}\) \(\def \Bzero {{\bf 0}}\) \(\def \BA {{\bf A}}\) \(\def \BB {{\bf B}}\) \(\def \BC {{\bf C}}\) \(\def \BH {{\bf H}}\) \(\def \BI {{\bf I}}\) \(\def \BM {{\bf M}}\) \(\def \Bu {{\bf u}}\) \(\def \Bx {{\bf x}}\) \(\def \By {{\bf y}}\) \(\def \wh #1{\widehat {#1}}\) \(\def \BAr {\wh {\bf A}}\) \(\def \BBr {\wh {\bf B}}\) \(\def \Bbr {\wh {\bf b}}\) \(\def \BCr {\wh {\bf C}}\) \(\def \Bcr {\wh {\bf c}}\) \(\def \BDr {\wh {\bf D}}\) \(\def \Bdr {\wh {\bf d}}\) \(\def \BEr {\wh {\bf E}}\) \(\def \BMr {\wh {\bf M}}\) \(\def \BHr {\wh {\BH }}\) \(\def \Bxr {\wh {\bf x}}\) \(\def \BVr {\wh {\bf V}}\) \(\def \BWr {\wh {\bf W}}\) \(\def \yr {\wh {y}}\) \(\def \Byr {\wh {\bf y}}\) \(\def \Sys {\Sigma }\) \(\def \Sysred {\wh \Sys }\) \(\newcommand {\imunit }{{\dot {\imath \hspace *{-0.2em}\imath }}}\)

\(\CH _2\) optimal model reduction of linear systems with quadratic outputs: from rational function interpolation to data-driven modeling

Sean Reiter , Ion Victor Gosea, Igor Pontes Duff, Serkan Gugercin

Abstract

\(\mathcal {H}_2\) optimal reduction of linear dynamical systems represents a long-lasting, worthwhile problem in system-theoretical model order reduction. In this short note, we propose extensions of first-order necessary conditions, both based on system Gramians and transfer functions, to the \(\mathcal {H}_2\) problem for the class of linear systems with quadratic outputs.

1 Introduction

Model-order reduction (MOR) refers to the procedure by which one approximates a large-scale dynamical system, modeled by systems of ordinary differential equations, with a comparatively lower-order surrogate model which can be used as a cheap-to-evaluate surrogate in downstream computational tasks, such as optimization or control. In order to be an effective surrogate, the computed reduced-order model (ROM) should recover the dominant input-to-output response characteristics of the original complex system, as well as preserve qualitative features like internal structures. We refer to [1, 2] for more details on system-theoretical MOR, since this category is of interest to us.

The primary consideration of this work is the development of methods for the MOR of linear dynamical systems which contain quadratic output functions, or linear quadratic-output (LQO) systems. Here, we develop extensions of classical MOR approaches applicable solely to systems with linear dynamics and linear outputs. Our contributions are threefold: First, we consider the \(\CH _2\) optimal model reduction problem, and derive first-order necessary conditions for \(\CH _2\) optimality based on rational transfer function interpolation. These provide a natural extension of the well-known interpolation-based \(\CH _2\) optimality framework of Meier and Leunberger  [7, 8] for linear model reduction. Based on the developed theoretical optimality framework, we propose an extension of the well-known iterative rational Krylov algorithm (IRKA)  [7] for linear \(\CH _2\) optimal model reduction. Finally, we show how to compute \(\CH _2\) optimal reduced models using only evaluations of the linear- and quadratic-output transfer functions.

2 Transfer functions, norms and MOR of LQO systems

In this work, we consider large-scale dynamical systems with linear dynamics and outputs which are (up to) quadratic functions of the state vector. In state-space, such systems are formulated as

\begin{align} \label {eq:lqosys} \Sys : \quad \begin{cases} \dot {\Bx }(t)=\BA \Bx (t)+\BB \Bu (t), \\ \By (t)={\BC \Bx (t)} + {\BM \left (\Bx (t)\otimes \Bx (t)\right )}, \end {cases} \end{align} where \(\Bx (t)\in \Rn \) is the system’s internal state, \(\Bu (t) \in \Rm \) are the control inputs, and \(\By (t)\in \Rp \) are the observed outputs. The system matrices satisfy \(\BA \in \Rnn \), \(\BB \in \Rnm \), \(\BC \in \Rpn \) and \(\BM \in \R ^{p\times n^2}\). We assume that the system in (1) is asymptotically stable, i.e., the eigenvalues of \(\BA \) to be in the left-half plane. Systems that consider quadratic observables as quantities of interest arise in a variety of applications, and particularly whenever one is interested in observing quantities computed as the product of time or frequency-domain components of the state  [4, 11].

The frequency-domain response of system (1) is fully specified by \(2\) rational transfer functions [4, 6]:

\begin{align} \label {eq:tfs} \BH _1(s)=\BC (s\BI -\BA )^{-1}\BB ~~~\mbox {and}~~~\BH _2(s,z)=\BM \left ((s \BI -\BA )^{-1}\BB \otimes (z\BI -\BA )^{-1}\BB \right ). \end{align} The first function \(\BH _1(s)\) is the typical transfer function of a linear-output system, and describes the transfer from input \(\Bu (t)\) to output \(\By _1(t):= \BC \Bx (t)\), which is linear in \(\Bx (t)\). The second bivariate function \(\BH _2(s,z)\) the transfer from input \(\Bu (t)\) to output \(\By _2(t):= \BM \left (\Bx (t)\otimes \Bx (t)\right )\), which is quadratic in \(\Bx (t)\). The \(\CH _2\) norm for systems of the form (1) can be defined via these transfer functions as  [4, 5]

\begin{align} \label {eq:h2} \|\Sigma \|_{\CH _2}^2:= \frac {1}{2\pi }\int _{-\infty }^\infty \|{\BH _1(\imunit \omega )}\|_{\textsf {F}}^2 \, d\omega + \frac {1}{(2\pi )^2}\int _{-\infty }^\infty \int _{-\infty }^\infty \|{\BH _2(\imunit \omega _1,\imunit \omega _2)}\|_{\textsf {F}}^2 \, d\omega _1 \, d\omega _2. \end{align} We note that when \(\BM = \mathbf {0}\), it implies that \(\BH (s,z) = 0\), and the formula above simplifies to the first integral term only, which is the standard formula as used in [7].

In practical applications, the state dimension \(n\) can be rather large, e.g., in the order of the millions, and any repeated action involving the full-order model (FOM) (1) becomes prohibitively expensive. Model reduction seeks to remedy this problem with the construction of cheap-to-evaluate surrogate models having the same form as (1), but described by a comparatively much smaller number of differential equations. Mathematically, this amounts to computing a system

\begin{align} \label {eq:lqosysred} \Sysred : \quad \begin{cases} \dot {\Bxr }(t)=\BAr \Bxr (t)+\BBr \Bu (t), \\ \Byr (t)={\BCr \Bxr (t)} + {\BMr \left (\Bxr (t)\otimes \Bxr (t)\right )}, \end {cases} \end{align} with a significantly reduced dimension \(1\leq r\ll n\). \(\Bxr (t) \in \Rr \) contains the reduced-order state variables, and \(\Byr (t) \in \Rp \) are the approximateds output. The reduced-order matrix operators satisfy \(\BAr \in \R ^{r\times r}\), \(\BBr \in \R ^{r\times m}\), \(\BCr \in \R ^{p\times r}\), and \(\BMr \in \R ^{p\times r^2}\). In order to be an effective surrogate, the reduced model (4) should replicate the input-to-output response characteristics of the large-scale system (1). In order words, for a given tolerance \(\tau > 0\), the output deviation should satisfy \(\|\By -\Byr \|\leq \tau \|\Bu \|\) in an appropriate norm for a range of admissible inputs \(\Bu \).

Suppose that one is interested in controlling the \(\CL ^p_\infty \), or “worst case” deviation in the output \(\|\By -\Byr \|_{\CL ^p_\infty }:=\sup _{t\geq 0}\|\By (t)-\Byr (t)\|_\infty \). Significantly, one can show following error bound  [4]:

\begin{equation} \label {eq:bound} \|{\By -\Byr }\|_{\CL _\infty ^p} \leq {\|{\Sys -\Sysred }\|_{\CH _2}} \big (\|\Bu \|_{\CL ^m_2}^2+\|\Bu \otimes \Bu \|_{\CL ^{m^2}_2}^2\big )^{1/2}. \end{equation}

In other words, the \(\CH _2\) model error bounds the \(\CL ^p_\infty \) output error. Based on the bound (5), we consider the \(\CH _2\) optimal model reduction problem for the LQO system class (1). Given the system in (1), find a ROM such that

\begin{equation} \label {eq:H2optimalMOR} \min _{\dim ({\Sys })=r} \mathcal {J}(\Sysred ),~~~\mathcal {J}(\Sysred ):=\|\Sys -{\Sysred }\|_{\CH _2}^2. \end{equation}

3 One main result

We follow here the results in [9]. To simplify the approximation problem, we assume the approximate system in (4) has simple poles. Then, the reduced-order linear- and quadratic-output transfer functions can be expressed in pole-residue form:

\begin{equation} \label {eq:poleResidue} \widehat {\BH }_1(s)=\sum _{i=1}^r\frac {\bbcred _i\bbbred _i^{\trans }}{s-\lambda _i}~~~\mbox {and}~~~\widehat {\BH }_2(s,z)=\sum _{j=1}^r\sum _{k=1}^r\frac {\bbmred _{j,k}\left (\bbbred _j\otimes \bbbred _k\right )^{\trans }}{(s-\lambda _j)(z-\lambda _k)}, \end{equation}

where \(\bbbred _k\in \C ^{m}\), and \(\bbcred _k\in \C ^p,\bbmred _{j,k}\in \C ^{p}\), for all \(j,k=1,\ldots ,r\). We define \(\bbcred _i\bbbred _i^{\trans }\in \Cpm \) and \(\bbmred _{j,k}(\bbbred _i\otimes \bbbred _i)^{\trans }\in \C ^{p\times m^2}\) to be the residues of \(\widehat {\BH }_1(s)\) and \(\widehat {\BH }_2(s,z)\) corresponding to \(\lambda _i\) and \((\lambda _j,\lambda _k)\), respectively. We are able to show that

\begin{equation} \label {eq:h2error} \|\Sys -\Sysred \|_{\CH _2}^2=\|\Sys \|_{\CH _2}^2 -2\left ( \sum _{i=1}^r {\bbcred _i^{\trans }} \BH _1({-\lambda _i}){\bbbred _i} +\sum _{j=1}^r\sum _{k=1}^r {\bbmred _{j,k}^{\trans }} \BH _2({-\lambda _j},-{\lambda _k})\left ({\bbbred _j}\otimes {\bbbred _k}\right ) \right )+ {\|\Sysred \|_{\CH _2}^2}. \end{equation}

This makes the \(\CH _2\) optimal model reduction problem tractable by minimally parameterizing the ROM in (4) in terms of the transfer function poles and residues.

  • Theorem 1 Suppose that \(\Sysred \) has simple poles \({\lambda _1,\ldots ,\lambda _r}\in \C _-\), and is a local minimizer of the squared \(\CH _2\) error \(\|\Sys -\Sysred \|_{\CH _2}^2\). Then, for all \(i,j,k=1,\ldots ,r\), it holds that

    \begin{align*} \Bzero &=\left (\BH _1({-\lambda _i})-\BHr _1({-\lambda _i})\right ){\bbbred _i},\\ \Bzero &= \left (\BH _2({-\lambda _j},{-\lambda _k})-\BHr _2({-\lambda _j},{-\lambda _k})\right )\left ({\bbbred _{j}}\otimes {\bbbred _k}\right ),\\ \Bzero &= {\bbcred _k^{\trans }}\left (\BH _1({-\lambda _k})-\BHr _1({-\lambda _k})\right ) +\sum _{\ell =1}^r{\bbmred _{k,\ell }^{\trans }} \left (\BH _2({-\lambda _k},{-\lambda _\ell })-\BHr _2({-\lambda _k},{-\lambda _\ell })\right )\left (\BI _m\otimes {\bbbred _\ell }\right ) \\[-2ex] &~~~~~~~~~+\sum _{\ell =1}^r{\bbmred _{\ell ,k}^{\trans }}\left (\BH _2({-\lambda _\ell },{-\lambda _k})-\BHr _2({-\lambda _\ell },{-\lambda _k})\right )\left ({\bbbred _\ell }\otimes \BI _m\right ),\\ 0&= {\bbcred _k^{\trans }}\left (\frac {d}{ds}\BH _1({-\lambda _k})-\frac {d}{ds}\BHr _1({-\lambda _k})\right ){\bbbred _k} +\sum _{\ell =1}^r{\bbmred _{k,\ell }^{\trans }} \left (\frac {\partial }{\partial s_1}\BH _2({-\lambda _k},{-\lambda _\ell })-\frac {\partial }{\partial s_1}\BHr _2({-\lambda _k},{-\lambda _\ell })\right )\left ({\bbbred _k}\otimes {\bbbred _\ell }\right ) \\[-1ex] &~~~~~~~~~+\sum _{\ell =1}^r{\bbmred _{\ell ,k}^{\trans }}\left (\frac {\partial }{\partial s_2}\BH _2({-\lambda _\ell },{-\lambda _k})-\frac {\partial }{\partial s_2}\BHr _2({-\lambda _\ell },{-\lambda _k})\right )\left ({\bbbred _\ell }\otimes {\bbbred _k}\right ). \end{align*}

In other words, tangential interpolation is a necessary condition for \(\CH _2\) optimality. We also note that when \(\BM = \mathbf {0}\), it implies that the formulae above simplify accordingly to the standard interpolation-based FONC’s for classical linear systems, as in [7, 8].

4 Summary of all proposed results

Based on the bound in (5), we have considered the \(\CH _2\) optimal model reduction problem for the class of systems in (1). We went about this in two different ways, corresponding to two types of FONCs, namely for the first one mentioned earlier introduced in [8]. Then, for the second, e.g., the Gramian-based FONCs as introduced in [13], we analyzed it in [10].

Our contributions to the interpolation-based formulation are threefold:

Due to space limitations, we are not able to go into a detailed analysis of the results concerning Gramian-based FONCs, as presented in [10]. In short, we derive gradients of the squared \(\CH _2\) system error with respect to the system matrices of the LQO-ROM as parameters. The stationary points of these gradients directly yield Gramian-based FONCs for \(\CH _2\) optimality. These results generalize the analogous Gramian-based FONCs for linear \(\CH _2\) optimal model  [12, 13] to the LQO setting. We also show that a \(\CH _2\) optimal LQO-ROM is necessarily defined by Petrov-Galerkin projection. The relevant projection matrices are obtained as solutions to a pair of Sylvester equations.

References