Time-limited $\mathcal{H}_{2}$-optimal Model Order Reduction of Linear Systems with Quadratic Outputs

This research work studies a special class of nonlinear dynamical systems characterized by weak nonlinearity. These systems retain linear time-invariant (LTI) state equations but incorporate quadratic nonlinear terms in their output equations, thus named linear quadratic output (LQO) systems [1]. They naturally arise in scenarios requiring the observation of quantities that involve the product of state components, either in the time or frequency domain. These systems are valuable in quantifying energy or power of dynamical systems, for example, in assessing a system’s internal energy [2] or the cost function in optimal quadratic control problems [3]. Moreover, they are used to measure deviations of state coordinates from a reference point, such as for calculating the root mean squared displacement of spatial coordinates around an excitation point or for estimating the variance of a random variable in stochastic modeling [4].

Consider an LQO system defined by the following state and output equations:

where $x(t)\in\mathbb{R}^{N\times 1}$ denotes the state vector, $u(t)\in\mathbb{R}^{m\times 1}$ represents inputs, $y(t)\in\mathbb{R}^{p\times 1}$ is the output vector. The matrices $A\in\mathbb{R}^{N\times N}$, $B\in\mathbb{R}^{N\times m}$, $C\in\mathbb{R}^{p\times N}$, and $M_{i}\in\mathbb{R}^{N\times N}$ consititute the state-space realization $(A,B,C,M_{1},\cdots,M_{p})$. The state equation is identical in form to that of a standard LTI system. However, the output equation introduces a nonlinearity through the quadratic terms $x(t)^{T}M_{i}x(t)$, which distinguishes the LQO system from a standard LTI system.

Accurately modeling complex physical phenomena frequently necessitates dynamical systems of exceptionally high order, often in the range of several thousand or more. Due to this very high order $N$, simulating and analyzing the model (1) becomes computationally intensive and impractical. Consequently, it becomes necessary to approximate (1) with a reduced-order model (ROM) of significantly lower order $n$ (where $n\ll N)$, which simplifies simulation and analysis [5]. Model order reduction (MOR) refers to the process of constructing a ROM while ensuring that the essential features and characteristics of the original model are preserved [6].

We define the $n^{th}$-order ROM of the system $H$ as $\hat{H}$, which is described by the state and output equations shown in (2)

Here, $\hat{A}$, $\hat{B}$, $\hat{C}$, and $\hat{M}_{i}$ are matrices derived from the original system matrices $A$, $B$, $C$, and $M_{i}$, respectively, through the Petrov-Galerkin projection condition $\hat{W}^{T}\hat{V}=I$. Specifically, $\hat{A}=\hat{W}^{T}A\hat{V}$, $\hat{B}=\hat{W}^{T}B$, $\hat{C}=C\hat{V}$, and $\hat{M_{i}}=\hat{V}^{T}M_{i}\hat{V}$. The projection matrices $\hat{V}\in\mathbb{R}^{N\times n}$ and $\hat{W}\in\mathbb{R}^{N\times n}$ are used to project the original system $H$ onto a reduced subspace, resulting in the ROM $\hat{H}$. Different MOR techniques vary in their approach to constructing $\hat{V}$ and $\hat{W}$, which is dependent on the specific characteristics of $H$ that need to be preserved in $\hat{H}$ [7, 8, 9]. Throughout this paper, it is assumed that both $A$ and $\hat{A}$ are Hurwitz matrices.

The Balanced Truncation (BT) method, developed in 1981, is a widely utilized MOR technique [10]. This method selectively retains states that significantly contribute to the energy transfer between inputs and outputs, while discarding those with minimal influence, as determined by their Hankel singular values. One key advantage of BT is its capability to estimate errors a priori before constructing the ROM [11]. Furthermore, BT ensures that the stability of the original model is maintained. Initially developed for standard LTI systems, BT has significantly expanded its applicability to encompass various system types, including descriptor systems [12], second-order systems [13], linear time-varying systems [14], parametric systems [15], nonlinear systems [16], and bilinear systems [17], among others. Additionally, BT has been adapted to preserve specific system properties, such as positive realness [18], bounded realness [19], passivity [20], and special structural characteristics [21]. For a comprehensive understanding of the diverse BT algorithm family, refer to the survey [22]. BT has been extended to LQO systems in [23, 24, 25]. Among these algorithms, only the one presented in [25] maintains the LQO structure in the ROM.

