# On Scaling Robust Feedback Control and State Estimation Problems
in Power Networks***This work is supported by the National Science Foundation
under Grants ECCS 2151571 and CMMI 2152450.

MirSaleh Bahavarnia<EMAIL_ADDRESS>Muhammad Nadeem
<EMAIL_ADDRESS>Ahmad F. Taha<EMAIL_ADDRESS>
###### Abstract

Many mainstream robust control/estimation algorithms for power networks are
designed using the Lyapunov theory as it provides performance guarantees for
linear/nonlinear models of uncertain power networks but comes at the expense
of scalability and sensitivity. In particular, Lyapunov-based approaches rely
on forming semi-definite programs (SDPs) that are (i) not scalable and (ii)
extremely sensitive to the choice of the bounding scalar that ensures the
strict feasibility of the linear matrix inequalities (LMIs). This paper
addresses these two issues by employing a celebrated non-Lyapunov approach
(NLA) from the control theory literature. In lieu of linearized models of
power grids, we focus on (the more representative) nonlinear differential
algebraic equation (DAE) models and showcase the simplicity, scalability, and
parameter-resiliency of NLA. For some power systems, the approach is nearly
fifty times faster than solving SDPs via standard solvers with almost no
impact on the performance. The case studies also demonstrate that NLA can be
applied to more realistic scenarios in which (i) only partial state data is
available and (ii) sparsity structures are imposed on the feedback gain. The
paper also showcases that virtually no degradation in state estimation quality
is experienced when applying NLA.

###### keywords:

Robust control , power systems , differential algebraic equations ,
decentralized control , dynamic state estimation , linear matrix inequalities
, Lyapunov theory

[inst1]organization=Department of Civil and Environmental Engineering,
Vanderbilt University,addressline=2201 West End Avenue, city=Nashville,
postcode=37235, state=TN, country=USA

[inst2]organization=Department of Electrical and Computer Engineering,
Vanderbilt University,addressline=2201 West End Avenue, city=Nashville,
postcode=37235, state=TN, country=USA

##  1 Introduction

In current and future power systems, many feedback control and state
estimation algorithms will jointly be deployed to perform network-driven,
real-time feedback control. On one hand, feedback control algorithms aim to
robustly stabilize the frequency oscillations and improve the transient
stability in the presence of unknown uncertainties [1, 2, 3, 4]. On the other
hand, state estimation algorithms enable system operators to estimate the
unmeasurable dynamic and algebraic states in the presence of unknown
uncertainties and abrupt changes [5, 6, 7, 8].

Many studies including the aforementioned robust control/estimation
algorithms, have used Lyapunov theory, resulting in convex SDPs/LMI
formulations. Lyapunov-based approaches are significantly powerful and have
extensively been applied to design various robust controllers and estimators
for diverse engineering applications [9]. One of the main Lyapunov-based
approaches has been built upon a well-known lemma, namely, Kalman-Popov-
Yakubovich (KYP) Lemma (alternatively known as the Bounded Real Lemma) [9,
10], to cast $\mathcal{H}_{\infty}$ control/estimation problems as convex SDPs
via incorporating the LMIs. In power systems, such convex SDP/LMI formulations
have widely been utilized to design robust $\mathcal{H}_{\infty}$
control/estimation algorithms. The central idea behind the
$\mathcal{H}_{\infty}$ control/estimation [11] is that by minimizing the
$\mathcal{H}_{\infty}$ norm, one minimizes the effects of the unknown
uncertainties and abrupt changes on the desired performance of the
control/estimation algorithm.

In [12, 7, 13], casting the control/estimation problems for DAE-modeled power
systems as convex SDPs/LMI formulations via KYP Lemma, the authors propose
dense robust $\mathcal{H}_{\infty}$ controllers and estimators, respectively.
However, the Lyapunov-based control/estimation algorithms in the
aforementioned research works are highly costly in terms of computational
time. For instance, computing the Lyapunov-based dense $\mathcal{H}_{\infty}$
controller in [12] via MOSEK [14] SDP solver implemented in YALMIP [15] can
take more than $4$ hours (on a personal computer). In the presence of abrupt
model changes/uncertainties, it is not practically acceptable to utilize such
highly time-consuming $\mathcal{H}_{\infty}$ syntheses.

Although the Lyapunov-based approaches provide a solid theoretical framework
to analyze the control/estimation problems in terms of stability and
performance, similar to any other LMI-based approaches, they face two main
issues, namely, computational scalability and computational sensitivity.

The former issue occurs when a Lyapunov-based approach casts a
control/estimation problem as a convex SDP via incorporating the LMIs. It is
evident that dealing with large-scale LMIs and solving the corresponding SDPs
are generally challenging tasks in terms of scalability [16]. The authors in
[16] have surveyed some alternative scalable (yet possibly inaccurate) SDP
solvers like SDPNAL+ [17, 18] which is particularly developed based on a
majorized semi-smooth Newton-CG augmented Lagrangian method. Although SDPNAL+
has significantly improved the scalability in a class of robust traffic
density estimation problems [19], it has unfortunately been unable to even
propose a feasible solution for the $\mathcal{H}_{\infty}$ control/estimation
problem of a medium-size power system.

The latter issue occurs due to the fine-tuning of $\epsilon I$-shifting (where
$\epsilon$ is the bounding scalar that ensures the strict feasibility of the
LMIs and $I$ is the identity matrix) in LMIs as there is no systematic way in
the literature to appropriately choose/tune $\epsilon$ in general. Note that
the appropriate choice of $\epsilon$ positively affects the feasibility and
sub-optimality of the solution at the cost of excessive computational time.
The $\epsilon I$-shifting idea is commonly utilized to implement the strict
LMIs as the non-strict LMIs because SDP solvers cannot deal with the strict
LMIs.

One might think that utilizing the generalized algebraic Riccati equations
(GAREs) could be effective in overcoming the issues mentioned above.
Unfortunately, although casting the control/estimation problems for a DAE-
modeled power system as GAREs is theoretically possible (thus increasing the
scalability significantly) [20, 21, 22], solving the cast GAREs practically
faces numerical issues/difficulties due to checking a rank equality constraint
and an indefinite matrix coefficient appearing in an algebraic quadratic
equation (AQE). In fact, such an AQE is efficiently solvable if that matrix
coefficient is not indefinite which is not the case for a DAE-modeled power
system in general. Then, casting the $\mathcal{H}_{\infty}$ control/estimation
problems as GAREs to bypass the aforementioned scalability and sensitivity
issues in the case of DAE-modeled power systems, unfortunately, becomes
impractical.

To that end, in light of the aforementioned limitations and to effectively
overcome the arising issues by the Lyapunov-based approaches (i.e.,
computational scalability and computational sensitivity), considering the NDAE
state-space model of a power system, we derive the equivalent nonlinear
ordinary differential equation (ENODE) state-space model without imposing any
simplifying specific form on the controller or the estimator. It is noteworthy
that the continuous-time algebraic Riccati equation (CARE)-based approaches
are not applicable to the derived ENODE state-space model as the dependency of
the state-space matrices on the controller or the estimator is not affine
(e.g., $A+BK$ or $A+LC$) unlike the simplified state-space model derived in
[12]. To design $\mathcal{H}_{\infty}$ controllers and estimators for given
NDAE state-space models, motivated by the aforementioned twofold issues, we
employ a celebrated non-Lyapunov approach (NLA) from the control theory
literature [23] and utilize off-the-shelf computational tools along with the
derived ENODE state-space model to bypass the computational scalability and
computational sensitivity issues arising in Lyapunov-based
$\mathcal{H}_{\infty}$ controller and estimator designs for DAE-modeled power
systems. The proposed $\mathcal{H}_{\infty}$ synthesis (NLA) in [23] is mainly
built upon non-convex non-smooth optimization techniques taking advantage of
the Clarke sub-differential [24] computation of the $\mathcal{H}_{\infty}$
norm. Such a structured $\mathcal{H}_{\infty}$ synthesis has found extensive
control applications including but not limited to proportional-integral-
derivative (PID) control and combination of proportional-integral (PI) and SOF
control of distillation column, aircraft autopilot, and helicopter as
highlighted by [25, 26, 27].

Note that NLA can alternatively be applied to any $\mathcal{H}_{\infty}$
controller and estimator design problem that can potentially be cast as a
convex SDP/LMI formulation, e.g., [28, 13, 1, 7, 12, 29]. Moreover, compared
to such Lyapunov-based syntheses in the literature, we take into account a
more general problem setup consisting of structured controllers and estimators
instead of solely dense counterparts. In other words, thanks to the
specialized features of NLA, we are able to propose structured
$\mathcal{H}_{\infty}$ controllers and estimators.

The paper’s contributions are as follows:

  * 1.

We design structured $\mathcal{H}_{\infty}$ controllers and estimators for
DAE-modeled power systems via employing a celebrated NLA from the control
theory literature [23]. Specifically, for the controller design, the proposed
controller is applicable to a general scenario in which we deal with
partially-accessed states and noisy measurement outputs. Moreover, for the
estimator design, the proposed estimator can deal with load demand and
renewable disturbances, non-Gaussian measurement noise, and unknown control
inputs.

  * 2.

NLA is noticeably faster than the Lyapunov-based counterpart while attaining
the same $\mathcal{H}_{\infty}$ performance (computational scalability).
Moreover, there is no need to deal with the fine-tuning of the bounding scalar
that ensures the strict feasibility of the LMIs anymore as we have
successfully bypassed it (computational sensitivity).

  * 3.

We investigate the conditions that can potentially affect the superiority of
NLA over its Lyapunov-based counterpart in terms of computational scalability.
It turns out that having a relatively large number of inputs (for the
controller design) and outputs (for the estimator design) can negatively
affect such a superiority. Beyond a specific level, such superiority can even
be reversed.

  * 4.

Utilizing a notion of block-sparsity, we study the impact of various sparsity
structures on the $\mathcal{H}_{\infty}$ performance of the proposed
controller design. Also, we empirically explore the potential effects of
output matrix selection on the $\mathcal{H}_{\infty}$ performance and how we
can take advantage of such empirical information to select an optimal output
matrix. Specifically, we can detect the states that have a central role in the
stabilization of the system.

Notations. Uppercase and lowercase letters denote matrices and vectors,
respectively. We denote the set of real-valued $r$ by $c$ matrices by
$\mathbb{R}^{r\times c}$. Also, we denote the $r$ by $c$ identity and all-ones
matrices with $I_{r\times c}$ and $\mathbf{1}_{r\times c}$, respectively. To
represent an $n$-dimensional identity matrix, we utilize $I_{n}$. For a matrix
$M$, symbols $M^{T}$, $M^{\perp}$, $\bar{\sigma}[M]$, and $\|M\|_{0}$ denote
its transpose, orthogonal complement, largest singular value, and number of
non-zero elements, respectively. We denote positive semi-definiteness and
positive/negative definiteness by $\succeq 0$ and $\succ 0$/$\prec 0$,
respectively. For a square matrix $M$, we define $\mathbf{He}(M):=M+M^{T}$,
and symbol $\mathbf{sa}(M)$ denotes the spectral abscissa of $M$, i.e., the
maximum real part of the eigenvalues of $M$. Symbol $\ast$ represents the
corresponding symmetric components of a symmetric matrix. We use $j$ and
$\Re(s)$ to showcase the $\sqrt{-1}$ and the real part of complex-valued
scalar $s$, respectively. We denote the matrix element-wise Hadamard product
and the matrix Kronecker product with $\odot$ and $\otimes$, respectively. For
a vector $v$, symbol $\|v\|$ denotes the Euclidean norm of $v$. Given a vector
$h(t)$ defined for $t\in[0,\infty)$, its $\mathcal{L}_{2}$-norm is defined as
$\sqrt{\int_{0}^{\infty}\|h(t)\|^{2}dt}$. To save space, we omit the
dependency on time, i.e., $(t)$ in denoting the time-dependent quantities. The
set of subgradients of a function $f$ at the point $x$ is defined as the
subdifferential of $f$ at $x$ and is denoted by $\partial f(x)$. We represent
the function composition operator by $\circ$. The standard basis vector
$\mathbf{e}_{i}$ represents a column vector with $1$ at the $i^{th}$ position
and $0$ at the rest of the positions.

The rest of the paper is organized as follows: Section 2 consisting of
preliminaries are three-fold: (i) Section 2.1 details the NDAE state-space
model of the power system, (ii) Section 2.2 elaborates on the problem
statement, and (iii) Section 2.3 briefly reviews dense $\mathcal{H}_{\infty}$
controller and estimator designs built upon the Lyapunov-based approach for
DAE-modeled power systems. Section 3 contains the central result of the paper
which is the implementation of NLA to design structured $\mathcal{H}_{\infty}$
controllers and estimators for DAE-modeled power systems. Section 4 via
various benchmarks thoroughly verifies the effectiveness of NLA in comparison
with its Lyapunov-based counterpart. Finally, Section 5 ends the paper with a
few concluding remarks.

##  2 Preliminaries

In Section 2.1, we include the details of the NDAE state-space model of the
power system. We formally state the structured $\mathcal{H}_{\infty}$
controller problem for DAE-modeled power systems in Section 2.2. The
$\mathcal{H}_{\infty}$ estimator problem can also be stated in a similar
fashion. Section 2.3 briefly sheds light on the Lyapunov-based approach
utilized by [12, 8] to design dense $\mathcal{H}_{\infty}$ controllers and
estimators for DAE-modeled power systems.

###  2.1 Power system NDAE state-space model

Here, we include the fourth-order dynamics of the modeled synchronous
generators. To this end, we borrow the notations from [3]. Let us consider a
power system comprised of $N_{b}$ buses, modeled by a graph
$(\mathcal{N},\mathcal{E})$ where $\mathcal{N}$ and $\mathcal{E}$ respectively
denote the set of nodes and edges. The nodes are of three main types: (i)
traditional synchronous generators, (ii) renewable energy resources (RERs),
and (iii) load buses. Note that
$\mathcal{N}=\mathcal{G}\cup\mathcal{R}\cup\mathcal{L}\cup\mathcal{U}$ holds
where $\mathcal{G}$, $\mathcal{R}$, $\mathcal{L}$, and $\mathcal{U}$
respectively collect $G$ generator buses, the buses containing $R$ renewables,
$L$ load buses, and $U$ non-unit buses.

