Loewner Framework

Note: This page has not been verified by our editors.

The Loewner Framework is a frequency-domain system identification, system realization and model reduction method, named after Charles Loewner, who introduced the Loewner Matrix which is essential to this method based on ideas from .

Introduction

Model order reduction is commonly used to simulate and control complex physical processes. This is done by replacing the large-scale system model described in terms of differential equations by a system of much lower dimension that has similar response characteristics. In many situations, a model (the collection of differential equations written in matrix format) for the underlying dynamical process is not explicitly available. Hence, the system matrix realization can not be directly obtained. Instead, a black box is provided that produces input-output measurements. These situations include VLSI modeling from chips or real time simulation of multi-body dynamics with constraints.

We study reduction methods based on moment matching, that is, matching of the coefficients of power series expansions of the transfer function at selected points in the complex plane. More precisely, values of the underlying transfer function, together with the values of its derivatives are interpolated. Hence, the main underlying problem is rational interpolation. Different methods based on Krylov iteration, as well as the Arnoldi or the Lanczos procedures, and multi-point (rational) versions thereof.

Interpolatory model reduction methods compute ROMs (reduced order models) by means of interpolating the response of the original system at selected points (usually the transfer function). The Loewner framework (introduced in ) falls into this category. It is a data-driven method that directly constructs a ROM using only measured data (samples of the transfer function). Additionally, by compressing the (usually large) data set, it extracts the dominant features and eliminates the inherent redundancies. For a recent survey on the Loewner framework for linear systems, see . For the nonlinear case, this method requires appropriate definition of transfer functions, but the main procedure does not differ too much from the one applied for the linear case. Some of the extensions include the bilinear systems case in , quadratic-bilinear systems and linear switched systems . The dissertation represents a comprehensive collection of these methods.

The motivation behind developing generalized versions of the Loewner framework to new systems is that such classes of systems can be viewed as a bridge between linear and nonlinear systems. For example, one can always write an approximation of a nonlinear system by means of a bilinear system (Carleman linearization ). Furthermore, for certain types of nonlinear systems, one can always find an equivalent quadratic-bilinear model without performing any approximation (McCormick relaxation ) whatsoever. Linear switched systems can be viewed as a special class of hybrid systems that are used to model coupled or distributed processes characterized by both discrete and continuous dynamics.

Loewner Framework for Linear Systems

We analyze linear time-invariant continuous-time systems described by the following equations:

$\Sigma : \begin{cases} E \dot{x}(t) = Ax(t) + Bu(t) \\ \;\;\; y(t) = Cx(t) + Du(t) \end{cases} \quad(1)$

where $A, E \in \mathbb{R}^{n \times n}$, $C \in \R^{p \times n}$, $B \in \mathbb{R}^{n \times m}$, $D \in \mathbb{R}^{p \times m}$ and $x(t) \in \R^n$, $u(t) \in R^m$, $y(t) \in \R^p$ are state, input and output vector respectively. The transfer function corresponding to $\Sigma:(A,B,C,D,E)$ is computed as $H(s) = C(sE-A)^{-1}B + D$.

The Loewner framework can construct a reduced-order system by means of measured data. In this context, data represent sample values of the underlying transfer function $H(s)$ at selected interpolation points (in some cases on the imaginary axis, when approximating the frequency response of the original system). A Loewner matrix pencil $\lambda\mathbb{L}-\mathbb{L}_s$ is composed of the Loewner matrix $\mathbb{L}$ and the shifted Loewner matrix $\mathbb{L}_s$. This pencil represents a surrogate of the original matrix pencil $\lambda E - A$.

The Loewner framework originates from rational approximation theory (see ) and can be used to recover or reduce the original dynamical system by means of interpolating the transfer function. For simplicity, assume that the data set contains an even number of measurements denoted with $2N$. The next step is to partition the data measurements set into two disjoint data subsets, i.e. the left subset and the right subset as described in the following

Left set: $(\mu_i,\ell_i,v_i), i = 1,\dots,N$, Right set: $(\lambda_j,r_j,w_j), j=1,\dots,N. \quad (2)$

More precisely, the data in (2) is composed of $N$ left interpolation points: $\{\mu_i\}_{i=1}^N \subset \mathbb{C}$, $N$ left sample values: $\{v_i\}_{i=1}^N \subset \mathbb{C}^m$ and $N$ left tangential directions: $\{\ell_i\}_{i=1}^N \subset \mathbb{C}^p$. Additionally, we have $N$ distinct right interpolation points: $\{\lambda_i\}_{i=1}^N \subset \mathbb{C}$ with $N$ right sample values: $\{w_i\}_{i=1}^N \subset \mathbb{C}^p$ and right tangential directions: $\{r_i\}_{i=1}^N \subset \mathbb{C}^m$. We seek a linear system $\widehat{\Sigma}$ such that the associated transfer function $\widehat{H}(s)$ is a tangential interpolant to the original one, $H(s)$. More precisely, the following interpolation conditions should then hold:

