Structure-Preserving Model Reduction for Nonlinear Power Grid Network
^††thanks: This work was supported in parts by National Science Foundation under Grant No. DMS-1923221

Power networks are complex and large-scale systems in operation. Simulation of these high dimensional power system models are computationally expensive and demand unmanageable levels of storage in dynamic simulation or trajectory sensitivity analysis. Hence, reducing the order of these models is of a great importance especially for real time applications, stability analysis, and control design [12].

Numerous model order reduction (MOR) techniques have been applied to the linearized power networks models including balanced truncation [20, 11, 32, 36, 18], Krylov subspace (interpolatory) methods [8, 38, 39, 35], proper orthogonal decomposition (POD) [41], singular perturbation theory [27, 13], sparse representation of the system matrices [19], and clustering analysis [10, 37, 16, 23]. For a comparative study of MOR techniques in power systems, see, e.g., [7, 1, 12]. Recently, development of MOR techniques for nonlinear systems including power networks has received great attention. In [28], POD is employed to reduce the hybrid, nonlinear model of a power network for the study of cascading failures. Trajectory piece-wise linearization (TPWL) method to approximate the nonlinear term in swing models together with POD to find a reduced order model for the swing dynamics was employed in [21]. Balanced truncation based on empirical controllability and observability covariances of nonlinear power systems is proposed in [30]. The empirical Gramians and balanced realization were also used in [42, 43] to build a nonlinear power system model suitable for controller design. In [29], real-time phasor measurements are used to cluster generators based on similar behaviors using recursive spectral bipartitioning and spectral clustering. The reduced model is formed by representative generators chosen from each cluster. A model reduction approach that adaptively switches between a hybrid system with selective linearization and a linear system based on the variations of the system state or size of a disturbance was developed in [26]. A structure-preserving, aggregation-based model order reduction framework based on synchronization properties of power systems is introduced in [24]. In [40] an SVD-based approach for real-time dynamic model reduction with the ability of preserving some nonlinearity in the reduced model is presented. The recent work [33] presents the idea of transforming the swing equations into a quadratic system, lifts nonlinear swing equations to a quadratic one, and employs quadratic balanced truncation model order reduction technique. A non-intrusive data-driven modeling framework for power network dynamics using Lift and Learn [31] has been studied in [34].

In this paper, we focus on designing a structure-preserving and input-independent MOR framework for a special class of power grid networks. We obtain the reduction bases using the idea of lifting the nonlinear dynamics to quadratic models and employing $\mathcal{H}_{2}$-quasi-optimal model reduction for quadratic systems (Q-IRKA) [6]. The reduction basis for the second-order nonlinear power network grids model is obtained with the help of the reduction bases obtained by Q-IRKA and their particular subspace structure. This reduction basis is then applied on the second-order nonlinear power network grids model to form a structure-preserving reduced model.

The rest of the paper is organized as follows: Section II briefly introduces the nonlinear dynamical model of the power grid network followed by a quick recap of projection-based structure-preserving model reduction in this setting in Section III. Section IV illustrates how to lift the second-order dynamics of the power grid network to quadratic dynamics and presents an equilibrium analysis. In Section V, we describe the model reduction for quadratic systems via the $\mathcal{H}_{2}$-based model reduction, Q-IRKA, followed by the required modification for the quadratized power grid model before employing $\mathcal{H}_{2}$-based model reduction. In Section VI, we present our structured-preserving model reduction approach for power grid networks by proving and exploiting the specific structure arising in Q-IRKA for these dynamics. Section VII illustrates the feasibility of our approach via numerical examples and by comparing our approach with two other methods, followed by conclusions in Section VIII.

There are three most commonly used models for describing a dynamical model of a power grid network: synchronous motor ($SM$), effective network ($EN$) and structure-preserving ($SP$). Each model is described as a network of $n$ coupled oscillators (generators or loads) whose dynamics at the $i$th node is expressed by the nonlinear second-order equation [25]

where $\delta_{i}$ is the angle of rotation for the $i$th oscillator; $J_{i}$ and $D_{i}$ are inertia and damping constants, respectively; $\omega_{R}$ is the angular frequency for the system; $K_{ij}\geq 0$ is dynamical coupling between oscillator $i$ and $j$; and $\gamma_{ij}$ is the phase shift in this coupling. Constants $B_{i}$, $K_{ij}$ and $\gamma_{ij}$ are computed by solving the power flow equations and applying Kron reduction [25]. The three models mainly differ by the way they model loads. The EN model assumes loads as constant impedances instead of oscillators. The SP model represents loads as first-order oscillators ($J_{i}=0$) while the SM model represent loads as synchronous motors expressed as second-order oscillators.