A fourth-order dynamics of synchronous generators can be modeled as [3, 30,
31]

$\displaystyle\dot{\delta}_{i}$ $\displaystyle=\omega_{i}-\omega_{0},$ (1a)
$\displaystyle\begin{split}M_{i}\dot{\omega}_{i}&=T_{\mathrm{M}i}-P_{\mathrm{G}i}-D_{i}(\omega_{i}-\omega_{0}),\end{split}$
(1b) $\displaystyle T^{\prime}_{\mathrm{d0}i}\dot{E}^{\prime}_{i}$
$\displaystyle=-\tfrac{x_{\mathrm{d}i}}{x^{\prime}_{\mathrm{d}i}}E^{\prime}_{i}+\tfrac{x_{\mathrm{d}i}-x^{\prime}_{\mathrm{d}i}}{x^{\prime}_{\mathrm{d}i}}v_{i}\cos(\delta_{i}-\theta_{i})+E_{\mathrm{fd}i},$
(1c) $\displaystyle T_{\mathrm{CH}i}\dot{T}_{Mi}$
$\displaystyle=-T_{\mathrm{M}i}-\tfrac{1}{R_{\mathrm{D}i}}(\omega_{i}-\omega_{0})+T_{\mathrm{r}i},$
(1d)

where $(\delta_{i},\omega_{i},E^{\prime}_{i},T_{Mi})$ and $(E_{fdi},T_{ri})$
respectively represent the generator’s internal states and the generator’s
inputs. The generator’s internal states, namely,
$(\delta_{i},\omega_{i},E^{\prime}_{i},T_{Mi})$, its supplied power
$(P_{Gi},Q_{Gi})$, and terminal voltage $v_{i}$ are related to each other via
the following two algebraic equations [3, 30]:

$\displaystyle\begin{split}P_{\mathrm{G}i}&=\tfrac{\sin(\delta_{i}-\theta_{i})}{x^{\prime}_{\mathrm{d}i}}(E^{\prime}_{i}v_{i}-\tfrac{(x_{\mathrm{q}i}-x^{\prime}_{\mathrm{d}i})v_{i}^{2}\cos(\delta_{i}-\theta_{i})}{x_{\mathrm{q}i}}),\end{split}$
(2a)
$\displaystyle\begin{split}Q_{\mathrm{G}i}&=\tfrac{E^{\prime}_{i}v_{i}\cos(\delta_{i}-\theta_{i})-v_{i}^{2}\cos^{2}(\delta_{i}-\theta_{i})}{x^{\prime}_{\mathrm{d}i}}-\tfrac{v_{i}^{2}\sin^{2}(\delta_{i}-\theta_{i})}{x_{\mathrm{q}i}}.\end{split}$
(2b)

The power flow (PF)/balance equations governing the power transfer among
generators, RERs, and loads are as follows [3, 31]:

$\displaystyle\begin{split}P_{\mathrm{G}i}+P_{\mathrm{R}i}+P_{\mathrm{L}i}&=\sum_{j=1}^{N_{b}}v_{i}v_{j}\left(G_{ij}\cos\theta_{ij}+B_{ij}\sin\theta_{ij}\right),\end{split}$
(3a)
$\displaystyle\begin{split}Q_{\mathrm{G}i}+Q_{\mathrm{R}i}+Q_{\mathrm{L}i}&=\sum_{j=1}^{N_{b}}v_{i}v_{j}\left(G_{ij}\sin\theta_{ij}-B_{ij}\cos\theta_{ij}\right),\end{split}$
(3b)

where $\theta_{ij}:=\theta_{i}-\theta_{j}$. The pair $(G_{ij},B_{ij})$
respectively denotes the conductance and susceptance between buses $i$ and
$j$. Also, the pairs $(P_{Ri},Q_{Ri})$ and $(P_{Li},Q_{Li})$ respectively
represent the active and reactive powers generated by the renewables and the
active and reactive powers consumed by the loads.

Let us consider the following NDAE state-space model:

$\displaystyle E\dot{x}$ $\displaystyle=Ax+Bu+B_{w}w+h(x,u,w),$ (4a)
$\displaystyle u$ $\displaystyle=Fy,\leavevmode\nobreak\ y=C_{y}x+D_{y}w,$
(4b)

where $E$ denotes a singular matrix encoding the algebraic equations with all-
zeros rows,
$x=\begin{bmatrix}x_{d}^{T}&x_{a}^{T}\end{bmatrix}^{T}\in\mathbb{R}^{n_{x}}$,
$x_{d}\in\mathbb{R}^{n_{d}}$, $x_{a}\in\mathbb{R}^{n_{a}}$,
$u\in\mathbb{R}^{n_{u}}$, $w\in\mathbb{R}^{n_{w}}$, and
$y\in\mathbb{R}^{n_{y}}$ respectively represent the state, dynamic state,
algebraic state, control input, disturbance input, and measurement output
vectors, $F$ denotes the static output feedback (SOF) controller matrix, and
$A$, $B$, and $B_{w}$ are extracted via Jacobian-based linearization while
$h(x,u,w)$ encompasses the linearization error vector associated with the
nonlinearities. Here, we emphasize the point that we have not utilized a
purely linearized model by overlooking the linearization error vector
associated with the nonlinearities.

To construct the NDAE state-space representation of the power system for the
controller design, let us define

$\displaystyle x_{d}$
$\displaystyle=\begin{bmatrix}\delta^{T}&\omega^{T}&E^{\prime
T}&T_{M}^{T}\end{bmatrix}^{T},\leavevmode\nobreak\
x_{a}=\begin{bmatrix}a^{T}&\tilde{v}^{T}\end{bmatrix}^{T},\leavevmode\nobreak\
a:=\begin{bmatrix}P_{G}^{T}&Q_{G}^{T},\end{bmatrix}^{T},$
$\displaystyle\tilde{v}$
$\displaystyle:=\begin{bmatrix}v^{T}&\theta^{T}\end{bmatrix}^{T},\leavevmode\nobreak\
u=\begin{bmatrix}E_{fd}^{T}&T_{r}^{T}\end{bmatrix}^{T},\leavevmode\nobreak\
w=\begin{bmatrix}P_{R}^{T}&Q_{R}^{T}&P_{L}^{T}&Q_{L}^{T}\end{bmatrix}^{T},$

where

$\displaystyle\delta$
$\displaystyle:=\\{\delta_{i}\\}_{i\in\mathcal{G}},\leavevmode\nobreak\
\omega:=\\{\omega_{i}\\}_{i\in\mathcal{G}},\leavevmode\nobreak\
E^{\prime}:=\\{E^{\prime}_{i}\\}_{i\in\mathcal{G}},\leavevmode\nobreak\
T_{M}:=\\{T_{Mi}\\}_{i\in\mathcal{G}},$ $\displaystyle P_{G}$
$\displaystyle:=\\{P_{Gi}\\}_{i\in\mathcal{G}},\leavevmode\nobreak\
Q_{G}:=\\{Q_{Gi}\\}_{i\in\mathcal{G}},\leavevmode\nobreak\
v:=\\{v_{i}\\}_{i\in\mathcal{N}},\leavevmode\nobreak\
\theta:=\\{\theta_{i}\\}_{i\in\mathcal{N}},$ $\displaystyle E_{fd}$
$\displaystyle:=\\{E_{fdi}\\}_{i\in\mathcal{G}},\leavevmode\nobreak\
T_{r}:=\\{T_{ri}\\}_{i\in\mathcal{G}},\leavevmode\nobreak\
P_{R}:=\\{P_{Ri}\\}_{i\in\mathcal{R}},\leavevmode\nobreak\
Q_{R}:=\\{Q_{Ri}\\}_{i\in\mathcal{R}},$ $\displaystyle P_{L}$
$\displaystyle:=\\{P_{Li}\\}_{i\in\mathcal{L}},\leavevmode\nobreak\
Q_{L}:=\\{Q_{Li}\\}_{i\in\mathcal{L}}.$

Based on the above vector representation and (1), (2), and (3), the NDAE
state-space model (4a) for the power system can be constructed.

###  2.2 Problem statement

Assuming the $\mathcal{L}_{2}$-norm boundedness of $h(x,u,w)$ and considering
it as $B_{h}w_{h}$, the NDAE state-space model (4) reduces to the following
linear DAE (LDAE) form:

$\displaystyle E\dot{x}$
$\displaystyle=(A+BFC_{y})x+\hat{B}_{w}\tilde{w},\leavevmode\nobreak\
\tilde{w}:=\begin{bmatrix}w^{T}&w_{h}^{T}\end{bmatrix}^{T},$ (5a)
$\displaystyle\hat{B}_{w}$
$\displaystyle:=\begin{bmatrix}B_{w}+BFD_{y}&B_{h}\end{bmatrix}.$ (5b)

With that in mind, to design an SOF $\mathcal{H}_{\infty}$ controller, we
define the following performance output vector:

$\displaystyle z$
$\displaystyle=Cx+Du+\hat{D}_{w}\tilde{w},\leavevmode\nobreak\
\hat{D}_{w}:=\begin{bmatrix}D_{w}&0\end{bmatrix}.$ (6)

The main idea is to design an SOF $\mathcal{H}_{\infty}$ controller $F$ in
(4b) such that it satisfies $\mathcal{H}_{\infty}$ performance criterion,
namely, $\left\lVert
z\right\rVert^{2}_{\mathcal{L}_{2}}<\mu^{2}\left\lVert{\tilde{w}}\right\rVert^{2}_{\mathcal{L}_{2}}$
[32] ($\mu$: the $\mathcal{H}_{\infty}$ value). In the sequel, we elaborate on
structurally-constrained $\mathcal{H}_{\infty}$ SOF controller design which
frequently appears in more realistic scenarios. Further details about the
design of such an SOF controller $F$ are later given in Section 3.

The sparsity structure of the SOF controller $F$ is of high significance in
decreasing the communication burden among generators. To name a few common
sparsity structures, we list the following two crucial sparsity structures:
(i) decentralized, and (ii) distributed. In the context of power systems, each
generator has access only to its own state information for a decentralized
sparsity structure while in the case of a distributed sparsity structure, each
generator has access to the state information of a few numbers of the other
generators. To encode the sparsity structures, we simply use the matrix
element-wise Hadamard product as follows:

$\displaystyle F$ $\displaystyle=F\odot S,$ (7)

where $S\in\\{0,1\\}^{n_{u}\times n_{y}}$ represents an imposed binary
sparsity structure. To that end, to design the structured SOF (SSOF)
$\mathcal{H}_{\infty}$ controller for the NDAE state-space model (4), we
define the following problem statement:

###### Problem 1

Given the NDAE state-space model (4) and its extracted LDAE form (5) along
with the performance output vector (6) and an imposed controller sparsity
structure (7), design an SSOF $\mathcal{H}_{\infty}$ controller.

###  2.3 Lyapunov-based approach

This section briefly details the Lyapunov-based approach utilized by [12, 8]
to design dense $\mathcal{H}_{\infty}$ controllers and estimators for DAE-
modeled power systems.

####  2.3.1 Dense H-infinity controller design

Considering (6) and LDAE state-space model (5) with $C_{y}=I_{n_{x}},D_{y}=0$,
dense static state feedback (DSSF) $\mathcal{H}_{\infty}$ controller is
proposed in [12] for which a convex SDP is cast and solved (as shown in A). To
derive the LMIs of the cast convex SDP therein, the $\mathcal{H}_{\infty}$
inequality (8) (based on KYP Lemma) is utilized as a cornerstone

$\displaystyle\dot{V}(x)+z^{T}z-\lambda\tilde{w}^{T}\tilde{w}<0,$ (8)

where $V(x)=x^{T}E^{T}Px$ denotes a quadratic Lyapunov candidate function for
which $E^{T}P=P^{T}E\succeq 0$ holds. Note that to implement the LMIs with
strict definiteness, one has to utilize $\epsilon I$-shifting and
appropriately choose/tune $\epsilon$ to get a well-performing controller
design which is a time-consuming process.

####  2.3.2 Dense H-infinity estimator design

Considering the following NDAE state-space error dynamics model [8]:

$\displaystyle E\dot{e}$ $\displaystyle=(A+LC_{y})e+B_{w}w+\Delta f,$ (9)

where $e:=\Delta x=x-\hat{x}$ represents the error between the estimated
$\hat{x}$ and the true $x$ values of state variables, $L$ denotes the
estimator, and $\|\Delta f\|$ is assumed to be Lipschitz bounded by
$\alpha\|\Delta x\|$ ($\alpha$: the Lipschitz constant), a dense
$\mathcal{H}_{\infty}$ estimator is proposed in [8] for which a convex SDP is
cast and solved (as shown in B). To derive the LMIs of the cast convex SDP
therein, the $\mathcal{H}_{\infty}$ inequality (10) (similarly, based on KYP
Lemma) and the $\mathcal{S}$-procedure lemma [33] are utilized as cornerstones

$\displaystyle\dot{V}(e)+e^{T}e-\lambda w^{T}w<0,$ (10)

where $V(e)=e^{T}E^{T}Pe$ denotes a quadratic Lyapunov candidate function for
which $E^{T}P=P^{T}E\succeq 0$ holds. Likewise, the fine-tuning of $\epsilon
I$-shifting in LMIs is inevitable.

In the next section, we present the proposed methodology to design dense and
structured $\mathcal{H}_{\infty}$ controllers and estimators that ensure the
power system is stable in the sense of $\mathcal{H}_{\infty}$.

##  3 Structured H-infinity Designs for DAE-Modeled Power Systems Via NLA