$\begin{array}{rcl} \ell_i^* \widehat{H}(\mu_i) &=& \ell_i^*H(\mu_i) = v_i^* \quad \text{for} \;\; i = 1,\dots,N, \\ \widehat{H}(\lambda_j)r_j &=& H(\lambda_j) r_j = w_j \quad \text{for} \;\; j = 1,\dots,N. \qquad(3) \end{array}$

Next, arrange the measured data into matrix format as follows

$\begin{array}{rcl} M &=& \operatorname{diag}([\mu_1,\mu_2,\dots,\mu_N]) \in \mathbb{C}^{N \times N}, \quad \Lambda = \operatorname{diag}([\lambda_1,\lambda_2,\dots,\lambda_N]) \in \mathbb{C}^{N \times N}, \\ L^* &=& [\ell_1 \; \ell_2 \; \dots \; \ell_N] \in \mathbb{C}^{p \times N}, \quad R = [r_1 \; r_2 \; \dots \; r_N] \in \mathbb{C}^{m \times N}, \\ V^* &=& [v_1 \; v_2 \; \dots \; v_N] \in \mathbb{C}^{m \times N}, \quad W = [w_1 \; w_2 \; \dots \; w_N] \in \mathbb{C}^{p \times N}. \end{array} \quad(4)$

The Loewner matrix is constructed in the following way

$\mathbb{L} = \begin{bmatrix} \frac{v_1^\intercal r_1 - \ell_1^\intercal w_1}{\mu_1 - \lambda_1} & \dots & \frac{v_1^\intercal r_N - \ell_1^\intercal w_N}{\mu_1 - \lambda_N} \\ \vdots & \ddots & \vdots \\ \frac{v_N^\intercal r_1 - \ell_N^\intercal w_1}{\mu_N - \lambda_1} & \dots & \frac{v_N^\intercal r_N - \ell_N^\intercal w_N}{\mu_N - \lambda_N} \end{bmatrix} \in \mathbb{C}^{N \times N} \quad(5)$

and is the solution of the following Sylvester equation $M\mathbb{L} - \mathbb{L}\Lambda = VR - LW$. The shifted Loewner matrix can be constructed in the following way:

$\mathbb{L}_s = \begin{bmatrix} \frac{v_1^\intercal r_1 \mu_1 - \ell_1^\intercal w_1 \lambda_1}{\mu_1 - \lambda_1} & \dots & \frac{v_1^\intercal r_N \mu_1 - \ell_1^\intercal w_N \lambda_N}{\mu_1 - \lambda_N} \\ \vdots & \ddots & \vdots \\ \frac{v_N^\intercal r_1 \mu_N - \ell_N^\intercal w_1 \lambda_1}{\mu_N - \lambda_1} & \dots & \frac{v_N^\intercal r_N \mu_N - \ell_N^\intercal w_N \lambda_N}{\mu_N - \lambda_N} \end{bmatrix} \in \mathbb{C}^{N \times N}, \quad(6)$

and is the solution of the following Sylvester equation $M\mathbb{L}_s - \mathbb{L}_s\Lambda = MVR - LW\Lambda$. After rearranging the data into matrix format, we construct with basically no computational cost a Loewner linear model described by the following matrices

$\Sigma_L = (-\mathbb{L}_s, V, W, 0, -\mathbb{L})$.

The $D$ term does not play an important role in the Loewner framework since it is incorporated in the Loewner matrix $\mathbb{L}$ (see ). Consequently, this term is taken to be $0$ in (7). Note that, in the SISO case we consider only one input and one output $(m = p = 1)$. Then one can take $R = L^\intercal = [1 \; 1 \; \dots \; 1] \in \mathbb{R}^N$ and explicitly write

$\mathbb{L}(i,j) = \frac{H(\mu_i) - H(\lambda_j)}{\mu_i - \lambda_j}, \quad \mathbb{L}_s(i,j) = \frac{\mu_i H(\mu_i) - \lambda_j H(\lambda_j)}{\mu_i - \lambda_j}, \quad V(i,1) = H(\mu_i), \quad W(1,j) = H(\lambda_j). \quad (8)$

The Loewner model can completely recover the dynamics of the original system if enough measurements are provided. We can hence distinguish two cases:

1. The right amount of data is available - in this case one can exactly recover the original system;
2. A redundant amount of data is provided - in this case one can approximate the original system.