Define the state-vector $\delta=[\delta_{1}~{}\delta_{2}~{}\dots~{}\delta_{n}]^{T}\in\mathbb{R}^{n}$. Then, the dynamics of a network of $n$ coupled oscillators, given in (1) for the $i$th node, can be described in the system form

where the matrices $\mathcal{M}$, $\mathcal{D}{\in\mathbb{R}}^{n\times n}$, and $\mathcal{B}\in\mathbb{R}^{n}$ are defined as

where $m_{i}=\frac{2J_{i}}{\omega_{R}}$ and $d_{i}=\frac{D_{i}}{\omega_{R}}$, $u(t)=1$ is input to the dynamics; and the nonlinear map $f:\mathbb{R}^{n}\to\mathbb{R}^{n}$ is defined such that its $i$th component for $i=1,2,\ldots,n$ is

The variable $y(t)\in\mathbb{R}^{p}$ in (3) corresponds to the output (quantity of interest) of the network dynamics, described by the state-to-output matrix $\mathcal{C}\in\mathbb{R}^{p\times n}$.

Large scale power network dynamics, i.e., those with large $n$, are computationally expensive to simulate and control. With the aid of the model reduction, we are able to replace these high dimensional models by lower dimensional models that can closely mimic the input-output behavior of the original high dimensional systems. To be more specific, we seek to develop a nonlinear structure-preserving model reduction approach such that the reduced model preserves the physically-meaningful second-order structure. In other words, we directly reduce the second-order dynamics (2)-(3) to obtain a reduced second-order system with dimension $r\ll n$ and of the form

where $\delta_{r}\in\mathbb{R}^{r}$ is the reduced state; $\mathcal{M}_{r},\mathcal{D}_{r}\in\mathbb{R}^{r\times r}$; $\mathcal{B}_{r}\in\mathbb{R}^{r\times 1}$; $f_{r}:\mathbb{R}^{r}\to\mathbb{R}^{r}$; and $\mathcal{C}_{r}\in\mathbb{R}^{p\times r}$. The goal of this structure-preserving model reduction process is that $y_{r}(t)$ approximates the original output $y(t)$ with high fidelity.

Given the full-order dynamics (2)-(3), the most common approach for constructing the reduced model (8)-(9) is the Petrov-Galerkin projection framework: Construct two model reduction bases $\mathcal{W},\mathcal{V}\in\mathbb{R}^{n\times r}$ such that $\delta(t)\approx\mathcal{V}\delta_{r}(t)$. Insert $\mathcal{V}\delta_{r}(t)$ into (2) and enforce a Petrov-Galerkin condition on the residual to obtain the reduced matrices in (8) and (9) as

To preserve the symmetry and positive definiteness of $\mathcal{M}$ and $\mathcal{D}$ in (8), we employ a Galerkin projection by setting $\mathcal{W}=\mathcal{V}$. However, unlike the usual MOR via Galerkin projection where the subspace $\mathcal{V}$ (and thus $\mathcal{W}$) only contains information from the dynamics equation (2) and ignores the output equation (3), our one-sided Galerkin projection will contain information from both (2) and (3) by taking advantage of the special structure arising in the network dynamics. These issues are discussed in detail in Section VI. The term $f_{r}(\delta_{r})$ in (12) is usually further approximated to resolve the lifting bottleneck [9]. Due to the page limit, we skip this detail here.

The accuracy of the structure-preserving reduced model (8) depends on the proper choice of $\mathcal{V}$. For general nonlinearities, POD is usually the method of choice. However, for special classes of nonlinear systems such as bilinear and quadratic-bilinear systems, one can indeed use well-established system theoretical methods; see, e.g., [6, 5, 3]. A large class of nonlinear systems, such as those with smooth nonlinearities, e.g., exponential, trigonometric, can indeed be represented as quadratic-bilinear systems by introducing new variables (in a lifting map) [22], [14], [4], [33], [31]. Converting a nonlinear model to its equivalent quadratic form provides an opportunity to perform input-independent model reduction techniques where a system-theoretic norm is defined [5, 6]. Although this transformation is not unique, it is exact in many cases (as in here), i.e., there is no approximation error [5].