In Section 3.1, we limit our attention to the derivation of preliminary
casting for structured $\mathcal{H}_{\infty}$ controllers for DAE-modeled
power systems. Similar derivation for the structured $\mathcal{H}_{\infty}$
estimator design can be obtained as presented in Section 3.3. Section 3.2 is
allocated to elaborate on the utilized non-Lyapunov approach (NLA) and the
thorough details of its implementation for the controller design. In Section
3.3, NLA is similarly implemented to obtain the estimator design. Procedures
Procedure 1 and Procedure 2 in Sections 3.2 and 3.3, contain the systematic
NLA of designing structured $\mathcal{H}_{\infty}$ controllers and estimators
specialized for the DAE-modeled power systems, respectively. Note that as a
notation rule, the subscripts $d$ and $a$ refer to the corresponding dynamic
and algebraic components, respectively. Also, for the sake of simplicity in
derivations, without loss of generality, we assume that
$E=\begin{bmatrix}I_{n_{d}}&0\\\ 0&0\end{bmatrix}$ holds.

###  3.1 Preliminary casting for the controller design

Let us consider

$\displaystyle A=\begin{bmatrix}A_{dd}&A_{da}\\\
A_{ad}&A_{aa}\end{bmatrix},\leavevmode\nobreak\ B=\begin{bmatrix}B_{d}\\\
B_{a}\end{bmatrix},\leavevmode\nobreak\
C_{y}=\begin{bmatrix}C_{yd}&C_{ya}\end{bmatrix},$
$\displaystyle\mathcal{A}_{ij}(F):=A_{ij}+B_{i}FC_{yj},\leavevmode\nobreak\
i,j\in\\{d,a\\},$ $\displaystyle
B_{w}=\begin{bmatrix}B_{wd}^{T}&B_{wa}^{T}\end{bmatrix}^{T},\leavevmode\nobreak\
h=\begin{bmatrix}h_{d}^{T}&h_{a}^{T}\end{bmatrix}^{T},$ $\displaystyle
B_{w}^{i}(F):=B_{wi}+B_{i}FD_{y},\leavevmode\nobreak\ i\in\\{d,a\\}.$

According to the algebraic equations of the NDAE state-space model (4), we
have

$\displaystyle\mathcal{A}_{aa}(F)x_{a}=$
$\displaystyle-(\mathcal{A}_{ad}(F)x_{d}+B_{w}^{a}(F)w+h_{a}).$ (11)

Assuming the invertibility of $\mathcal{A}_{aa}(F)$ (less restrictive than the
invertibility of $A_{aa}$ assumed in [12]) and utilizing (11), we can
eliminate $x_{a}$ from the NDAE state-space model (4) and get the following
ENODE state-space model of the NDAE state-space model (4):

$\displaystyle\dot{x}_{d}=\tilde{A}(F)x_{d}+\bar{B}_{w}(F)w+\tilde{h}(F),$
(12a)
$\displaystyle\tilde{A}(F):=\mathcal{A}_{dd}(F)-\mathcal{A}_{da}(F)\mathcal{A}_{aa}(F)^{-1}\mathcal{A}_{ad}(F),$
(12b)
$\displaystyle\bar{B}_{w}(F):=B_{w}^{d}(F)-\mathcal{A}_{da}(F)\mathcal{A}_{aa}(F)^{-1}B_{w}^{a}(F),$
(12c)
$\displaystyle\tilde{h}(F):=h_{d}-\mathcal{A}_{da}(F)\mathcal{A}_{aa}(F)^{-1}h_{a}.$
(12d)

Assuming the $\mathcal{L}_{2}$-norm boundedness of $h(x,u,w)$ and considering
$h(x,u,w)$ as $B_{h}w_{h}$, the ENODE state-space model (12) boils down to the
following linear ODE (LODE) form:

$\displaystyle\boxed{\dot{x}_{d}=\tilde{A}(F)x_{d}+\tilde{B}_{w}(F)\tilde{w},}$
(13a)
$\displaystyle\tilde{B}_{w}(F):=\begin{bmatrix}\bar{B}_{w}(F)&\bar{B}_{h}(F)\end{bmatrix},\leavevmode\nobreak\
\bar{B}_{h}(F):=B_{hd}-\mathcal{A}_{da}(F)\mathcal{A}_{aa}(F)^{-1}B_{ha},$
(13b) $\displaystyle
B_{h}=\begin{bmatrix}B_{hd}^{T}&B_{ha}^{T}\end{bmatrix}^{T},$ (13c)

and defining

$\displaystyle\check{C}(F):=C+DFC_{y},\leavevmode\nobreak\
\check{D}_{w}(F):=D_{w}+DFD_{y},$

the performance output vector (6) boils down to

$\displaystyle\boxed{z=\tilde{C}(F)x_{d}+\tilde{D}_{w}(F)\tilde{w},}$ (14a)
$\displaystyle\tilde{C}(F):=\check{C}(F)\begin{bmatrix}I_{n_{d}}\\\
-\mathcal{A}_{aa}(F)^{-1}\mathcal{A}_{ad}(F)\end{bmatrix},\leavevmode\nobreak\
\tilde{D}_{w}(F):=\begin{bmatrix}\bar{D}_{w}(F)&\bar{D}_{h}(F)\end{bmatrix},$
(14b)
$\displaystyle\bar{D}_{w}(F):=\check{D}_{w}(F)+\check{C}(F)\begin{bmatrix}0\\\
-\mathcal{A}_{aa}(F)^{-1}B_{w}^{a}(F)\end{bmatrix},$ (14c)
$\displaystyle\bar{D}_{h}(F):=\check{C}(F)\begin{bmatrix}0\\\
-\mathcal{A}_{aa}(F)^{-1}B_{ha}\end{bmatrix}.$ (14d)

For LODE (13) and the reduced performance output vector (14), we define the
transfer function $T_{z\tilde{w}}(s)$ from the disturbance input vector
$\tilde{w}$ to the performance output vector $z$ in (14) as follows:

$\displaystyle T_{z\tilde{w}}(s)$
$\displaystyle:=\tilde{C}(F)(sI_{n_{d}}-\tilde{A}(F))^{-1}\tilde{B}_{w}(F)+\tilde{D}_{w}(F),$

and denote its $\mathcal{H}_{\infty}$ norm by
$\|T_{z\tilde{w}}(s)\|_{\mathcal{H}_{\infty}}$ which is defined as follows
[11]:

$\displaystyle\|T_{z\tilde{w}}(s)\|_{\mathcal{H}_{\infty}}:=\underset{\Re(s)>0}{\sup}\leavevmode\nobreak\
\bar{\sigma}[T_{z\tilde{w}}(s)]\overset{\mathrm{for\leavevmode\nobreak\
a\leavevmode\nobreak\ stable}\leavevmode\nobreak\
T_{z\tilde{w}}(s)}{=}\|T_{z\tilde{w}}(s)\|_{\mathcal{L}_{\infty}}=\underset{\omega\in\mathbb{R}}{\sup}\leavevmode\nobreak\
\bar{\sigma}[T_{z\tilde{w}}(j\omega)].$

Note that for an unstable $T_{z\tilde{w}}(s)$ (the closed-loop system with a
destabilizing controller $F$, i.e., $\mathbf{sa}(\tilde{A}(F))\nless 0$)
$\|T_{z\tilde{w}}(s)\|_{\mathcal{H}_{\infty}}=\infty$ holds. To compute
$\|T_{z\tilde{w}}(s)\|_{\mathcal{H}_{\infty}}$, we utilize the following
MATLAB built-in functions:

$\displaystyle\boxed{\texttt{sys}=\texttt{ss}(\tilde{A}(F),\tilde{B}_{w}(F),\tilde{C}(F),\tilde{D}_{w}(F)),}$
(15a)
$\displaystyle\|T_{z\tilde{w}}(s)\|_{\mathcal{H}_{\infty}}=\texttt{hinfnorm}(\texttt{sys}).$
(15b)

###  3.2 NLA and its implementation for the controller design

In this section, we propose a solution to Problem 1 on the basis of NLA. Given
the NDAE state-space model (4), assuming the $\mathcal{L}_{2}$-norm
boundedness of the linearization error vector associated with the
nonlinearities, we extract the LODE state-space model (13) along with the
reduced performance output vector (14) and utilize NLA to design an SSOF
$\mathcal{H}_{\infty}$ controller subject to an imposed sparsity structure
encoded by (7).

As NLA, we utilize the MATLAB built-in function hinfstruct which has been
developed based on a structured $\mathcal{H}_{\infty}$ synthesis originally
proposed by [23]. The proposed synthesis in [23] is mainly built upon non-
convex non-smooth optimization techniques taking advantage of the Clarke sub-
differential [24] computation of the $\mathcal{H}_{\infty}$ norm. Such a
structured $\mathcal{H}_{\infty}$ synthesis has found extensive control
applications including but not limited to proportional-integral-derivative
(PID) control and combination of proportional-integral (PI) and SOF control of
distillation column, aircraft autopilot, and helicopter as highlighted by [25,
26, 27]. Note that hinfstruct computes the $\mathcal{H}_{\infty}$ norm based
on efficient algorithms developed by [34, 35, 36].