A locally optimal solution in the $\mathcal{H}_{2}$ norm is achievable using several efficient algorithms. The necessary conditions for this local optimum, known as Wilson’s conditions [26], involve interpolation of the original system at selected points. The iterative rational Krylov algorithm (IRKA) is a well-known algorithm that can achieve this local optimum [27]. A more general algorithm, which does not assume the original system having simple poles, is also available [28]. This algorithm is based on Sylvester equations and its numerical properties have been improved in [29]. Recently, the $\mathcal{H}_{2}$ MOR problem for LQO systems has been addressed, and an algorithm based on Sylvester equations has been proposed as a solution [30].

BT aims to approximate the original model dynamics over the entire time horizon. However, practical systems and simulations often operate within finite time intervals due to real-world constraints. For example, in interconnected power systems, low-frequency oscillations typically last around 15 seconds and are mitigated by power system stabilizers and damping controllers [31]. The first $15$ seconds are crucial for small-signal stability assessments [32]. Similarly, in finite-time optimal control problems, the system’s behavior within a specific time frame is of utmost importance [33]. This need has led to the development of time-limited MOR, which focuses on achieving good accuracy within a specified time interval rather than the entire time domain [34].

The time-limited MOR problem focuses on developing a ROM that guarantees a small deviation between the outputs of the original model $y(t)$ and the ROM’s output $y_{n}(t)$ for a given input signal $u(t)$, but only within a specified limited time interval $[0,\tau]$ seconds. In other words, the goal is to ensure that the approximation error $\|y(t)-y_{n}(t)\|$ remains minimal within this specified time frame [35].

To address the time-limited MOR problem, BT has been adapted into the time-limited BT (TLBT) algorithm [34]. Although TLBT does not fully retain all characteristics of traditional BT, such as stability guarantees and a priori error bounds, it effectively addresses the time-limited MOR scenario. The computational aspects of TLBT and strategies for handling large-scale systems have been explored in [36]. Furthermore, TLBT has been extended to handle descriptor systems [36], second-order systems [37], and bilinear systems [38], enhancing its applicability. Recently, TLBT has also been extended to LQO systems in [39].

In [40], a new norm called the time-limited $\mathcal{H}_{2}$ norm is defined, and the conditions needed to find a local optimum for this norm are derived. An algorithm called time-limited IRKA (TLIRKA) is then proposed based on IRKA to construct a local optimum, but it does not meet any necessary conditions. In [41], new optimality conditions based on interpolation are derived, and a nonlinear optimization algorithm is proposed to construct the local optimum, whose applicability is limited to medium-scale systems. This paper studies the time-limited $\mathcal{H}_{2}$-optimal MOR problem for LQO systems.

The key contributions of this research work are manifold. First, it defines the time-limited $\mathcal{H}_{2}$ norm ($\mathcal{H}_{2,\tau}$ norm) for LQO systems and demonstrates its computation using time-limited system Gramians defined in [39]. Second, it derives the necessary condition for achieving a local optimum of $||H-\hat{H}||_{\mathcal{H}_{2,\tau}}^{2}$. Third, these conditions are then contrasted with those related to the standard $\mathcal{H}_{2}$-optimal MOR [30]. Notably, it is shown that Petrov-Galerkin projection cannot, in general, attain a local optimum in the time-limited setting. Fourth, a stationary point algorithm, rooted in the Petrov-Galerkin projection, is proposed. Upon convergence, this algorithm satisfies three of the four necessary conditions for optimality. An illustrative numerical example is presented to demonstrate the accuracy of the proposed algorithm within the specified time interval.

In this section, we will briefly discuss two of the most relevant MOR algorithms for LQO systems in the context of the problem under consideration. The first algorithm is the $\mathcal{H}_{2}$-optimal MOR method [30], and the second is the TLBT [39].

This section introduces the time-limited $\mathcal{H}_{2}$ norm and its connection to the time-limited observability Gramian. Subsequently, necessary conditions for the local optimality of the squared time-limited $\mathcal{H}_{2}$ norm of the error are derived. Based on these optimality conditions, a Petrov-Galerkin projection-based iterative algorithm is proposed that satisfies three of the four optimality conditions. The challenges associated with fulfilling the remaining optimality condition within the Petrov-Galerkin projection framework are also discussed. Finally, the computational efficiency of the proposed algorithm is briefly analyzed.