Second-order model (1) inherently contains a quadratic nonlinearity in $f(\delta)$ and the nonlinear system (2) can be lifted to a quadratic dynamic as we will see next in Section IV. However, as stated in section (III), we seek to obtain a reduced second-order model, not a reduced quadratic model. Therefore, by exploiting the inherent structure of the model reduction bases obtained for the quadratic model, we form the final reduction base $\mathcal{V}$ to construct the final structure preserving reduced model (12). This is achieved in Sections V–VI.

As we mentioned in Section III for systems with quadratic nonlinearities, effective and in some cases optimal, systems-theoretic MOR methods already exist. Therefore, in obtaining high-fidelity reduced models, we can significantly benefit from a quadratic representation of the nonlinear swing dynamics.

Towards this goal, recall the original second-order dynamics (1). Expand the nonlinearity $\sin(\delta_{i}-\delta_{j}-\gamma_{ij})$ using the trigonometric identity and introducing $s:=\sin$ and $c:=\cos$

Note that equation (13) is quadratic in $\sin(\delta_{i})$ and $\cos(\delta_{i})$. This directly hints at choosing $\sin(\delta)$ and $\cos(\delta)$ as auxiliary variables in a new state vector to convert the nonlinear dynamics to a quadratic one. This is precisely what has been done in [33] to find a quadratic representation for (2). However, since the transformation, the so-called lifting, is not unique and we have a slightly different form in the resulting state-space form, we include its derivation here for completeness.

Now that we have represented the original nonlinear swing dynamics as an equivalent quadratic dynamics, we will investigate the MOR problem for those systems in more details.

Consider the quadratic dynamical system

where $E,A\in\mathbb{R}^{N\times N}$, $H\in\mathbb{R}^{N\times N^{2}}$, $B\in\mathbb{R}^{N\times l}$, and $C\in\mathbb{R}^{p\times N}$. One can view $H$ as the mode-1 matricization of a third-order tensor $\mathcal{H}\in\mathbb{R}^{N\times N\times N}$. In other words, let $\mathcal{H}_{i}\in\mathbb{R}^{N\times N}$ for $i=1,\dots,N$ be the frontal slices of $\mathcal{H}$. Then,

where the notation ${\mathcal{H}^{(1)}}$ denotes the mode-1 matricization.

For the cases where the state-space dimension $N$ is large, the goal is to find a reduced quadratic system

where ${E}_{r},{{A}}_{r}\in\mathbb{R}^{r_{q}\times r_{q}}$ ${H}_{r}\in\mathbb{R}^{r_{q}\times r_{q}^{2}}$, ${C}_{r}\in\mathbb{R}^{p\times r_{q}}$ and ${{B}}_{r}\in\mathbb{R}^{r_{q}\times l}$ with $r_{q}\ll N$ such that $y(t)\approx y_{r}(t)$ for a wide range of inputs. We construct the reduce system (52) via projection: Construct two model reduction bases $V$, $W\in\mathbb{R}^{N\times r_{q}}$ such that the reduced matrices in (52) are given by

The quality of the reduced system (52) depends on the choice of the model reduction bases $V$ and $W$. There are various model reduction techniques to determine theses bases including trajectory-based method such as POD (see, e.g., [15]); and input-independent systems theoretic methods for quadratic systems such as balanced truncation [5], [33], interpolatory projections [4], [14] and $\mathcal{H}_{2}$-quasi-optimal model reduction [6]. In this paper, we employ the $\mathcal{H}_{2}$-based approach in [6].

Our aim is to perform structure-preserving model reduction of the second-order dynamics (2)-(3) using the reduction bases $\mathcal{V},\mathcal{W}\in\mathbb{R}^{n\times r}$ to obtain the reduced second-order model (8)-(12). To keep the symmetry of the underlying dynamics, we perform a Galerkin projection by choosing $\mathcal{W}=\mathcal{V}$. Converting the nonlinear second-order dynamics to the quadratic form has allowed us to use (optimal) systems-theoretic techniques to apply, such as Q-IRKA. Although we cannot simply use the output of Q-IRKA as the reduced model (since it does not preserve the second-order structure), we still would like to employ the high-fidelity model reduction bases $V$ and $W$ resulting from Q-IRKA in constructing $\mathcal{V}$. In this section, by analyzing and exploiting the structure of the underlying power network dynamics, we show how to efficiently construct $\mathcal{V}$ using $V$ and $W$ from Q-IRKA.