Procedure Procedure 1 summarizes NLA to design an SSOF $\mathcal{H}_{\infty}$
controller. In Procedure Procedure 1, $F_{0}$ denotes the initialization for
$F$ and the MATLAB built-in function realp$(F,\textrm{`}F_{0}\textrm{'})$
creates a real-valued free parameter initialized by $F_{0}$. To encode the
imposed sparsity structure $S$ in (7), we utilize F.Free(i,j) = false for all
zero elements of $S$. It is remarkable that extra design specifications can be
incorporated via options embedded in hinfstruct. For example, the maximum
closed-loop natural frequency and the minimum decay rate for closed-loop poles
can be incorporated via options.MaxFrequency and options.MinDecay in
hinfstruct, respectively. Also, to avoid high-gain controller designs or
unwanted fast dynamics, one may set options.MaxFrequency to a finite value.

2 input: $A$, $B$, $B_{w}$, $B_{h}$, $C_{y}$, $D_{y}$, $C$, $D$, $D_{w}$, $S$.

3set $F_{0}=0$.

4set $F$ = realp($F$,‘$F_{0}$’).

5for  _$i=1:n_{u}$_ do

6 for  _$j=1:n_{y}$_ do

7 if  _$S(i,j)=0$_ then

8 $\texttt{F.Free}(i,j)=\texttt{false}$.

9

10

11

12construct $\tilde{A}(F)$, $\tilde{B}_{w}(F)$, $\tilde{C}(F)$,
$\tilde{D}_{w}(F)$ in (13) and (14).

13construct sys via (15).

14utilize hinfstruct(sys) to obtain $F^{\ast}$.

15output: $F^{\ast}$.

Procedure 1 SSOF $\mathcal{H}_{\infty}$ Controller Design

######  Remark 1 (On hinfstruct)

Borrowing from [25], we present a high-level description of the structured
$\mathcal{H}_{\infty}$ solver implemented as a building block of hinfstruct.
The detailed explanation of the theoretical aspects of hinfstruct can be found
in [25, 23].
Structured $\mathcal{H}_{\infty}$ synthesis deals with the following semi-
infinite, non-convex, and non-smooth optimization problem:

$\displaystyle\underset{K}{\min}\leavevmode\nobreak\
\|T_{z\tilde{w}}(s)\|_{\mathcal{H}_{\infty}}\iff\underset{\kappa}{\min}\leavevmode\nobreak\
\underset{\omega\in[0,\infty]}{\max}\leavevmode\nobreak\
\bar{\sigma}[T_{z\tilde{w}}(j\omega)],$ (16)

where $K$ denotes the structured controller/estimator and $\kappa$ represents
all the low-level tunable free parameters. It is noteworthy that the function
in the right-hand-side of (16) is the composition of the convex non-smooth
function $\max_{\omega}\circ\leavevmode\nobreak\ \bar{\sigma}(.)$ with the
non-convex differentiable mapping $\kappa\longrightarrow
T_{z\tilde{w}}(j\omega)$. Fortunately, such composite functions are Clarke
regular [24] meaning that a complete description of the Clarke subdifferential
is available. For the sake of simplicity in representation, (16) can be
rewritten as

$\displaystyle\underset{\kappa}{\min}\leavevmode\nobreak\
[f_{\infty}(\kappa):=\underset{\omega\in[0,\infty]}{\max}\leavevmode\nobreak\
f(\omega,\kappa)].$ (17)

Clarke regularity ensures that critical point $\kappa^{\ast}$ (usually local
minima) are characterized via $0\in\partial f(\kappa)$.
To solve (17), the developers of hinfstruct [25, 23] construct a tangent model
around the current iterate $\kappa$ that constitutes quadratic first-order
local approximation of the original problem. An adequate descent direction
$\eta$ is then computed by solving the following convex quadratic optimization
problem:

$\displaystyle\underset{\eta}{\min}\leavevmode\nobreak\ $
$\displaystyle\hat{f}_{\infty}(\kappa+\eta),$ (18a)
$\displaystyle\hat{f}_{\infty}(\kappa+\eta):=[\underset{\omega\in\Omega_{f}}{\max}\leavevmode\nobreak\
$ $\displaystyle
f(\omega,\kappa)]-f_{\infty}(\kappa)+\phi_{\omega}^{T}\eta+\frac{1}{2}\eta^{T}Q\eta,$
(18b)

where $\Omega_{f}$ denotes some finite set of frequencies, and
$\phi_{\omega}\in\partial f(\omega,\kappa)$ represents a subgradient of
$f(\omega,\kappa)$. A minimal requirement to implement this scheme is that
$\Omega_{f}$ contains the active frequencies $\omega_{\mathrm{active}}$
achieving the peak value in (17), namely,
$f_{\infty}(\kappa)=f(\omega_{\mathrm{active}},\kappa)$. Such a simple
requirement is sufficient for the convergence of the algorithm. Nevertheless,
by adding a few extra well-chosen frequencies, one can often improve the
quality of the tangent model (18) and take longer steps at each iteration.
It is noteworthy that the algorithm implemented in hinfstruct [23], starting
from an initial guess $K_{\mathrm{ig}}$ searches for a stabilizing starting
point $K_{0}$ via minimizing the $a$-shifted $\mathcal{H}_{\infty}$ norm [37]
of the closed-loop control system where the shift $a>0$ is usually kept fixed
at the initial $a_{0}>0$ for which the $a$-shifted $\mathcal{H}_{\infty}$ norm
of the closed-loop control system with $K_{\mathrm{ig}}$ is finite. The
authors in [23] highlight that the reason for not utilizing the spectral
abscissa $\mathbf{sa}$ minimization (to obtain $K_{0}$) is that minimizing the
$a$-shifted $\mathcal{H}_{\infty}$ norm is relatively more compatible with
their proposed structured $\mathcal{H}_{\infty}$ synthesis.

###  3.3 Casting and NLA implementation for the estimator design

To design structured $\mathcal{H}_{\infty}$ estimators for DAE-modeled power
systems, a problem similar to as stated in Problem 1 can accordingly be
stated. NLA to propose a structured $\mathcal{H}_{\infty}$ estimator is quite
the same as the one proposed for an SSOF $\mathcal{H}_{\infty}$ controller. In
the sequel, we present a brief set of similar key formulas. To that end, let
us consider the following Luenberger-type estimator dynamics:

$\displaystyle\begin{split}E\dot{\hat{x}}&={A}{\hat{x}}+L\left(y-\hat{y}\right)+{{B}}{{{u}}}+B_{w}w_{0}+h(\hat{x},u,w),\end{split}$
(19a) $\displaystyle\begin{split}\hat{y}&=C\hat{x},\end{split}$ (19b)

where $\hat{x}$ and $\hat{y}$ denote the estimated states and outputs,
respectively and $w_{0}$ represents the steady-state value of $w$ which
contains load demands and renewable generations. Now, defining the error as
$e:=x-\hat{x}$, then, the model of error dynamics can be computed as given in
(9). The main idea is designing $L$ such that the error model (9)
asymptotically converges to zero. To that end, considering

$\displaystyle
e=\begin{bmatrix}e_{d}^{T}&e_{a}^{T}\end{bmatrix}^{T},\leavevmode\nobreak\
\mathcal{A}_{ij}(L):=A_{ij}+L_{i}C_{yj},\leavevmode\nobreak\ i,j\in\\{d,a\\},$
$\displaystyle\Delta f=\begin{bmatrix}{\Delta f}_{d}^{T}&{\Delta
f}_{a}^{T}\end{bmatrix}^{T},\leavevmode\nobreak\
L=\begin{bmatrix}L_{d}^{T}&L_{a}^{T}\end{bmatrix}^{T},\leavevmode\nobreak\
C_{y}=\begin{bmatrix}C_{yd}&C_{ya}\end{bmatrix},$

and assuming the invertibility of $\mathcal{A}_{aa}(L)$ to eliminate $e_{a}$
from the NDAE state-space model (9) via

$\displaystyle\mathcal{A}_{aa}(L)e_{a}$
$\displaystyle=-(\mathcal{A}_{ad}(L)e_{d}+B_{wa}w+{\Delta f}_{a}),$

we get the following ENODE state-space model of the NDAE state-space model
(9):

$\displaystyle\dot{e}_{d}=\tilde{A}(L)e_{d}+\bar{B}_{w}(L)w+\Delta\tilde{f}(L),$
(20a)
$\displaystyle\tilde{A}(L):=\mathcal{A}_{dd}(L)-\mathcal{A}_{da}(L)\mathcal{A}_{aa}(L)^{-1}\mathcal{A}_{ad}(L),$
(20b)
$\displaystyle\bar{B}_{w}(L):=B_{wd}-\mathcal{A}_{da}(L)\mathcal{A}_{aa}(L)^{-1}B_{wa},$
(20c) $\displaystyle\Delta\tilde{f}(L):={\Delta
f}_{d}-\mathcal{A}_{da}(L)\mathcal{A}_{aa}(L)^{-1}{\Delta f}_{a}.$ (20d)

Assuming the $\mathcal{L}_{2}$-norm boundedness of $\Delta f$ and considering
it as $B_{\Delta f}w_{\Delta f}$, the ENODE state-space model (20) reduces to
the following LODE form:

$\displaystyle\boxed{\dot{e}_{d}=\tilde{A}(L)e_{d}+\tilde{B}_{w}(L)\tilde{w},}$
(21a)
$\displaystyle\tilde{B}_{w}(L):=\begin{bmatrix}\bar{B}_{w}(L)&\bar{B}_{\Delta
f}(L)\end{bmatrix},$ (21b) $\displaystyle\bar{B}_{\Delta f}(L):=B_{\Delta
fd}-\mathcal{A}_{da}(L)\mathcal{A}_{aa}(L)^{-1}B_{\Delta
fa},\leavevmode\nobreak\ B_{\Delta f}=\begin{bmatrix}B_{\Delta
fd}^{T}&B_{\Delta fa}^{T}\end{bmatrix}^{T},$ (21c)

and the performance output vector $z=e$ reduces to

$\displaystyle\boxed{z=\tilde{C}(L)e_{d}+\tilde{D}_{w}(L)\tilde{w},}$ (22a)
$\displaystyle\tilde{C}(L):=\begin{bmatrix}I_{n_{d}}\\\
-\mathcal{A}_{aa}(L)^{-1}\mathcal{A}_{ad}(L)\end{bmatrix},\leavevmode\nobreak\
\tilde{D}_{w}(L):=\begin{bmatrix}\bar{D}_{w}(L)&\bar{D}_{\Delta
f}(L)\end{bmatrix},$ (22b) $\displaystyle\bar{D}_{w}(L):=\begin{bmatrix}0\\\
-\mathcal{A}_{aa}(L)^{-1}B_{wa}\end{bmatrix},\leavevmode\nobreak\
\bar{D}_{\Delta f}(L):=\begin{bmatrix}0\\\ -\mathcal{A}_{aa}(L)^{-1}B_{\Delta
fa}\end{bmatrix}.$ (22c)

For LODE (21) and the reduced performance output vector (22), we have

$\displaystyle\boxed{\texttt{sys}=\texttt{ss}(\tilde{A}(L),\tilde{B}_{w}(L),\tilde{C}(L),\tilde{D}_{w}(L)),}$
(23a)
$\displaystyle\|T_{z\tilde{w}}(s)\|_{\mathcal{H}_{\infty}}=\texttt{hinfnorm}(\texttt{sys}).$

Note that we can impose a sparsity structure to the estimator $L$ via
$L=L\odot S$ similar to as in (7). With a slight modification to Procedure
Procedure 1, we obtain a similar Procedure to design a structured
$\mathcal{H}_{\infty}$ estimator, namely, Procedure Procedure 2, where $L_{0}$
denotes the initialization for $L$.

2

3input: $A$, $C_{y}$, $B_{w}$, $B_{\Delta f}$, $S$.

4set $L_{0}=0$.

5set $L$ = realp($L$,‘$L_{0}$’).

6for  _$i=1:n_{x}$_ do

7 for  _$j=1:n_{y}$_ do

8 if  _$S(i,j)=0$_ then

9 $\texttt{L.Free}(i,j)=\texttt{false}$.

10

11

12

13construct $\tilde{A}(L)$, $\tilde{B}_{w}(L)$, $\tilde{C}(L)$,
$\tilde{D}_{w}(L)$ in (13) and (14).

14construct sys via (23).

15utilize hinfstruct(sys) to obtain $L^{\ast}$.

16output: $L^{\ast}$.

Procedure 2 Structured $\mathcal{H}_{\infty}$ Estimator Design

##  4 Numerical Simulations

In this section, we assess the effectiveness of NLA in designing structured
controllers and estimators for the DAE-modeled power systems. We design
various types of controllers and estimators with different sparsity structures
such as centralized (dense), decentralized, and distributed architectures.
Notice that in the proposed methodology, this can easily be achieved by
appropriately selecting the sparsity structure $S$. With that in mind, we
assess the performance of the proposed methodology on four IEEE test systems
with the following specifications: (i) IEEE $9$-bus with
$(n_{x},n_{d},n_{a},n_{u},n_{w})=(36,12,24,6,18)$, (ii) IEEE $14$-bus with
$(n_{x},n_{d},n_{a},n_{u},n_{w})=(58,20,38,10,28)$, (iii) IEEE $39$-bus with
$(n_{x},n_{d},n_{a},n_{u},n_{w})=(138,40,98,20,78)$, and (iv) IEEE $57$-bus
with $(n_{x},n_{d},n_{a},n_{u},n_{w})=(156,28,128,14,114)$. To compute the
proposed structured $\mathcal{H}_{\infty}$ controllers and estimators, we run
Procedure Procedure 1 and Procedure Procedure 2, respectively. All the
numerical experiments have been run in MATLAB R$2022$b on a MacBook Pro with a
$3.1$ GHz Intel Core i$5$ and memory $8$ GB $2133$ MHz. For solving the convex
SDPs, we have utilized CVX [38] to implement the LMIs in MOSEK [14] SDP
solver. We simply set $D_{w}=0$ for all the numerical experiments.

Through the comprehensive case studies in this section, we seek answers to the
following questions:

  * 1.

Q1: How efficient is NLA in comparison with the Lyapunov-based approach in
terms of computation scalability and computation sensitivity?

  * 2.

Q2: Under what conditions NLA can potentially face a scalability issue leading
to even underperformance compared to the Lyapunov-based approach?

  * 3.

Q3: Can NLA be utilized to propose $\mathcal{H}_{\infty}$ controller designs
subject to the imposed structural/sparsity constraints with partially-accessed
states in the presence of noisy measurement outputs?

  * 4.

Q4: How can output matrix selection affect the $\mathcal{H}_{\infty}$
performance and the computational time corresponding to the SSOF controller
design built upon NLA?

  * 5.

Q5: What advantages can be achieved via the proposed structured
$\mathcal{H}_{\infty}$ estimator design built upon NLA?

(a)  (b)

(c)  (d)

Figure 1: (Section 4.1) The generator frequencies for $9$-bus (top-left),
$14$-bus (top-right), $39$-bus (bottom-left), and $57$-bus (bottom-right) test
systems, for moderate disturbance in load demand and renewable power
generation.

(a)  (b)

(c)  (d)

Figure 2: (Section 4.1) Generator $3$ frequency and power for $9$-bus (top)
and $57$-bus (bottom) test systems, under large disturbance in load demand and
renewable power generation.

###  4.1 Dense controller design with fully-accessed states and noiseless
measurement states

In this section, we design a full-state DSSF $\mathcal{H}_{\infty}$ controller
with no noise in the measurement states. This can be achieved by setting
$S=\mathbf{1}_{n_{u}\times n_{y}}$, $C_{y}=I_{n_{x}}$, and $D_{y}=0$.
Furthermore, we consider $B_{h}=0.1I_{n_{x}}$. To showcase the superiority of
the proposed methodology, we compare the performance of the proposed
methodology with the conventional LQR-type controller presented in [39] and
the Lyapunov-based controllers given in [12, 3] for the various IEEE test
systems.

We observe from Tab. 1 that as the size of the system increases, the
computational time for calculating the controller matrix $F$ using the
Lyapunov-based approach [12] increases significantly while for NLA, the
computational time is much less and is more than $47$
($\frac{5361.67-111.31}{111.31}\approx 47.17$) times faster on the largest
IEEE test system, i.e., IEEE $57$-bus. We can also see that the proposed
methodology yields exactly the same $\mathcal{H}_{\infty}$ performance value
as compared to [12]. This means that the performance of the proposed
controller is not degraded and hence provides the same stability improvement
on power systems after disturbances.

Also, we observe that unlike the Lyapunov-based approach [12], NLA controller
is not a high-gain controller which makes it more suitable for practical
purposes. It is noteworthy that the corresponding values of the Lyapunov-based
approach [12] reported in Tab. 1 are obtained after fine-tuning of $\epsilon
I$-shifting of the corresponding LMIs (setting $\epsilon=10^{-3.3}$) which is
a time-consuming process by itself. Nevertheless, we have not taken into
account the time of fine-tuning in computational times reflected in Tab. 1.
Moreover, through the numerical experiments, we observe that although the
approach proposed in [12] is capable of designing the dense
$\mathcal{H}_{\infty}$ controllers with the same $\mathcal{H}_{\infty}$
performance level as the presented NLA, there is a noticeable gap between the
value of $\lambda:=\mu^{2}$ in the corresponding convex SDP and the true
$\mathcal{H}_{\infty}$-squared value i.e.,
$\|T_{z\tilde{w}}(s)\|_{\mathcal{H}_{\infty}}^{2}$ evaluated by hinfnorm.
Therefore, we have reported the corresponding true $\mathcal{H}_{\infty}$
values in Tab. 1.

That being said, we now assess the performance of NLA controller in
stabilizing the power system under transient conditions. To that end, large
disturbances in load and renewable (which are modeled as negative loads) have
been added and the simulations are performed as follows: Initially, the system
load demand and power generation are exactly equal, and thus the system rests
in equilibrium conditions. Then, immediately after $t>0$, there is an abrupt
change in load demand and renewable power generation, and their new values are
given as follows:

$\displaystyle
P^{\prime}_{L}+Q^{\prime}_{L}=(1+\Delta_{L})(P^{0}_{L}+Q^{0}_{L}),\leavevmode\nobreak\
P^{\prime}_{R}+Q^{\prime}_{R}=(1+\Delta_{R})(P^{0}_{L}+Q^{0}_{L}),$

where $P^{0}_{L}$, $Q^{0}_{L}$, $P^{0}_{R}$, and $Q^{0}_{R}$ denote the
initial values of active/reactive load demand and renewable generation, while
$P^{\prime}_{L}$, $Q^{\prime}_{L}$, $P^{\prime}_{R}$, and $Q^{\prime}_{R}$
denote their values after the disturbances, and $\Delta_{L}$ and $\Delta_{R}$
denote the severity of the disturbances in load and renewable power
generation, respectively. We apply different severity of disturbances by
choosing different values of $\Delta_{L}$, and $\Delta_{R}$. To that end, we
run two simulation studies for all the IEEE test systems. For $9$-bus test
system, we choose $(\Delta_{L},\Delta_{R})=(0.5,-0.3)$ for the first case
study and $(\Delta_{L},\Delta_{R})=(0.9,-0.5)$ for the second case study.
Similarly, for $14$-bus and $39$-bus test systems, we respectively select
$(\Delta_{L},\Delta_{R})=(0.05,-0.03)$ and
$(\Delta_{L},\Delta_{R})=(0.15,-0.15)$, while for $57$-bus system, we use
$(\Delta_{L},\Delta_{R})=(0.01,-0.01)$ for the first case study and
$(\Delta_{L},\Delta_{R})=(0.05,-0.05)$ for the second case study.

The results are presented in Figs. 1 and 2. Notice that in all the
visualizations, we refer to the controller proposed in [3] as LRFC (load and
renewable following control) and to the one proposed in [12] as Lyapunov-based
(abbreviated as Lyap) while to the controller proposed in the current paper as
NLA (non-Lyap). From Fig. 1, we observe that for a moderate disturbance in
load and renewable, the LQR-type controller proposed in [39] is unable to
stabilize the system while the Lyapunov-based controllers [12, 3] (note that
LRFC is also classified as a Lyapunov-based controller) and NLA controller can
successfully keep the system stable and synchronized.

Furthermore, from Fig. 2, we observe that as the severity of the disturbances
in load and renewable increases, only the controller proposed in [12] and the
presented controller in this paper are able to keep the system synchronized.
This corroborates the result presented in Tab. 1 from which, we see that the
Lyapunov-based approach and presented NLA $\mathcal{H}_{\infty}$ attain almost
the same $\mathcal{H}_{\infty}$ value, thus, providing almost the same
performance. However, the main benefit of the presented approach over [12] is
that it is far less computationally cumbersome. This holds because unlike [12]
no convex SDP needs to be cast and solved to compute the controller matrix
$F$. Moreover, the reason behind the fact that controllers in [3, 39] are
unable to keep the power system synchronized, is that the disturbance in load
and renewable are not modeled in the controller architecture. In both of these
studies, the disturbances are assumed to be zero in the controller design thus
making them less robust.

Table 1: (Section 4.1) The $\mathcal{H}_{\infty}$ values and computational times for the Lyapunov-based approach [12] and NLA corresponding to DSSF $\mathcal{H}_{\infty}$ controller with $S=\mathbf{1}_{n_{u}\times n_{y}}$, $C_{y}=I_{n_{x}}$, $D_{y}=0$, and $B_{h}=0.1I_{n_{x}}$ for the IEEE test systems. Approach | $\|T_{z\tilde{w}}(s)\|_{\mathcal{H}_{\infty}}$ | Computational Time | $(n_{u},n_{y})$
---|---|---|---
non-Lyap $9$-bus | $3.1856$ | $2.09$ s | $(6,36)$
Lyap $9$-bus | $3.1857$ |  $2.72$ s | $(6,36)$
non-Lyap $14$-bus | $7.4260$ | $6.05$ s | $(10,58)$
Lyap $14$-bus | $7.4274$ |  $19.41$ s | $(10,58)$
non-Lyap $39$-bus | $9.2018$ | $259.56$ s | $(20,138)$
Lyap $39$-bus | $9.1701$ |  $1611.45$ s | $(20,138)$
non-Lyap $57$-bus | $23.5716$ | $111.31$ s | $(14,156)$
Lyap $57$-bus | $23.5731$ |  $5361.67$ s | $(14,156)$

###  4.2 Structured controller design with partially-accessed states and noisy
measurement outputs

In the previous section, we designed a full-state feedback controller, meaning
that the control action taken by each generator also depends on all the rest
of the generators and network states. However, this could be unrealistic and
require a dense and reliable communication network. To that end, here, we
design an SSOF $\mathcal{H}_{\infty}$ controller and simultaneously consider
the scenario in which, we have noisy measurement outputs. In the proposed
methodology, various sparsity structures (such as decentralized and
distributed) on the controller matrix $F$ can be imposed using $S$ as
discussed in Procedure Procedure 1. Here, we design two types of controller
architecture: (i) completely decentralized design, meaning that each generator
utilizes only its own state information, and (ii) distributed design in which
the control action taken by each generator depends on a few other state
information. Furthermore, to consider noisy output measurements, we set
$D_{y}=0.1I_{n_{y},n_{w}}$.

Tab. 2 represents the $\mathcal{H}_{\infty}$ values and computational times
for NLA corresponding to SSOF $\mathcal{H}_{\infty}$ controllers for different
sparsity structures for the IEEE test systems. Such different sparsity
structures are three-fold as follows:

  * 1.

Centralized: in this case, we set $S$ as the following dense sparsity
structure:

$\displaystyle S$ $\displaystyle=\mathbf{1}_{n_{u}\times n_{y}}.$

  * 2.

Distributed: in this case, we set $S$ as the following randomly generated
zero-one sparsity structure:

$\displaystyle S$
$\displaystyle=\texttt{randi}(\begin{bmatrix}0&1\end{bmatrix},n_{u},n_{y}).$

where $\texttt{randi}(\begin{bmatrix}0&1\end{bmatrix},n_{u},n_{y})$ is a
MATLAB built-in function that generates a $n_{u}$ by $n_{y}$ matrix consisting
of uniformly distributed pseudorandom zero-one values.

  * 3.

Decentralized: in this case, noting that $(n_{u},n_{d})=(2N,4N)$ holds, we
choose $C_{y}=\begin{bmatrix}I_{n_{d}}&0\end{bmatrix}$ and set $S$ as the
following decentralized sparsity structure:

$\displaystyle S$ $\displaystyle=\begin{bmatrix}I_{N}&I_{N}&I_{N}&I_{N}\\\
I_{N}&I_{N}&I_{N}&I_{N}\end{bmatrix}=\mathbf{1}_{2\times 4}\otimes I_{N}.$

Notice that the smaller the $\mathcal{H}_{\infty}$ value is, the better the
controller performance is. We see from Tab. 2 that as sparsity of the
controller is promoted the corresponding $\mathcal{H}_{\infty}$ value
increases which makes sense because the controller loses some state
information. In other words, there exists a fundamental trade-off between the
$\mathcal{H}_{\infty}$ performance and the controller sparsity. Fig. 3 depicts
the randomly generated zero-one sparsity structure of the Distributed $39$-bus
$F$ reported in Tab. 2.

Figure 3: (Section 4.2) The randomly generated zero-one sparsity structure of
the Distributed $39$-bus $F$ reported in Tab. 2. The blue dots represent the
nonzero values of $F$.

The dynamic performance of the system is assessed under load and renewable
disturbance as discussed in the previous section. Note that the Lyapunov-based
approaches proposed in [12, 3, 39] are not capable of designing an SSOF
$\mathcal{H}_{\infty}$ controller, thus their results are not included here.
The frequencies of all the generators are presented in Fig. 4. We see that
after a large abrupt disturbance in load and renewable, NLA controller can
successfully keep the system transiently stable and synchronized. We can also
see that the non-Lyap controller can also bring the system back to its nominal
frequency after the large disturbance.

Table 2: (Section 4.2) The $\mathcal{H}_{\infty}$ values and computational times for NLA corresponding to SSOF $\mathcal{H}_{\infty}$ controllers with different choices of $S$: Centralized (Cen.), Distributed (Dis.), and Decentralized (Dec.), $C_{y}=I_{n_{y},n_{x}}$, $D_{y}=0.1I_{n_{y},n_{w}}$, $n_{y}=n_{d}$, and $B_{h}=0.1I_{n_{x}}$ for the IEEE test systems. Sparsity Structure | $\|T_{z\tilde{w}}(s)\|_{\mathcal{H}_{\infty}}$ | Computational Time | $(n_{u},n_{y})$
---|---|---|---
Cen. $9$-bus | $3.2314$ |  $1.69$ s | $(6,12)$
Dis. $9$-bus | $3.3007$ |  $4.86$ s | $(6,12)$
Dec. $9$-bus | $3.2784$ |  $4.11$ s | $(6,12)$
Cen. $14$-bus | $7.3964$ |  $2.04$ s | $(10,20)$
Dis. $14$-bus | $7.4446$ |  $5.40$ s | $(10,20)$
Dec. $14$-bus | $7.4274$ |  $1.78$ s | $(10,20)$
Cen. $39$-bus | $9.2120$ |  $42.59$ s | $(20,40)$
Dis. $39$-bus | $9.9968$ |  $145.53$ s | $(20,40)$
Dec. $39$-bus | $9.8467$ |  $82.38$ s | $(20,40)$
Cen. $57$-bus | $23.6002$ |  $35.33$ s | $(14,28)$
Dis. $57$-bus | $23.6168$ |  $51.04$ s | $(14,28)$
Dec. $57$-bus | $23.6180$ |  $56.93$ s | $(14,28)$

(a)  (b)

(c)  (d)

Figure 4: (Section 4.2) The dynamic performance of the proposed SSOF
$\mathcal{H}_{\infty}$ controller under large disturbance in load demand and
renewable power generation and noisy measurement states; Generator $1$
frequency and rotor angle for $9$-bus (top) and $57$-bus (bottom) test
systems.

###### Remark 2

Note that unlike the [12], we take into account a more general form for the
nonlinearity $h$ by choosing $B_{h}\propto I_{n_{x}}$ instead of $B_{h}\propto
B_{w}$ as in the former one, the number of columns of $B_{h}$ is $n_{x}$ while
in the latter one, the number of columns of $B_{h}$ is $n_{w}<n_{x}$. Then,
such a more general choice of $B_{h}$ makes the Lyapunov-based approach even
more inefficient in terms of scalability as it makes the LMI in the convex SDP
even larger by adding extra rows and columns.

Defining the following block matrix notations for $i,j\in\\{1,\dots,N\\}$:

$\displaystyle\hat{u}$
$\displaystyle=\begin{bmatrix}\hat{u}_{1}^{T}&\cdots&\hat{u}_{N}^{T}\end{bmatrix}^{T},\leavevmode\nobreak\
\hat{u}_{i}=\begin{bmatrix}u_{i}&u_{i+N}\end{bmatrix}^{T},\leavevmode\nobreak\
\hat{x}=\begin{bmatrix}\hat{x}_{1}^{T}&\cdots&\hat{x}_{N}^{T}\end{bmatrix}^{T},$
$\displaystyle\hat{x}_{j}$
$\displaystyle=\begin{bmatrix}x_{j}&x_{j+N}&x_{j+2N}&x_{j+3N}\end{bmatrix}^{T},$
$\displaystyle\hat{F}_{i,j}$
$\displaystyle=\begin{bmatrix}F_{i,j}&F_{i,j+N}&F_{i,j+2N}&F_{i,j+3N}\\\
F_{i+N,j}&F_{i+N,j+N}&F_{i+N,j+2N}&F_{i+N,j+3N}\end{bmatrix},$

we get $\hat{u}=\hat{F}\hat{x}$. Then, we can define the block-sparsity
structure of $\hat{F}$ as follows:

$\displaystyle\hat{S}_{i,j}$
$\displaystyle=\begin{cases}0,&\mathrm{if}\leavevmode\nobreak\
\hat{F}_{i,j}=0,\\\ 1&\mathrm{if}\leavevmode\nobreak\ \hat{F}_{i,j}\neq
0.\end{cases}$

Based on the definition of this block-sparsity structure, the decentralized
sparsity structure $S=\mathbf{1}_{2\times 4}\otimes I_{N}$ would attain the
block-sparsity structure $\hat{S}=I_{N}$. It is noteworthy that since
$\hat{S}_{i,i}=1$ requires $\hat{F}_{i,i}\neq 0$, then a total number of
$(2^{2\times 4}-1)^{N}=255^{N}$ sparsity structures $S$ can lead to the block-
sparsity structure $\hat{S}=I_{N}$. In this sense, the sparsity structure
$S=\mathbf{1}_{2\times 4}\otimes I_{N}$ is the densest decentralized sparsity
structure. Moreover, there exist a total number of $(2\times 4)^{N}=8^{N}$
sparsest decentralized sparsity structures (for each $i\in\\{1,\dots,N\\}$,
$\hat{F}_{i,i}$ has only $1$ non-zero element). Although investigating the
$255^{N}$ possibilities leading to $\hat{S}=I_{N}$ is computationally
impossible, we may consider the following subset of those sparsity structures:

$\displaystyle S=\tilde{S}\otimes I_{N},$ (24)

where $\tilde{S}\in\\{0,1\\}^{2\times 4}\backslash\\{0\\}$. Note that in this
specific subset, the sparsity structure of $\hat{F}$ would be
$I_{N}\otimes\tilde{S}$. Fig. 5 illustrates the $\mathcal{H}_{\infty}$ values
for all $255$ possible selections of $\tilde{S}$. As Fig. 5 shows, for the
IEEE test systems, the corresponding optimal selections of $\tilde{S}$ are as
follows:

$\displaystyle\tilde{S}^{\mathrm{opt}}_{9-\mathrm{bus}}$
$\displaystyle=\begin{bmatrix}0&1&1&1\\\
1&1&1&1\end{bmatrix},\leavevmode\nobreak\
\tilde{S}^{\mathrm{opt}}_{14-\mathrm{bus}}=\begin{bmatrix}0&1&1&1\\\
1&1&0&1\end{bmatrix},\leavevmode\nobreak\
\tilde{S}^{\mathrm{opt}}_{39-\mathrm{bus}}=\begin{bmatrix}1&0&1&0\\\
0&1&1&1\end{bmatrix},$
$\displaystyle\tilde{S}^{\mathrm{opt}}_{57-\mathrm{bus}}$
$\displaystyle=\begin{bmatrix}1&1&1&1\\\ 1&1&1&1\end{bmatrix}.$

with the $\mathcal{H}_{\infty}$ values of $3.2729$, $7.4251$, $9.6947$, and
$23.6180$, respectively. Comparing these values with the corresponding values
reflected in Tab. 2, we observe that except for the $57$-bus test system,
imposing the sparsity structure on $F$ can help the non-Lyapunov solver to
obtain a better solution in terms of $\mathcal{H}_{\infty}$ optimality. Note
that ideally the best solution should be associated with the densest sparsity
structure, i.e., $\tilde{S}=\mathbf{1}_{2\times 4}$, however, due to the sub-
optimal nature of the solutions of non-Lyapunov solver in practice
(highlighted by Remark 1), it may not be the case. Nevertheless, it suggests
that exploring the sparsity-constrained solutions may be beneficial to us in
practice via potentially detecting better optimal solutions in terms of
$\mathcal{H}_{\infty}$ optimality. Fig. 6 visualizes the relationship between
the sparsity level of $\tilde{S}$ and the average computational time (s). We
observe that by increasing the size of the test system, computing the sparsest
solutions takes less time compared to the densest solutions. We may interpret
it as the positive effect of imposing the sparsity structure on computational
time improvement for the larger test systems in which we deal with a larger
search space for solutions. In other words, imposing the sparsity structure
significantly reduces the search space in the case of large test systems.

(a)  (b)

(c)  (d)

Figure 5: (Section 4.2) The $\mathcal{H}_{\infty}$ values for all $255$
possible selections of $\tilde{S}$ for $9$-bus (top-left), $14$-bus (top-
right), $39$-bus (bottom-left), and $57$-bus (bottom-right) test systems. The
selections leading to the instability have been excluded.

(a)  (b)

(c)  (d)

Figure 6: (Section 4.2) The average computational time (ACT) (s) versus
$\|\tilde{S}\|_{0}$ for all $255$ possible selections of $\tilde{S}$ for
$9$-bus (top-left), $14$-bus (top-right), $39$-bus (bottom-left), and $57$-bus
(bottom-right) test systems.

###  4.3 Output matrix selection

####  4.3.1 The additionally-accessed algebraic states

In the previous section, we chose the output matrix $C_{y}$ as
$C_{y}=\begin{bmatrix}I_{n_{d}}&0\end{bmatrix}$. In other words, we assumed
that we only have access to the information of dynamic states, i.e.,
$x_{d}=\begin{bmatrix}\delta^{T}&\omega^{T}&E^{\prime
T}&T_{M}^{T}\end{bmatrix}^{T}$ which includes the generators’ internal states.
Here, we select a relatively general output matrix $C_{y}$ as
$C_{y}=\begin{bmatrix}I_{6N}&0\end{bmatrix}$ by additionally taking into
account the $2N$ first elements of algebraic states, i.e.,
$a=\begin{bmatrix}P_{G}^{T}&Q_{G}^{T}\end{bmatrix}^{T}$ which includes the
generators’ supplied power states.

Considering the following decentralized sparsity structure:

$\displaystyle S$
$\displaystyle=\begin{bmatrix}I_{N}&I_{N}&I_{N}&I_{N}&I_{N}&I_{N}\\\
I_{N}&I_{N}&I_{N}&I_{N}&I_{N}&I_{N}\end{bmatrix}=\mathbf{1}_{2\times 6}\otimes
I_{N},$

we obtain the results reflected in Tab. 3.

Table 3: (Section 4.3) The $\mathcal{H}_{\infty}$ values and computational times for NLA corresponding to SSOF $\mathcal{H}_{\infty}$ controllers with different choices of Decentralized (Dec.) $S$: $(n_{u},n_{y})=(2N,4N)$ and $(n_{u},n_{y})=(2N,6N)$, $C_{y}=I_{n_{y},n_{x}}$, $D_{y}=0.1I_{n_{y},n_{w}}$, and $B_{h}=0.1I_{n_{x}}$ for the IEEE test systems. Sparsity Structure | $\|T_{z\tilde{w}}(s)\|_{\mathcal{H}_{\infty}}$ | Computational Time | $(n_{u},n_{y})$
---|---|---|---
Dec. $9$-bus | $3.2784$ |  $4.11$ s | $(6,12)$
Dec. $9$-bus | $3.1946$ |  $2.68$ s | $(6,18)$
Dec. $14$-bus | $7.4274$ |  $1.78$ s | $(10,20)$
Dec. $14$-bus | $7.4273$ |  $1.65$ s | $(10,30)$
Dec. $39$-bus | $9.8467$ |  $82.38$ s | $(20,40)$
Dec. $39$-bus | $9.4424$ |  $32.17$ s | $(20,60)$
Dec. $57$-bus | $23.6180$ |  $56.93$ s | $(14,28)$
Dec. $57$-bus | $23.5797$ |  $21.21$ s | $(14,42)$

Tab. 3 demonstrates that utilizing the information of generators’ supplied
power states, i.e., $2N$ extra state measurements, both $\mathcal{H}_{\infty}$
norm values and computational times can significantly be improved.

####  4.3.2 The partially-accessed dynamic states

We can consider the scenario in which, we want to optimally select and place
$n_{y}=p$ states out of $n_{d}=4N$ dynamic states in the
$\mathcal{H}_{\infty}$ sense. To this end, there exist
$\frac{n_{d}!}{n_{y}!(n_{d}-n_{y})!}$ possibilities in the case of
$S=\mathbf{1}_{n_{u}\times n_{y}}$. Note that according to the existence of
permutation matrix transformations detailed later on, we do not take into
account the $n_{y}!$ permutations of each selection when
$S=\mathbf{1}_{n_{u}\times n_{y}}$ holds. Then, denoting the positions of non-
zero elements of $C_{y}$ as $(1,j_{1})$ to $(n_{y},j_{n_{y}})$ at which each
element of $C_{y}$ is set to $1$, without loss of generality, we can assume
that the ascending order $j_{1}<\dots<j_{n_{y}}$ holds for the case of
$S=\mathbf{1}_{n_{u}\times n_{y}}$. Assuming that $F$ and $F^{\pi}$ denote the
solutions corresponding to $(j_{1},\dots,j_{n_{y}})$ and
$(j_{\pi(1)},\dots,j_{\pi(n_{y})})$, respectively where
$\pi:\\{1,\dots,n_{y}\\}\longrightarrow\\{1,\dots,n_{y}\\}$ represents a
permutation mapping, we have

$\displaystyle F^{\pi}$ $\displaystyle=FP_{\pi}^{T},\leavevmode\nobreak\
P_{\pi}^{T}=\begin{bmatrix}\mathbf{e}_{\pi(1)}&\cdots&\mathbf{e}_{\pi(n_{y})}\end{bmatrix}.$

However, for the case of $S=I_{n_{u}\times n_{y}}$, we need to take into
account all the permutations as the decentralized sparsity structure on the
controller is not preserved under the permutation matrix transformation. Then,
the total number of possible selections in such a scenario would be
$\frac{n_{d}!}{(n_{d}-n_{y})!}$. Since investigating all
$\frac{n_{d}!}{(n_{d}-n_{y})!}$ possibilities is highly time-consuming,
similar to the case of $S=\mathbf{1}_{n_{u}\times n_{y}}$, we only investigate
the $\frac{n_{d}!}{n_{y}!(n_{d}-n_{y})!}$ possibilities in ascending order.
Setting $n_{y}=p=2N$, Fig. 7 depicts the $\mathcal{H}_{\infty}$ values for all
$\frac{n_{d}!}{n_{y}!(n_{d}-n_{y})!}$ possible selections of $C_{y}$ for
$9$-bus test system including both centralized $S=\mathbf{1}_{n_{u}\times
n_{y}}$ and decentralized $S=I_{n_{u}\times n_{y}}$ scenarios. As Fig. 7
demonstrates, the optimal selection for the centralized scenario is associated
with $i=453$, i.e., $(j_{1},\dots,j_{6})=(1,6,8,9,10,12)$ and the
corresponding $\mathcal{H}_{\infty}$ value is $3.2579$. Also, the optimal
selection for the decentralized scenario is associated with $i=19$, i.e.,
$(j_{1},\dots,j_{6})=(1,2,3,4,8,9)$ and the corresponding
$\mathcal{H}_{\infty}$ value is $5.6785$. Fig. 7 interestingly shows that the
selections after $i=840$, i.e., the selections with $j_{1}>3$ cannot stabilize
the system. In other words, if none of $(x_{1},x_{2},x_{3})$ is selected,
there is no stabilizing $F$. The physical interpretation of such an
observation is that the state information of $\delta$ of at least one of the
generators is required to stabilize the system.

In Fig. 8, for $9$-bus test system and $S=I_{n_{u}\times n_{y}}$, we explore
the $n_{y}!$ permutations of the optimal selection among the
$\frac{n_{d}!}{n_{y}!(n_{d}-n_{y})!}$ possibilities in ascending order. Also,
we explore all $\frac{n_{d}!}{n_{y}!(n_{d}-n_{y})!}$ possible selections of
$C_{y}$ subject to the permutation order associated with the optimal
selection. Interestingly, both explorations suggest
$(j_{1},\dots,j_{6})=(3,8,9,4,2,1)$ with $\mathcal{H}_{\infty}$ value of
$4.1145$ as the optimal placement for $C_{y}$. Then, we realize that exploring
the $n_{y}!$ permutations of the optimal selection for $C_{y}$, i.e., $720$
possibilities associated with $6$ optimally selected states
$\\{1,2,3,4,8,9\\}$, significantly improves the $\mathcal{H}_{\infty}$ value
from $5.6785$ to $4.1145$.

(a)  (b)

Figure 7: (Section 4.3) The $\mathcal{H}_{\infty}$ values for all
$\frac{n_{d}!}{n_{y}!(n_{d}-n_{y})!}$ possible selections of $C_{y}$ for
$9$-bus test system: Centralized, $S=\mathbf{1}_{n_{u}\times n_{y}}$ (left),
and Decentralized, $S=I_{n_{u}\times n_{y}}$ (right). The selections leading
to the instability have been excluded.

(a)  (b)

Figure 8: (Section 4.3) The $\mathcal{H}_{\infty}$ values for the $n_{y}!$
permutations of the optimal selection of $C_{y}$ among the
$\frac{n_{d}!}{n_{y}!(n_{d}-n_{y})!}$ possibilities in ascending order for
$9$-bus test system, Decentralized, $S=I_{n_{u}\times n_{y}}$ (left), the
$\mathcal{H}_{\infty}$ values for all $\frac{n_{d}!}{n_{y}!(n_{d}-n_{y})!}$
possible selections of $C_{y}$ subject to the permutation order associated
with the optimal selection for $9$-bus test system, Decentralized,
$S=I_{n_{u}\times n_{y}}$ (right (the selections leading to the instability
have been excluded)).

###  4.4 Dense estimator design

In this section, we design a dense $\mathcal{H}_{\infty}$ estimator via
setting $S=\mathbf{1}_{n_{x}\times n_{y}}$. Considering the IEEE test systems,
we utilize a convex SDP similar to the one proposed in [8] and Procedure
Procedure 2 to obtain the Lyapunov-based and NLA dense $\mathcal{H}_{\infty}$
estimators, respectively. In order to make a fair comparison between the
Lyapunov-based approach and NLA, we slightly modify the convex SDP proposed by
[8] provided that $\Delta f$ is $\mathcal{L}_{2}$-norm bounded (instead of
Lipschitz boundedness assumption utilized in [8]). Then, we utilize a modified
version of (10) as follows:

$\displaystyle{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\dot{V}(e)+e^{T}e-\lambda\tilde{w}^{T}\tilde{w}<0.}$

Tab. 4 lists the $\mathcal{H}_{\infty}$ values and computational times for the
Lyapunov-based approach and NLA corresponding to the dense
$\mathcal{H}_{\infty}$ estimator with $S=\mathbf{1}_{n_{x}\times n_{y}}$ and
$B_{\Delta f}=0.1I_{n_{x}}$ for the IEEE test systems. As Tab. 4 indicates,
NLA proposes better solutions in terms of $\mathcal{H}_{\infty}$ values
compared to its Lyapunov counterpart. However, due to the extremely large
number of variables in NLA (e.g., $12480$ for the $57$-bus test system) and
the simpler (sparser) structure of the LMIs in the Lyapunov-based approach,
NLA is highly time-consuming, in the case of estimator design, unlike the
controller design case.

Tab. 5 represents the $\mathcal{H}_{\infty}$ values and computational times
for the Lyapunov-based approach and NLA corresponding to the dense
$\mathcal{H}_{\infty}$ estimator with $S=\mathbf{1}_{n_{x}\times n_{y}}$,
$B_{\Delta f}=0.1I_{n_{x}}$, $C_{y}=I_{n_{y},n_{x}}$, and $n_{y}=2N$ for the
IEEE test systems. Comparing the results of Tab. 5 with the results of Tab. 4
reveals the fact that the scalability of NLA can significantly be improved for
the case of having less number of output measurements (notice that the
scalability superiority of the Lyapunov-based approach compared to NLA gets
reversed). For instance, for the $57$-bus test system, reducing the $n_{y}$
from $80$ to $14$, the computational time reduces from more than $8$ hours to
less than $1$ minute. Then, minimizing the number of output measurements can
play a significant role in terms of scalability improvement while attaining a
reasonable $\mathcal{H}_{\infty}$ performance degradation. As another
observation, we realize that the $\mathcal{H}_{\infty}$ performance of the
Lyapunov-based approach can significantly be improved for the case of having a
smaller number of output measurements. For instance, except for the $14$-bus
test system, the Lyapunov-based approach attains the same
$\mathcal{H}_{\infty}$ performance as NLA. Although the computational time is
also improved, it is not significant compared to the scalability improvement
in the case of NLA.

Table 4: (Section 4.4) The $\mathcal{H}_{\infty}$ values and computational times for the Lyapunov-based approach and NLA corresponding to the Dense $\mathcal{H}_{\infty}$ estimator with $S=\mathbf{1}_{n_{x}\times n_{y}}$ and $B_{\Delta f}=0.1I_{n_{x}}$ for the IEEE test systems. Approach | $\|T_{z\tilde{w}}(s)\|_{\mathcal{H}_{\infty}}$ | Computational Time | $(n_{x},n_{y})$
---|---|---|---
non-Lyap $9$-bus | $1.5100$ |  $4.76$ s | (36,14)
Lyap $9$-bus | $1.5100$ | $2.26$ s | $(36,14)$
non-Lyap $14$-bus | $1.4140$ |  $44.30$ s | $(58,28)$
Lyap $14$-bus | $2.6656$ | $7.58$ s | $(58,28)$
non-Lyap $39$-bus | $3.1129$ |  $6841.34$ s | $(138,58)$
Lyap $39$-bus | $4.7038$ | $943.74$ s | $(138,58)$
non-Lyap $57$-bus | $2.1155$ |  $29639.22$ s | $(156,80)$
Lyap $57$-bus | $3.6180$ |  $4113.78$ s | $(156,80)$
Table 5: (Section 4.4) The $\mathcal{H}_{\infty}$ values and computational times for the Lyapunov-based approach and NLA corresponding to the Dense $\mathcal{H}_{\infty}$ estimator with $S=\mathbf{1}_{n_{x}\times n_{y}}$, $B_{\Delta f}=0.1I_{n_{x}}$, $C_{y}=I_{n_{y},n_{x}}$, and $n_{y}=2N$ for the IEEE test systems. Approach | $\|T_{z\tilde{w}}(s)\|_{\mathcal{H}_{\infty}}$ | Computational Time | $(n_{x},n_{y})$
---|---|---|---
non-Lyap $9$-bus | $3.0921$ |  $0.85$ s | (36,6)
Lyap $9$-bus | $3.0921$ |  $6.06$ s | $(36,6)$
non-Lyap $14$-bus | $3.2893$ |  $1.95$ s | $(58,10)$
Lyap $14$-bus | $4.3332$ |  $7.25$ s | $(58,10)$
non-Lyap $39$-bus | $2.6402$ |  $45.25$ s | $(138,20)$
Lyap $39$-bus | $2.6402$ |  $556.77$ s | $(138,20)$
non-Lyap $57$-bus | $17.5631$ |  $17.28$ s | $(156,14)$
Lyap $57$-bus | $17.5631$ |  $1337.05$ s | $(156,14)$

###  4.5 Structured estimator design

In this section, we design a structured $\mathcal{H}_{\infty}$ estimator via
imposing the sparsity structure $S$ on the estimator. Considering the IEEE
test systems, we utilize Procedure Procedure 2 to obtain NLA structured
$\mathcal{H}_{\infty}$ estimators.

Tab. 6 summarizes the $\mathcal{H}_{\infty}$ values and computational times
for NLA corresponding to the dense $\mathcal{H}_{\infty}$ estimator with
$S=\mathbf{1}_{n_{x}\times n_{y}}$ and the structured $\mathcal{H}_{\infty}$
estimator with $S=I_{n_{x},n_{y}}$ under $B_{\Delta f}=0.1I_{n_{x}}$ for the
IEEE test systems. It shows that imposing the sparsity structure on $L$ can
significantly expedite the estimator computation at the cost of
$\mathcal{H}_{\infty}$ performance degradation. Such an observation is due to
the fact that imposing the sparsity structure on $L$ significantly reduces the
number of variables (e.g., the reduction from $12480$ to $80$ for the $57$-bus
test system leading to the reduction from more than $8$ hours to $5$ minutes).
Comparing the results of Tab. 6 with the results of Tab. 4, suggests that for
the large-scale power systems, utilizing the structured $\mathcal{H}_{\infty}$
estimator design (obtained by NLA) can significantly save computational time
compared to the dense $\mathcal{H}_{\infty}$ estimator design (obtained by the
Lyapunov-based approach) while degrading the $\mathcal{H}_{\infty}$ by some
extent. Moreover, such performance degradation can potentially be minimized by
optimally choosing the output matrix $C_{y}$.

Table 6: (Section 4.5) The $\mathcal{H}_{\infty}$ values and computational times for NLA corresponding to structured $\mathcal{H}_{\infty}$ estimators with different choices of $S$: Dense ($S=\mathbf{1}_{n_{x}\times n_{y}}$) and Structured ($S=I_{n_{x},n_{y}}$) under $B_{\Delta f}=0.1I_{n_{x}}$ for the IEEE test systems. Sparsity Structure | $\|T_{z\tilde{w}}(s)\|_{\mathcal{H}_{\infty}}$ | Computational Time | $(n_{x},n_{y})$
---|---|---|---
Dense $9$-bus | $1.5100$ |  $4.76$ s | $(36,14)$
Structured $9$-bus | $3.4120$ |  $4.63$ s | $(36,14)$
Dense $14$-bus | $1.4140$ |  $44.30$ s | $(58,28)$
Structured $14$-bus | $8.7075$ |  $11.63$ s | $(58,28)$
Dense $39$-bus | $3.1129$ |  $6841.34$ s | $(138,58)$
Structured $39$-bus | $13.3505$ |  $223.11$ s | $(138,58)$
Dense $57$-bus | $2.1155$ |  $29639.22$ s | $(156,80)$
Structured $57$-bus | $5.5907$ |  $299.10$ s | $(156,80)$

###  4.6 Estimator performance under transient conditions

In this section, we discuss the performance of the estimator in estimating all
the states of the power system under various transient conditions and noisy
measurements. Notice that for all the case studies here, we assume that phasor
measurement units (PMUs) are already placed optimally in the power system and
thus, making it completely observable as given in [40, 41, 7]. It is
noteworthy that the proposed estimator has the following key advantages as
compared to the current literature on power system dynamic state estimation
(DSE):

  * 1.

The proposed estimator in this study can simultaneously estimate both dynamic
states (states of generators) and algebraic states (states of the network such
as voltages and currents). In the current literature on power system DSE, they
are usually estimated separately [42, 43] because of the complexity of
handling the complete power system NDAE models. Few studies have been carried
out to estimate them simultaneously [8, 44]. However, they use Lyapunov-based
approaches which are much more difficult to solve for a large-scale power
system model as discussed in this paper.

  * 2.

The proposed estimator does not require any statistical properties of the
disturbance/noise and can provide accurate estimation results as long as the
disturbance is norm-bounded. It can also seamlessly handle the situation when
the real-time control inputs are not known to the estimator and only steady-
state values are given.

  * 3.

The presented estimator also only requires a few measurements from PMUs placed
optimally such that the whole system is observable as compared to the
literature where it is commonly required that all the generator buses need to
be equipped with PMUs [45].

For all the case studies, the estimator dynamics (19) are initialized from
random initial conditions having $20\%$ maximum deviation from the steady-
state values of the power system. With that in mind, we consider the following
various scenarios:

####  4.6.1 Estimation under disturbance in load demand and renewable

Here, we consider that the minute or hour ahead prediction (or steady-state
values) of load $P_{L}$, $Q_{L}$ and renewable generation $P_{R}$, $Q_{R}$ are
known (this is realistic as system operators obtain and publish these
quantities regularly—see[46]) while the disturbance/uncertainty in them is not
known to the estimator. To that end, we add disturbance in load and renewable
as discussed in Section 4.1, and then we estimate the state of the power
system using both dense and structured estimators (with only steady-state
values of load and renewable known to them) as designed in the previous
sections. The results are shown in Fig. 9 where we see that the proposed
estimator can successfully estimate all the states of the power system with
appropriate accuracy. This can be verified from the error norm also given in
Fig. 10 where we observe that it asymptotically converges to zero for all test
networks.

(a)  (b)

(c)  (d)

Figure 9: (Section 4.6) Estimation results for $57$-bus test system under load
and renewable disturbance; Generator $1$ transient voltage and rotor angle and
Bus $7$ voltage and angle.

(a)  (b)

(c)  (d)

Figure 10: (Section 4.6) Estimation error norm for $9$-bus, $14$-bus,
$39$-bus, and $57$-bus test systems under disturbance in load demand and
renewable power generation.

####  4.6.2 Estimation under non-Gaussian measurement noise

Here, we consider the scenario when the measurement $y$ received by the
estimator contains some non-Gaussian noise. Notice that it is important to
check the performance of the estimator under non-Gaussian noise because the
measurement noise is usually non-Gaussian and assuming it to be Gaussian is a
serious simplification as discussed in [47]. We want to mention that the
vintage Kalman-based estimators (such as EKF, UKF, etc) that are commonly used
in power system DSE cannot handle non-Gaussian noise as they require some
statistical properties of the noise to work. This is one of the main
advantages of the proposed estimator as it can handle any type of bounded
noise because of the robust $\mathcal{H}_{\infty}$-notion used to design it.
To that end, Cauchy noise has been generated as $w_{mi}=a+b(\pi(R-0.5))$,
where $a=0$, $b=1\times 10^{-5}$ and $R$ is a random variable inside $(0,1)$
and has been added to the PMU measurements. The estimation results are shown
in Fig. 11 where we see that the proposed estimator can still accurately
estimate all the states of the power system.

(a)  (b)

(c)  (d)

Figure 11: (Section 4.6) Estimation error norm for $9$-bus (top-left),
$14$-bus (top-right), $39$-bus (bottom-left), and $57$-bus (bottom-right) test
systems under non-Gaussian measurement noise.

####  4.6.3 Estimation under unknown control inputs

Here, we consider the scenario when the real-time control input $u(t)$ is not
available to the estimator and only steady-state values of control input $u$
are known. Notice that this is important to consider because, in the brushless
excitation system of synchronous generators, it is difficult to measure field
current and voltage in real-time [48]. Thus, assuming that the real-time
information of control inputs is known to the estimator is unrealistic. To
that end only steady-state values of $u(t)$ are supplied to the estimator
dynamics (19) and state estimation is performed. The results are shown in Fig.
12 where we observe that the proposed estimator can still provide appropriate
estimation results. This can also be verified from Fig. 13 where we see that
the trajectory of error norm between actual $x$ and estimated states $\hat{x}$
successfully converges to zero.

(a)  (b)

(c)  (d)

Figure 12: (Section 4.6) Estimation results for $57$-bus test system under
unknown inputs; Generator $4$ speed, rotor angle, generated power, and
terminal bus voltage.

(a)  (b)

Figure 13: (Section 4.6) Estimation error norm for $14$-bus and $57$-bus test
systems under unknown control inputs.

##  5 Concluding Remarks

In this paper, we employ a celebrated non-Lyapunov approach (NLA) from the
control theory literature to design scaled-up structured
$\mathcal{H}_{\infty}$ controllers and estimators. Through extensive numerical
experiments on the IEEE test systems, we empirically validate that such NLA
significantly improves the scalability and sensitivity simultaneously compared
to its Lyapunov-based counterparts. Moreover, in the case of controller
design, it can handle the cases in which we deal with partially-accessed
states and noisy measurement outputs which are more realistic compared to the
fully-accessed states and noiseless measurement states considered by the
Lyapunov-based approach. Also, in the case of estimator design, the proposed
estimator can deal with load demand and renewable disturbances, non-Gaussian
measurement noise, and unknown control inputs. The employed NLA in this paper
can readily be applied to various robust controller/estimator algorithms in
power systems to overcome the computational scalability and computational
sensitivity issues arising from the Lyapunov-based LMI formulations.
Nonetheless, NLA also has some limitations. One limitation is that NLA becomes
impractical without a linear time-invariant (LTI) model. In other words, one
should have an LTI model or derive an approximate LTI model before employing
NLA. As another limitation, NLA cannot directly handle the continuous-time LTI
models with delay while it can handle the discrete-time LTI models with delay.
Furthermore, depending on the stability properties of the system, NLA can fail
to stabilize the closed-loop system and may need an appropriate selection of
an initialization for its built-in non-convex non-smooth optimization
algorithm. Fortunately, NLA can successfully stabilize the closed-loop system
in the case of power systems considered in this paper (by selecting $0$ as an
initialization). Thus, in the current paper, we do not face the last expressed
limitation of NLA.

The answers to the itemized questions in Section 4 (Q1–Q5) are as follows:

  * 1.

A1: The extensive numerical simulations (Tabs. 1 and 5) reveal that the
structured $\mathcal{H}_{\infty}$ designs can be computed via NLA in an
extremely reduced computation time compared to the Lyapunov-based approach
while attaining the same $\mathcal{H}_{\infty}$ performance. For instance, for
the $57$-bus test system, NLA controller and estimator are $47.17$ times and
$76.38$ times faster than the Lyapunov-based counterparts.

  * 2.

A2: For the dense controller design, as $\frac{n_{u}}{n_{x}}$ becomes larger,
the scalability quality of NLA deteriorates. Similarly, as
$\frac{n_{y}}{n_{x}}$ increases, the scalability quality of NLA degrades (as
illustrated by Tab. 4). As Tab. 4 shows, the superiority of NLA over the
Lyapunov-based counterpart has been reversed due to the large value of
$\frac{n_{y}}{n_{x}}$. In Tab. 5, the small value of $\frac{n_{y}}{n_{x}}$
ensures the superiority of NLA over the Lyapunov-based counterpart. An
intuitive justification for such an observation could be that the smaller the
number of variables, the smaller the search space for the solutions.

  * 3.

A3: Unlike the Lyapunov-based approach, NLA can be utilized to propose
$\mathcal{H}_{\infty}$ controller designs subject to the imposed
structural/sparsity constraints (e.g., distributed and decentralized) with
partially-accessed states in the presence of noisy measurement outputs (Tab. 2
and Fig. 4). Tab. 2 reflects a fundamental trade-off between the
$\mathcal{H}_{\infty}$ performance and the controller sparsity. Defining the
block-sparsity structure and comparing the results demonstrated by Fig. 5 and
Tab. 2, we realize that imposing the more sparse structures may enable us to
obtain the more sparse structured designs with better $\mathcal{H}_{\infty}$
performance. A possible justification behind such a counterintuitive
observation could be that the smaller the number of nonzero elements
(variables), the smaller the search space for the solutions. Similarly, NLA
can be utilized to propose structured $\mathcal{H}_{\infty}$ estimator
designs. Tab. 6 demonstrates that by imposing the sparsity structure to the
estimator, computational time can significantly be improved at the cost of
extra $\mathcal{H}_{\infty}$ performance degradation.

  * 4.

A4: Tab. 3 depicts that by taking advantage of the information of additional
algebraic states (generators’ supplied power states) in addition to the
dynamics states (generators’ internal states), both $\mathcal{H}_{\infty}$
performance and computational time can significantly be improved. Also, Figs.
7 and 8 visualize the process of optimal selection/placement of $n_{y}=p$
states out of $n_{d}=4N$ dynamic states in the $\mathcal{H}_{\infty}$ sense.
As an interesting observation, we realize that the state information of
$\delta$ of at least one of the generators is required to stabilize the
system.

  * 5.

A5: Figs. 9, 10, 11, 12, and 13 highlight that NLA can be utilized to propose
estimators capable of effectively dealing with load demand and renewable
disturbances, non-Gaussian measurement noise, and unknown control inputs. The
proposed estimator in this study can simultaneously estimate both dynamic
states (states of generators) and algebraic states (states of the network such
as voltages and currents). In the current literature on power system DSE, they
are usually estimated separately because of the complexity of handling the
complete power system NDAE models. The proposed estimator does not require any
statistical properties of the disturbance/noise and can provide accurate
estimation results as long as the disturbance is norm-bounded. The presented
estimator also only requires a few measurements from PMUs placed optimally
such that the whole system is observable as compared to the literature where
it is commonly required that all the generator buses need to be equipped with
PMUs.

## References

  * [1] A. Kanchanaharuthai, P. Ngamsom, Robust $\mathcal{H}_{\infty}$ load-frequency control for interconnected power systems with D-stability constraints via LMI approach, in: Proceedings of American Control Conference., IEEE, 2005, pp. 4387–4392.
  * [2] N. Chuang, Robust load-frequency control in interconnected power systems, IET Control Theory & Applications 10 (1) (2016) 67–75.
  * [3] S. A. Nugroho, A. F. Taha, Load-and renewable-following control of linearization-free differential algebraic equation power system models, IEEE Transactions on Control Systems Technology 31 (4) (2023) 1774–1786.
  * [4] M. Nadeem, M. Bahavarnia, A. F. Taha, On wide-area control of solar-integrated dae models of power grids, in: American Control Conference (ACC), IEEE, 2023, pp. 4495–4500.
  * [5] J. Zhao, Dynamic state estimation with model uncertainties using ${H}_{\infty}$ extended Kalman filter, IEEE Transactions on Power Systems 33 (1) (2017) 1099–1100.
  * [6] J. Zhao, L. Mili, A theoretical framework of robust H-infinity unscented Kalman filter and its application to power system dynamic state estimation, IEEE Transactions on Signal Processing 67 (10) (2019) 2734–2746.
  * [7] M. Nadeem, S. A. Nugroho, A. F. Taha, Dynamic state estimation of nonlinear differential algebraic equation models of power networks, IEEE Transactions on Power Systems 38 (3) (2022) 2539–2552.
  * [8] M. Nadeem, A. F. Taha, Robust dynamic state estimation of multi-machine power networks with solar farms and dynamics loads, in: 61st Conference on Decision and Control (CDC), IEEE, 2022, pp. 7174–7179.
  * [9] S. Boyd, L. El Ghaoui, E. Feron, V. Balakrishnan, Linear matrix inequalities in system and control theory, SIAM, 1994.
  * [10] S. Gusev, A. Likhtarnikov, Kalman-Popov-Yakubovich lemma and the S-procedure: A historical essay, Automation & Remote Control 67 (11) (2006).
  * [11] K. Zhou, J. C. Doyle, Essentials of robust control, Vol. 104, Prentice Hall, 1998\.
  * [12] M. Nadeem, M. Bahavarnia, A. F. Taha, Robust feedback control of power systems with solar plants and composite loads, IEEE Transactions on Power Systems (2023).  doi:10.1109/TPWRS.2023.3323222.
  * [13] R. Majumder, B. Chaudhuri, H. El-Zobaidi, B. Pal, I. Jaimoukha, LMI approach to normalised $\mathcal{H}_{\infty}$ loop-shaping design of power system damping controllers, IEE Proceedings-Generation, Transmission and Distribution 152 (6) (2005) 952–960.
  * [14] E. D. Andersen, K. D. Andersen, The mosek interior point optimizer for linear programming: an implementation of the homogeneous algorithm, High performance optimization (2000) 197–232.
  * [15] J. Lofberg, Yalmip: A toolbox for modeling and optimization in matlab, in: IEEE international conference on robotics and automation (IEEE Cat. No. 04CH37508), 2004, pp. 284–289.
  * [16] A. Majumdar, G. Hall, A. A. Ahmadi, Recent scalability improvements for semidefinite programming with applications in machine learning, control, and robotics, Annual Review of Control, Robotics, and Autonomous Systems 3 (2020) 331–360.
  * [17] L. Yang, D. Sun, K.-C. Toh, Sdpnal+: a majorized semismooth newton-cg augmented lagrangian method for semidefinite programming with nonnegative constraints, Mathematical Programming Computation 7 (3) (2015) 331–366.
  * [18] D. Sun, K.-C. Toh, Y. Yuan, X.-Y. Zhao, Sdpnal+: A matlab software for semidefinite programming with bound constraints (version 1.0), Optimization Methods and Software 35 (1) (2020) 87–115.
  * [19] S. A. Nugroho, A. F. Taha, C. G. Claudel, A control-theoretic approach for scalable and robust traffic density estimation using convex optimization, IEEE Transactions on Intelligent Transportation Systems 22 (1) (2021) 64–78.
  * [20] K. Takaba, N. Morihira, T. Katayama, ${H}_{\infty}$ control for descriptor systems-a J-spectral factorization approach, in: Proceedings of 33rd IEEE Conference on Decision and Control, Vol. 3, 1994, pp. 2251–2256.
  * [21] A. Kawamoto, T. Katayama, Standard ${H}_{\infty}$ control problem for descriptor system, in: Proceedings of the 36th IEEE Conference on Decision and Control, Vol. 4, 1997, pp. 4130–4133.
  * [22] Y. Feng, M. Yagoubi, Robust control of linear descriptor systems, Springer, 2017\.
  * [23] P. Apkarian, D. Noll, Nonsmooth ${H}^{\infty}$ synthesis, IEEE Transactions on Automatic Control 51 (1) (2006) 71–86.
  * [24] F. H. Clarke, Optimization and nonsmooth analysis, SIAM, 1990.
  * [25] P. Gahinet, P. Apkarian, Structured $\mathcal{H}_{\infty}$ synthesis in matlab, IFAC Proceedings Volumes 44 (1) (2011) 1435–1440.
  * [26] P. Gahinet, P. Apkarian, Decentralized and fixed-structure $\mathcal{H}_{\infty}$ control in matlab, in: 50th IEEE Conference on Decision and Control and European Control Conference, 2011, pp. 8205–8210.
  * [27] P. Gahinet, P. Apkarian, Frequency-domain tuning of fixed-structure control systems, in: Proceedings of UKACC International Conference on Control, IEEE, 2012, pp. 178–183.
  * [28] H. Werner, P. Korba, T. C. Yang, Robust tuning of power system stabilizers using LMI-techniques, IEEE Transactions on Control Systems Technology 11 (1) (2003) 147–152.
  * [29] A. Basati, J. M. Guerrero, J. C. Vasquez, A. Fakharian, K. H. Johansson, S. Golestan, Robust internal model-based voltage control for DC microgrids: An LMI based ${H}_{\infty}$ control, Sustainable Energy, Grids and Networks (2023) 101094.
  * [30] A. F. Taha, M. Bazrafshan, S. A. Nugroho, N. Gatsis, J. Qi, Robust control for renewable-integrated power networks considering input bound constraints and worst case uncertainty measure, IEEE Transactions on Control of Network Systems 6 (3) (2019) 1210–1222.
  * [31] P. W. Sauer, J. H. Chow, Power system dynamics and stability: with synchrophasor measurement and power system toolbox, John Wiley & Sons, 2017.
  * [32] I. Masubuchi, Y. Kamitane, A. Ohara, N. Suda, $\mathcal{H}_{\infty}$ control for descriptor systems: A matrix inequalities approach, Automatica 33 (4) (1997) 669–673.
  * [33] K. Derinkuyu, M. Ç. Pınar, On the S-procedure and some variants, Mathematical Methods of Operations Research 64 (2006) 55–77.
  * [34] S. Boyd, V. Balakrishnan, P. Kabamba, A bisection method for computing the ${H}_{\infty}$ norm of a transfer matrix and related problems, Mathematics of Control, Signals and Systems 2 (3) (1989) 207–219.
  * [35] S. Boyd, V. Balakrishnan, A regularity result for the singular values of a transfer matrix and a quadratically convergent algorithm for computing its ${L}_{\infty}$-norm, Systems & Control Letters 15 (1) (1990) 1–7.
  * [36] P. Gahinet, P. Apkarian, Numerical computation of the ${L}_{\infty}$ revisited, in: Proceedings of the 31st IEEE Conference on Decision and Control, 1992, pp. 2257–2258.
  * [37] S. P. Boyd, C. H. Barratt, Linear controller design: limits of performance, Vol. 7, Citeseer, 1991.
  * [38] M. Grant, S. Boyd, Cvx: Matlab software for disciplined convex programming, version 2.1 (2014).
  * [39] T. Sadamoto, A. Chakrabortty, T. Ishizaki, J.-i. Imura, Dynamic modeling, stability, and control of power systems with distributed energy resources: Handling faults using two control methods in tandem, IEEE Control Systems Magazine 39 (2) (2019) 34–65.
  * [40] S. Chakrabarti, E. Kyriakides, Optimal placement of phasor measurement units for power system observability, IEEE Transactions on Power Systems 23 (3) (2008) 1433–1440.  doi:10.1109/TPWRS.2008.922621.
  * [41] N. M. Manousakis, G. N. Korres, A weighted least squares algorithm for optimal pmu placement, IEEE Transactions on Power Systems 28 (3) (2013) 3499–3500.  doi:10.1109/TPWRS.2013.2242698.
  * [42] J. Zhang, G. Welch, G. Bishop, Z. Huang, A two-stage Kalman filter approach for robust and real-time power system state estimation, IEEE Transactions on Sustainable Energy 5 (2) (2014) 629–636.  doi:10.1109/TSTE.2013.2280246.
  * [43] A. Rouhani, A. Abur, Linear phasor estimator assisted dynamic state estimation, IEEE Transactions on Smart Grid 9 (1) (2018) 211–219.  doi:10.1109/TSG.2016.2548244.
  * [44] S. A. Nugroho, A. F. Taha, N. Gatsis, J. Zhao, Observers for differential algebraic equation models of power networks: Jointly estimating dynamic and algebraic states, IEEE Transactions on Control of Network Systems 9 (3) (2022) 1531–1543.  doi:10.1109/TCNS.2022.3158754.
  * [45] J. Zhao, Z. Zheng, S. Wang, R. Huang, T. Bi, L. Mili, Z. Huang, Correlation-aided robust decentralized dynamic state estimation of power systems with unknown control inputs, IEEE Transactions on Power Systems 35 (3) (2020) 2443–2451.  doi:10.1109/TPWRS.2019.2953256.
  * [46] CAISO, Real time demand curves (Nov 2021).
URL http://www.caiso.com/outlook/SystemStatus.html

  * [47] S. Wang, J. Zhao, Z. Huang, R. Diao, Assessing gaussian assumption of pmu measurement error using field data, IEEE Transactions on Power Delivery 33 (6) (2018) 3233–3236.  doi:10.1109/TPWRD.2017.2762927.
  * [48] G. Anagnostou, B. C. Pal, Derivative-free Kalman filtering based approaches to dynamic state estimation for power systems with unknown inputs, IEEE Transactions on Power Systems 33 (1) (2018) 116–130.  doi:10.1109/TPWRS.2017.2663107.

##  Appendix A The Lyapunov-based approach to design dense controllers

The dense static state feedback (DSSF) $\mathcal{H}_{\infty}$ controller
proposed by [12] is obtained as follows:

$\displaystyle\boxed{F=H(XE^{T}+E^{\perp}W)^{-1}},$ (25)

for which the following convex SDP is solved accordingly:

$\displaystyle\underset{\lambda,H,X,W}{\min}$
$\displaystyle\leavevmode\nobreak\ \lambda$ $\displaystyle\mathrm{s.t.}$
$\displaystyle\leavevmode\nobreak\ X\succ
0,\lambda>0,\begin{bmatrix}\mathbf{He}(A(XE^{T}+E^{\perp}W)+BH)&\ast&\ast\\\
\hat{B}_{w}^{T}&-\lambda I&\ast\\\
C(XE^{T}+E^{\perp}W)+DH&\hat{D}_{w}&-I\end{bmatrix}\prec 0.$

##  Appendix B The Lyapunov-based approach to design dense estimators

The dense $\mathcal{H}_{\infty}$ estimator proposed by [8] is obtained as
follows:

$\displaystyle\boxed{L=(N(XE+E^{T\perp}Y)^{-1})^{T}},$ (26)

for which the following convex SDP is solved accordingly:

$\displaystyle\underset{\epsilon,\lambda,N,X,Y}{\min}$
$\displaystyle\leavevmode\nobreak\ \lambda$ $\displaystyle\mathrm{s.t.}$
$\displaystyle\leavevmode\nobreak\ X\succ
0,\epsilon>0,\lambda>0,\begin{bmatrix}\Omega&\ast&\ast\\\
XE+E^{T\perp}Y&-\epsilon I&\ast\\\ B_{w}^{T}(XE+E^{T\perp}Y)&0&-\lambda
I\end{bmatrix}\prec 0.$

where
$\Omega=\mathbf{He}(A^{T}(XE+E^{T\perp}Y)+C_{y}^{T}N)+\epsilon\alpha^{2}I$.