In the first case, the matrix pencil $\lambda \mathbb{L} - \mathbb{L}_s$ is regular. Then, the transfer function of the Loewner model can be computed as follows $\widehat{H}(s) = -W(s\mathbb{L}-\mathbb{L}_s)^{-1} V$. This situation is encountered when the number of measurements is enough to recover the original system (one would require exactly $2n$ measurements to recover a system of order $n$). In the second case, the pencil $\lambda \mathbb{L} - \mathbb{L}_s$ is singular. Then, the definition of the transfer function can be adapted as in , i.e. $\widehat{H}(s) = -W(s\mathbb{L}-\mathbb{L}_s)^\# V$. Here, the $\#$ operator represents the Moore-Penrose or Drazin generalized inverse. In real world applications, the reason behind the singularity of the matrix pencil $\lambda \mathbb{L} - \mathbb{L}_s$ is the overfitting of the data that contains in most cases redundant or inaccurate measurements.

We focus on the second case, which is more relevant in practice. In this case, since the Loewner pencil is singular, we can use the singular value decomposition (SVD) to compress the data. The truncation order $r$ can be chosen by inspecting the singular value decay of the data. Deciding on a specific value of $r$ is usually done by making a trade-off between reduced-order model dimension and accuracy of fit. Consider the SVD of the Loewner matrix $\mathbb{L}$ as $\mathbb{L} = Y_1 S_1 X_1^*$ where $S_1, Y_1, X_1 \in \mathbb{C}^{N \times N}$. The singular values are nonnegative real scalars that give information on the numerical rank of matrices and are stored on the main diagonal of matrix $S_1$, as $S_1 = \operatorname{diag}([\sigma_1,\sigma_2,\dots,\sigma_N$ with $\sigma \geq \sigma_2\geq\dots\geq\sigma_N\geq 0$. the matrices $Y,X \in \mathbb{C}^{N \times r}$ are obtained by selecting the first $r$ columns of the matrices $Y_1$ and $X_1$, respectively. The reduced Loewner system is constructed by projecting the matrices $\mathbb{L}$, $\mathbb{L}_s$, $V$, $W$ with $Y^*$ and $X$ to the left and respectively, to the right, as

$\hat{\mathbb{L}} = Y^*\mathbb{L}X, \quad \hat{\mathbb{L}}_s = Y^* \mathbb{L}_s X, \quad \widehat{V} = Y^* V, \quad \widehat{W} = WX. \quad(9)$

The dimension $r$ of the the reduced system is typically chosen as follows; first set a tolerance value $\tau > 0$ and then select the truncation order $r$ as the maximum number of normalized singular values that are larger than $\tau$. More precisely, $\frac{\sigma_i}{\sigma_1} > \tau$ for $i \in \{1,2,\dots,r\}$. The projected Loewner model is given by

$\widehat{\Sigma}_L = (\hat{\mathbb{L}}_s,\widehat{V},\widehat{W},0,\hat{\mathbb{L}}). \quad (10)$

Given the SISO linear system as in (1), consider the tuples of left and right interpolation points: $[\mu_1,\mu_2], [\lambda_1,\lambda_2]$. The recovered matrices can be written in terms of data samples in the following way:

$\mathbb{L} = \begin{bmatrix} \frac{H(\mu_1)-H(\lambda_1)}{\mu_1-\lambda_1} & \frac{H(\mu_1)-H(\lambda_2)}{\mu_1-\lambda_2} \\ \frac{H(\mu_2)-H(\lambda_1)}{\mu_2-\lambda_1} & \frac{H(\mu_2)-H(\lambda_2)}{\mu_2-\lambda_2} \end{bmatrix}, \quad \mathbb{L}_s = \begin{bmatrix} \frac{\mu_1 H(\mu_1)-\lambda_1 H(\lambda_1)}{\mu_1-\lambda_1} & \frac{\mu_1 H(\mu_1)-\lambda_2 H(\lambda_2)}{\mu_1-\lambda_2} \\ \frac{\mu_2 H(\mu_2)-\lambda_1 H(\lambda_1)}{\mu_2-\lambda_1} & \frac{\mu_2 H(\mu_2)-\lambda_2 H(\lambda_2)}{\mu_2-\lambda_2} \end{bmatrix}, \quad V = \begin{bmatrix} H(\mu_1) \\ H(\mu_2) \end{bmatrix}, \quad W = \begin{bmatrix} H(\lambda_1) & H(\lambda_2) \end{bmatrix}.$

This reduced system matches $2k = 4$, for $k = 2$ moments of the original system, namely:

$\text{four of} \; H: \; H(\mu_1), \; H(\mu_2), \; H(\lambda_1), \; H(\lambda_2).$

Loewner Framework for Bilinear Systems

Next, we turn our attention to bilinear systems $\Sigma_B = (A,B,C,0,E,N)$ characterized by equations:

$\Sigma_B : \begin{cases} E\dot{x}(t) &= Ax(t) + Nx(t)u(t) + Bu(t), \\ y(t) &= Cx(t), \end{cases} \qquad (11)$

where $N \in \mathbb{R}^{n \times n}$ and the other quantities are exactly as in (1). For simplicity of exposition, we will assume single-input and single-output scenarios (so $m = p = 1$) and that the $D$ term is incorporated into the $E$ matrix (so $D = 0$). Bilinear systems can be equivalently rewritten as an infinite set of systems of the form:

$\begin{cases} E\dot{x}_1(t) &= Ax_1(t), \\ E\dot{x}_2(t) &= Ax_2(t) + Nx_1(t)u(t), \\ &\vdots \\ E\dot{x}_i(t) &= Ax_i(t) + Nx_{i-1}(t)u(t), \\ &\vdots \end{cases} \qquad (12)$

The external representation of the bilinear system $\Sigma_B$ can be expressed in terms of Volterra series, which describe the nonlinear mapping of admissible inputs $u$ to the output $y$. The solution of (11) can be written as $x(t) = \sum_{i=1}^\infty x_i(t)$. By considering $x_{i-1}(t)$ as a pseudo-input for the $i$th equation in (12), the frequency-domain behavior is described by a series of generalized transfer functions obtained by taking the multivariate Laplace transform of the degree $\ell$ regular kernel. More specifically, introduce the generalized multivariate transfer functions of system $\Sigma_B$, as follows (for $\ell = 1,2 \dots$)

$H_\ell(s_1,s_2,\dots,s_\ell) = C(s_1E-A)^{-1}N(s_2E-A)^{-1}N \dots N(s_\ell E-A)^{-1}B. \quad (13)$

The explicit derivation of these types of transfer functions is based on the so-called Volterra series representation. The characterization of bilinear systems by means of rational functions suggests that reduction of such systems can be performed by means of interpolatory methods, i.e. by extending the Loewner framework.

Let $\{\mu_1,\mu_2,\dots,\mu_N\}$ and $\{\lambda_1,\lambda_2,\dots,\lambda_N\}$ be two set of left and respectively, right interpolation points in the complex plane. Using samples of transfer functions in (13) evaluated at these points, we put together a realization of an order $N$ bilinear system. The main property is that this reduced system matches those original samples. Hence, for $i,j \in \{1,2,\dots,N\}$, one can explicitly write

$\begin{array}{rcl} \mathbb{L}(i,j) &=& \frac{H_{i+j-1}(\mu_1,\dots,\mu_i,\lambda_{j-1},\dots,\lambda_1) - H_{i+j-1}(\mu_1,\dots,\mu_{i-1},\lambda_j,\dots,\lambda_1)}{\mu_i-\lambda_j},\\ \mathbb{L}_s(i,j) &=& \frac{\mu_i H_{i+j-1}(\mu_1,\dots,\mu_i,\lambda_{j-1},\dots,\lambda_1) - \lambda_j H_{i+j-1}(\mu_1,\dots,\mu_{i-1},\lambda_j,\dots,\lambda_1)}{\mu_i-\lambda_j},\\ V(i,1) &=& H_j(\mu_1,\dots,\mu_{i-1},\mu_i), \\ W(1,j) &=& H_i(\lambda_j,\lambda_{j-1},\dots,\lambda_1), \\ \Psi(i,j) &=& H_{i+j}(\mu_1,\dots,\mu_i,\lambda_{j-1},\dots,\lambda_1). \end{array} \quad (14)$

Given a SISO bilinear system as in (11), consider the ordered tuples of left and right interpolation points: $[\{\mu_1\}, \{\mu_1,\mu_2\}],[\{\lambda_1\},\{\lambda_1,\lambda_2\}]$. The recovered matrices can be written in terms of data samples in the following way:

$\begin{array}{rcl} \mathbb{L} &=& \begin{bmatrix} \frac{H_1(\mu_1)-H_1(\lambda_1)}{\mu_1 - \lambda_1} & \frac{H_2(\mu_1,\lambda_1)-H_2(\lambda_2,\lambda_1)}{\mu_1 - \lambda_2} \\ \frac{H_2(\mu_1,\mu_2)-H_2(\mu_1,\lambda_1)}{\mu_2 - \lambda_1} & \frac{H_3(\mu_1,\mu_2,\lambda_1)-H_3(\mu_1,\lambda_2,\lambda_1)}{\mu_2 - \lambda_2} \end{bmatrix}, \\ \mathbb{L} &=& \begin{bmatrix} \frac{\mu_1 H_1(\mu_1)-\lambda_1 H_1(\lambda_1)}{\mu_1 - \lambda_1} & \frac{\mu_1 H_2(\mu_1,\lambda_1)-\lambda_2 H_2(\lambda_2,\lambda_1)}{\mu_1 - \lambda_2} \\ \frac{\mu_2 H_2(\mu_1,\mu_2)- \lambda_1 H_2(\mu_1,\lambda_1)}{\mu_2 - \lambda_1} & \frac{\mu_2 H_3(\mu_1,\mu_2,\lambda_1)- \lambda_2 H_3(\mu_1,\lambda_2,\lambda_1)}{\mu_2 - \lambda_2} \end{bmatrix}, \\ \Psi &=& \begin{bmatrix} H_2(\mu_1,\mu_2) & H_3(\mu_1,\lambda_2,\lambda_1) \\ H_3(\mu_1,\mu_2,\lambda_1) & H_4(\mu_1,\mu_2,\lambda_2,\lambda_1) \end{bmatrix}, \\ V &=& \begin{bmatrix} H_1(\mu_1) \\ H_2(\mu_1,\mu_2) \end{bmatrix} = \mathcal{O}B, \\ W &=& \begin{bmatrix} H_1(\lambda_1) & H_2(\lambda_2,\lambda_1) \end{bmatrix} = C\mathcal{R}. \end{array}$

This reduced system matches eight moments of the original system, namely:

$\begin{array}{rcl} \text{two of} \;\; H_1&:& H_1(\mu_1), H_1(\lambda_1) \\ \text{three of} \;\; H_2 &:& H_2(\mu_1,\mu_2), H_2(\mu_1,\lambda_1), H_2(\lambda_2,\lambda_1), \\ \text{two of} \;\; H_3 &:& H_3(\mu_1,\mu_2,\lambda_1), H_3(\mu_1,\lambda_2,\lambda_1) \\ \text{one of} \;\; H_4&:& H_4(\mu_1,\mu_2,\lambda_2,\lambda_1), \end{array}$

i.e. in total $2k+k^2 = 8$, for $k=2$, moments are matched using this procedure.

Loewner Framework for Quadratic-Bilinear Systems

Consider quadratic-bilinear systems $\Sigma_{QB}$ characterized by the following equations:

$\Sigma_{QB} = \begin{cases} E\dot{x}(t) &= Ax(t) + Q(x(t) \otimes x(t)) + Nx(t)u(t) + Bu(t), \\ y(t) &= Cx(t), \end{cases} \quad (15)$

where $Q \in \mathbb{R}^{n \times n^2}$ and the other quantities are exactly as in (1) and (11) for the SISO case ($m=p=1$). The multi-input case is technically more involved but it is based on the same ideas. It follows that the differential equation in (15) can be split into an infinite series of equations:

$\begin{cases} E\dot{x}_1(t) &= Ax_1(t) + Bu(t), \quad\ell = 1, \\ E\dot{x}_2(t) &= Ax_2(t) + Q(x)_1(t) \otimes x_1(t)) + Nx_1(t)u(t), \quad\ell = 2, \\ &\vdots \\ E\dot{x}_1(t) &= Ax_\ell(t) + Q(\sum_{k=1}^{\ell-1}x(t) \otimes x(t)) + Nx_{\ell - 1}(t)u(t), \quad\ell \geq 2, \\ &\vdots \end{cases} \quad (16)$

Assume zero initial conditions and integrate both sides of the differential equations in (16). By using the second equation in (15), we can find a closed form relation between the output $y(t)$ and the input $u(t)$ in terms of infinitely many multivariate time domain kernels. Again, by transforming these functions to the frequency domain, we are able to construct rational multivariate transfer functions. Some of these will be used in the proposed procedure. More precisely, the one corresponding to level zero is the linear transfer function $H_0^\epsilon(s) = C(sE-A)^{-1} B$, which can also be denoted with $H(s)$. Furthermore, the two functions corresponding to level one are written as $H_1^N(s_0,s_1) = C(s_0E-A)^{-1} N (S_1E-A)^{-1} B$ and

$H_1^Q(s_0,s_1,s_2) = C(s_0E-A)^{-1} Q ((s_1E-A)^{-1} B \otimes (s_2 E - A)^{-1} B),$

For the second level, proceed with the following example

$H_2^{N,Q}(s_0,s_1,s_2,s_3) = C(s_0E-A)^{-1} N (s_1 E - A)^{-1} Q ((s_2E-A)^{-1} B \otimes (s_3 E - A)^{-1} B),$

In total, we can write $2^\ell$ functions for each level $\ell \geq 0$, as explicitly presented in . Given a SISO quadratic-bilinear system $\Sigma_{QB}$ as in (15), choose $6$ interpolation points $\{\omega_1,\omega_2,\dots,\omega_6\}$. First partition this set into two disjoint sets corresponding to right and left points, as $\{\lambda_1,\lambda_2,\lambda_3\} \cup \{\mu_1,\mu_2,\mu_3\}$. Next, consider the multi-tuples $\boldsymbol{\lambda}$ and $\boldsymbol{\mu}$, composed of the tuples

$\begin{array}{rcl} \boldsymbol{\lambda} &=& \{ (\lambda_1), (\lambda_2,\lambda_1), (\lambda_3,\lambda_1,\lambda_1) \}, \\ \boldsymbol{\mu} &=& \{ (\mu_1), (\mu_1,\mu_2), (\mu_1,\lambda_1,\mu_3) \}. \end{array}$

The Loewner matrix $\mathbb{L} \in \mathbb{R}^{3 \times 3}$ is a divided difference matrix, which is written as,

$\mathbb{L} = \begin{bmatrix} \frac{H(\mu_1)-H(\lambda_1)}{\mu_1-\lambda_1} & \dots & \frac{H_1^Q(\mu_1,\lambda_1,\lambda_1)-H_1^Q(\lambda_3,\lambda_1,\lambda_1)}{\mu_1-\lambda_3} \\ \frac{H_1^N(\mu_1,\mu_2)-H_1^N(\mu_1,\lambda_1)}{\mu_2-\lambda_1} & \ddots & \frac{H_1^{N,Q}(\mu_1,\mu_2,\lambda_1,\lambda_1)-H_1^{N,Q}(\mu_1,\lambda_3,\lambda_1,\lambda_1)}{\mu_2-\lambda_3} \\ \frac{H_1^Q(\mu_1,\lambda_1,\mu_3)-H_1^Q(\mu_1,\lambda_1,\lambda_1)}{\mu_3-\lambda_1} & \dots & \frac{H_2^{Q,Q}(\mu_1,\lambda_1,\mu_3,\lambda_1,\lambda_1)-H_2^{Q,Q}(\mu_1,\lambda_1,\lambda_3,\lambda_1,\lambda_1)}{\mu_3-\lambda_3} \\ \end{bmatrix}$

Similarly, the shifted Loewner matrix $\mathbb{L} \in \mathbb{R}^{3 \times 3}$ is computed; also the vectors $V$ and $W$

$\begin{array}{rcl} V &=& \begin{bmatrix} H(\mu_1 \\ H_1^N(\mu_1,\mu_2) \\ H_1^Q(\mu_1,\lambda_1,\mu_3) \end{bmatrix} \\ W &=& \begin{bmatrix} H(\lambda_1) & H_1^N(\lambda_2,\lambda_1) & H_1^Q(\lambda_3,\lambda_1,\lambda_1) \end{bmatrix} \end{array}$

Next, construct the matrices $\Psi$ and $\Omega$, corresponding to the bilinear and quadratic dynamics,

$\begin{array}{rcl} \Psi &=& \begin{bmatrix} H_1^N(\mu_1,\lambda_1) & H_2^{N,N}(\mu_1,\lambda_2,\lambda_1) & H_2^{N,Q}(\mu_1,\lambda_3,\lambda_2,\lambda_1) \\ H_2^{N,N}(\mu_1,\mu_2,\lambda_1) & H_3^{N,N,N}(\mu_1,\mu_2,\lambda_2,\lambda_1) & H_3^{N,N,Q}(\mu_1,\mu_2,\lambda_3,\lambda_2,\lambda_1) \\ H_2^{Q,N}(\mu_1,\lambda_1,\mu_3,\lambda_1) & H_3^{Q,Q}(\mu_1,\lambda_1,\mu_3,\lambda_2,\lambda_1) & H_3^{Q,N,Q}(\mu_1.\lambda_1,\mu_2,\lambda_3,\lambda_2,\lambda_1) \end{bmatrix}, \\ \Omega &=& \begin{bmatrix} H_1^Q(\mu_1,\lambda_1,\lambda_1) & H_2^{Q,N}(\mu_1,\lambda_1,\lambda_2,\lambda_1) & H_2^{Q,Q}(\mu_1,\lambda_1,\lambda_3,\lambda_1,\lambda_1) & \dots \\ H_2^{N,Q}(\mu_1,\mu_2,\lambda_1,\lambda_1) & H_3^{N,Q,N}(\mu_1,\mu_2,\lambda_1,\lambda_2,\lambda_1) & H_3^{N,Q,Q}(\mu_1,\mu_2,\lambda_1,\lambda_3,\lambda_1,\lambda_1) & \ddots \\ H_2^{Q,Q}(\mu_1,\lambda_1,\mu_3,\lambda_1,\lambda_1) & H_3^{Q,Q,N}(\mu_1,\lambda_1,\mu_3,\lambda_1,\lambda_2,\lambda_1) & H_3^{Q,Q,Q}(\mu_1,\lambda_1,\mu_3,\dots,\lambda_3,\lambda_1,\lambda_1) & \dots \end{bmatrix} \end{array}$

The matrix $\Omega \in \mathbb{R}^{3 \times 9}$ represents the reduced $Q$ matrix. It follows that given the system $\Sigma_{QB}$ a reduced system $\hat{\Sigma}_{QB} := (-\mathbb{L},V,W,0,-\mathbb{L}_s,\Psi,\Omega)$ of order $k=3$ can be constructed. This reduced system matches at most $2k + k^2 + k^3 = 42$ moments of the original system, i.e.

$\begin{array}{rcl} \text{two of} \;\; H_0^\epsilon &:& H(\mu_1), H(\lambda_1), \\ \text{three of} \;\; H_1^N &:& H_1(\mu_1,\mu_2), H_1^N(\mu_1,\lambda_1), H_1^N(\lambda_2,\lambda_1), \\ \text{three of} \;\; H_1^Q &:& H_1^Q(\mu_1,\lambda_1,\mu_3), H_1^Q(\mu_1,\mu_2,\lambda_1), H_1^Q(\mu_1,\lambda_2,\lambda_1), \\ \text{two of} \;\; H_2^{N,Q} &:& H_2^{N,Q}(\mu_1,\lambda_3,\lambda_2,\lambda_1), H_2^{N,Q}(\mu_1,\mu_2,\lambda_1,\lambda_1). \end{array}$

Loewner Framework for Linear Switched Systems

A continuous time linear switched system (LSS) is a dynamical system described by equations

$\Sigma_{\text{LSS}} : \begin{cases} E_{\sigma(t)} \dot{x}(t) &= A_{\sigma(t)}x(t) + B_{\sigma(t)}u(t), \quad x(0) = x_0, \\ y(t) &= C_{\sigma(t)} x(t) \end{cases}, \qquad (17)$

where $\sigma(t)$ is the switching signal ($\sigma:\mathbb{R}\to\mathcal{M}$), $u$ is the input, $x$ is the state, $y$ is the output, and $\mathcal{M} = \{1,2,\dots,P\}$, $P>1$, is the set of discrete modes. The system matrices $E_q, A_q \in \mathbb{R}^{n_q \times n_q}$, $B_q \in \mathbb{R}^{n_q \times m}$, $C_q \in \mathbb{R}^{p \times n_q}$, where $q \in \mathcal{M}$, correspond to the linear system active in mode $q \in \mathcal{M}$, and $x_0 \in \mathbb{R}^{n_{q_s}}$. Here, $n_1,n_2,\dots,n_P,m$ and $p$ are positive integers and $q_s \in \mathcal{M}$ is the mode in which the system is initialized. We consider the $E_q$ matrices to be invertible and $m=p=1$. Furthermore, the transition from one mode to another is made via the so called switching or coupling matrices $K_{q,\tilde{q}} \in \mathbb{R}^{n_\tilde{q} \times n_q}$, where $q,\tilde{q} \in \mathcal{M}$, as $E_{q_{i+1}} x(T_i) = K_{q,\tilde{q}} \lim_{t \nearrow T_i} x(t)$.

The restriction of the switching signal $\sigma(t)$ to a finite interval of time $[0,T]$ is given as

$\sigma(t) = \begin{cases} q_1 &\text{if} \;\; t \in [0,T_1), \\ q_i &\text{if} \;\; t \in [T_{i-1},T_i), \;\; i>2. \end{cases}$

where $q_1,\dots,q_k \in \mathcal{M}$ and $0<t_1<t_2<\dots<t_k \in \mathbb{R}_+$, $T_i := t_1 + \dots + t_{i-1}+t_i$, $T_k := T$. The input-output mapping is characterized in frequency domain by a series of multivariate rational functions obtained by computing the Laplace transform of generalized time kernels:

$\begin{array}{rcl} H_{q_1}(s_1) &=& C_{q_1}(s_1E_1-A_1)^{-1} B_{q_1} \\ H_{q_1,q_2}(s_1,s_2) &=& C_{q_1}(s_1E_1-A_1)^{-1} K_{q_1,q_2} (s_2E_2-A_2)^{-1} B_{q_2} \\ &\dots& \end{array}$

In general, we can write the level $k$ generalized transfer function associated to the switching sequence $q_1,q_2,\dots,q_k$, and evaluated at the points $s_1,s_2,\dots,s_k)$ in the complex plane, as

$H_{q_1,q_2,\dots q_k}(s_1,s_2,\dots,s_k) = C_{q_1}(s_1E_1-A_1)^{-1} K_{q_1,q_2} (s_2E_2-A_2)^{-1} \dots K_{q_{k-1},q_k} (s_kE_k-A_k)^{-1} B_{q_k}.$

By using samples of such functions, we can directly construct (reduced) switched models that interpolate the original model, i.e. the extension of the Loewner framework to LSS. For example, given a system $\Sigma_{\text{LSS}}$ as in (17) with $P=2$ modes, consider the ordered tuples of left interpolation points: $\{ (\mu_1), (\mu_2,\mu_3)\}$, $\{ (\mu_2), (\mu_1,\mu_4)\}$ and right interpolation points $\{ (\lambda_1), (\lambda_3,\lambda_2)\}$, $\{ (\lambda_2), (\lambda_4,\lambda_1)\}$. The projected Loewner matrices can be written as

$\begin{array}{rcl} \mathbb{L}_1 &=& \begin{bmatrix} \frac{H_1(\mu_1)-H_1(\lambda_1)}{\mu_1-\lambda_1} & \frac{H_{1,2}(\mu_1,\lambda_2)-H_{1,2}(\lambda_3,\lambda_2)}{\mu_1-\lambda_3} \\ \frac{H_{2,1}(\mu_2,\mu_3)-H_{2,1}(\mu_2,\lambda_1)}{\mu_3-\lambda_1} & \frac{H_{2,1,2}(\mu_2,\mu_3,\lambda_2)-H_{2,1,2}(\mu_2,\lambda_3,\lambda_2)}{\mu_3-\lambda_3} \end{bmatrix}, \\ \mathbb{L}_2 &=& \begin{bmatrix} \frac{H_2(\mu_2)-H_2(\lambda_2)}{\mu_2-\lambda_2} & \frac{H_{2,1}(\mu_2,\lambda_1)-H_{2,1}(\lambda_4,\lambda_1)}{\mu_2-\lambda_4} \\ \frac{H_{1,2}(\mu_1,\mu_4)-H_{1,2}(\mu_1,\lambda_2)}{\mu_4-\lambda_2} & \frac{H_{1,2,1}(\mu_1,\mu_4,\lambda_4)-H_{1,2,1}(\mu_1,\lambda_4,\lambda_1)}{\mu_4-\lambda_4} \end{bmatrix}. \end{array}$

Then explicitly write the matrices $V_i$, $W_j$ and $\Xi_j$ matrices (for $i,j \in \{1,2\}$) as

$\begin{array}{rcl} V_1 &=& \begin{bmatrix} H_1(\mu_1) \\ H_{2,1}(\mu_2,\mu_3) \end{bmatrix}, \\ V_2 &=& \begin{bmatrix} H_2(\mu_2) \\ H_{1,2}(\mu_1,\mu_4) \end{bmatrix}, \\ W_1 &=& \begin{bmatrix} H_1(\lambda_1) & H_{1,2}(\lambda_3,\lambda_2) \end{bmatrix}, \\ W_2 &=& \begin{bmatrix} H_2(\lambda_2) & H_{2,1}(\lambda_4,\lambda_1) \end{bmatrix}, \\ \Xi_{1,2} &=& \begin{bmatrix} H_{2,1}(\mu_2,\lambda_1) & H_{2,1,2}(\mu_2,\lambda_3,\lambda_2) \\ H_{1,2,1}(\mu_1,\mu_4,\lambda_1) & H_{1,2,1,2}(\mu_1,\mu_4,\lambda_3,\lambda_2) \end{bmatrix}, \\ \Xi_{2,1} &=& \begin{bmatrix} H_{1,2}(\mu_1,\lambda_2) & H_{1,2,1}(\mu_1,\lambda_4,\lambda_1) \\ H_{2,1,2}(\mu_2,\mu_3,\lambda_2) & H_{2,1,2,1}(\mu_2,\mu_3,\lambda_4,\lambda_1) \end{bmatrix}, \end{array}$

where $\Xi_{i,j}$ represent the reduced coupling matrices $K_{i,j}$. It follows that a reduced LSS of order $k=2$, denoted with $\hat{\Sigma}_{LSS} = (-\mathbb{L}_i,V_i,W_i,0,-\mathbb{L}_{si},\Xi_{i,j})$, $i,j \in \{1,2\}$, can be obtained without computation matrix factorizations or solves. This reduced system matches $2(2k + k^2) = 16$ moments of the original system, namely:

$\begin{array}{rcl} \text{four of} \;\; H_1/H_2 &:& H_1(\mu_1), H_2(\mu_2), H_1(\lambda_1), H_2(\lambda_2) \\ \text{three of} \;\; H_{1,2} &:& H_{1,2}(\mu_1,\mu_4), H_{1,2}(\mu_1,\lambda_2), H_{1,2}(\lambda_3,\lambda_2) \\ &\dots& \\ \text{one of} \;\; H_{1,2,1,2} &:& H_{1,2,1,2}(\mu_1,\mu_4,\lambda_3,\lambda_2) \\ \text{one of} \;\; H_{2,1,2,1} &:& H_{2,1,2,1}(\mu_2,\mu_3,\lambda_4,\lambda_1). \end{array}$

Summary

In this note we have addressed the Loewner framework for model reduction of some classes of (nonlinear) dynamical systems. The underlying philosophy valid for each of the extensions is as follows: One needs to collect data, i.e. values of high-order transfer functions (either by performing measurements or by direct numerical simulations) and then extract the desired information directly from the model constructed. The strength of the proposed approaches consists in the fact that, given input-output data, one can construct with basically no computational cost, a model of the data in generalized state-space format.

Further Variants

Additionally, there are further variants of the Loewner framework for parametric systems , , singular systems and based on time-domain data .