Decentralized control methodology for multi-machine/multi-converter power
systems

Semester Thesis
Spring Semester 2018

_Author:_

Aidar Zhetessov
ID: 15-948-268

_Supervisor:_
Catalin Arghir
_Professor:_
Prof. Dr. F. Dörfler

June 4, 2018

###### Abstract

In this project we evaluate a framework for synchronization of mixed machine-
converter power grids. Synchronous machines are assumed to be actuated by
mechanical torque injections, while the converters by DC-side current
injections. As this approach is based on model-matching, the converter’s
modulation angle is driven by the DC-side voltage measurement, while its
modulation amplitude is assigned analogously to the electrical machine’s
excitation current.

In this way we provide extensions to the swing-equations model, retaining
physical interpretation, and design controllers that achieve various
objectives: frequency synchronization while stabilizing an angle configuration
and a bus voltage magnitude prescribed by an optimal power flow (OPF) set-
point. We further discuss decentralization issues related to clock drifts,
loopy graphs, model reduction, energy function selection and characterizations
of operating points. Finally, a numerical evaluation is based on experiments
from three- and two- bus systems.

###### Contents

  1. 1 Introduction
    1. 1.1 Background and motivation
    2. 1.2 Problem formulation
    3. 1.3 Matching-control review
    4. 1.4 Modeling and notation
    5. 1.5 Outline
  2. 2 Operating point characterization
    1. 2.1 Reduced model
    2. 2.2 Full model
  3. 3 Topology independence
    1. 3.1 Cycled topologies
    2. 3.2 Angle consensus in a cycled graph
    3. 3.3 Shifted angle consensus in a cycled graph
    4. 3.4 Remarks on full model, robustness and decentralization
  4. 4 Decentralization
    1. 4.1 Full model approaches
      1. 4.1.1 Feedback linearization
      2. 4.1.2 Incremental energy function
    2. 4.2 Reduced model decentralized control
    3. 4.3 Simplified full model decentralized control
    4. 4.4 Full model decentralized control
  5. 5 Conclusion
    1. 5.1 Summary
    2. 5.2 Outlook

###### List of Figures

  1. 1.1 Circuit diagram of a 3-phase inverter [1]
  2. 1.2 Matrix representation of reduced network model
  3. 1.3 Matrix representation of full network model
  4. 2.1 Schematics of decentralized control setup [2]
  5. 3.1 Reduced power system model with cycled graph of generators
  6. (a) Cycled graph model
  7. (b) Single synchronous generator model
  8. 3.2 Example of two stable load flow solutions [3]
  9. (a) Nominal power flow
  10. (b) Local power flow
  11. 3.3 Connection graph $\mathbf{E}$ and communication graph $\mathbf{B}$
  12. 3.4 Two-node gradient-creating normalized energy function $S_{B}$
  13. 3.5 Three-node gradient-creating normalized energy functions
  14. (a) 3-node normalized $S_{E}$
  15. (b) 3-node normalized $S_{B}$
  16. 3.6 Energy-function-candidate plots for various angle differences
  17. 4.1 Angular frequency and angle transients with exact controller
  18. (a) Angular frequency transient
  19. (b) Angles transient
  20. 4.2 Node voltage and line current transients with exact controller
  21. (a) $v_{\mathsf{dq}}-v_{\mathsf{dq}}^{*}$
  22. (b) $i_{\mathsf{t,dq}}-i_{\mathsf{t,dq}}^{*}$
  23. 4.3 Angular frequency and angle transients with approximate controller
  24. (a) Angular frequency transient
  25. (b) Angles transient
  26. 4.4 Node voltage and line current transients with approximate controller
  27. (a) $v_{\mathsf{dq}}-v_{\mathsf{dq}}^{*}$
  28. (b) $i_{\mathsf{t,dq}}-i_{\mathsf{t,dq}}^{*}$
  29. 4.5 Control action measures for generator 1
  30. 4.6 Control action measures for generator 2
  31. 4.7 Control action measures for generator 3

##  Chapter 1 Introduction

###  1.1 Background and motivation

Modern power system is undergoing significant changes with more renewable
energy sources (RES) penetrating conventional grid through DC/AC converters,
or simply inverters. In the light of this trend preserving voltage and
frequency stability in the grid is becoming a more and more difficult task due
to intermittent nature and inherent lack of rotational inertia in RES. On the
other hand, interfacing power electronics provide much more flexibility in
control that could be used to make up for mentioned RES-inherent deficiencies.

A broad spectrum of inverter control approaches has two main categories: grid-
forming and grid-following [4]. In the grid-following control schemes an
inverter, assuming stiff grid, simply attaches to the network phase through
PLL to inject power. Since a stable low-inertia grid cannot be formed only by
grid-following inverters [1], grid-forming control approaches, which provide
network stability and ancillary services, should be investigated, especially
for islanded micro-grids.

There are several grid-forming control approaches including, but not limited
to droop, virtual synchronous machine (VSM), virtual oscillator control (VOC)
and IoT/ICT based approaches. A brief, yet comprehensive overview of those
with further references is provided in [1]-[5]. Another promising approach
called matching control was recently presented in [1]. It utilizes the
similarities between synchronous generator (SG) and inverter dynamics under
certain modulation of the latter. This allows not merely emulating SG
behavior, but exactly match it, so that there is no observable difference from
grid side. This approach has several advantages such as elimination of costly
AC side measurements, associated redundant computations and consequent
actuation delays, incorporation of DC side to the dynamics, which ameliorates
the demanding time-scale separation condition, ease of DC side measurement and
remaining degree of freedom for ancillary services.

Within this local model-matching an overarching secondary control needs to be
designed to achieve network level objectives like frequency synchronization,
voltage control and optimal power flow (OPF) imposition. Energy-based (or
passivity-based) framework, inspired by [6]-[7] and similar works, was pursued
in previous works of IfA lab of ETH Zürich, and this thesis continues progress
in that direction.

###  1.2 Problem formulation

With this being said, the problem at hand could be formulated as follows:
first, given a prescribed OPF set-points for a network, within the local
model-matching, what are the references for every local generator/inverter? In
other words, given network-level objectives, what should every separate
generator/inverter do locally? And second, supposing that there is a way of
explicit network-level objective characterization through local references,
how to stabilize the desired set-points in a topology-independent [8] and
decentralized [9] manner? Research in these directions is the main focus of
this thesis.

###  1.3 Matching-control review

For in-depth discussion on matching control the reader is referred to [1].
Here only a brief review will be added to remind the main aspects related the
present work.

As it was mentioned, matching control builds on similarities between inverter
and SG dynamics. Specifically it was observed that under certain inverter
modulation strategy, in which the modulation reference depends on virtual
angle, being the integral of DC-side voltage, inverter and SG dynamics can be
exactly model-matched, thus eliminating the necessity of emulation, virtual
modeling etc. Mentioned model-match becomes more evident if one writes the
dynamic equations of both systems as well as associated figure:

Figure 1.1: Circuit diagram of a 3-phase inverter [1]

$\begin{split}C_{dc}\dot{v}_{dc}&=-G_{dc}v_{dc}+i_{dc}-\frac{1}{2}i_{\alpha\beta}^{\top}m_{\alpha\beta}\\\
L\dot{i}_{\alpha\beta}&=-Ri_{\alpha\beta}-v_{\alpha\beta}+\frac{1}{2}v_{dc}\\\
C\dot{v}_{\alpha\beta}&=-Gv_{\alpha\beta}-i_{l}+i_{\alpha\beta}\\\
m_{\alpha\beta}&=\mu\begin{bmatrix}-sin(\theta_{v})\\\
cos(\theta_{v})\end{bmatrix}\ \omega_{v}=\eta
v_{dc}\end{split}\Leftrightarrow\begin{split}\dot{\theta}&=\omega\\\
M\dot{\omega}&=-D\omega+\tau_{m}+i_{\alpha\beta}^{\top}L_{m}i_{f}\begin{bmatrix}-sin(\theta_{v})\\\
cos(\theta_{v})\end{bmatrix}\\\
L_{s}\dot{i}_{\alpha\beta}&=-Ri_{\alpha\beta}-v_{\alpha\beta}-\omega_{v}L_{m}i_{f}\begin{bmatrix}-sin(\theta_{v})\\\
cos(\theta_{v})\end{bmatrix}\\\
C\dot{v}_{\alpha\beta}&=-Gv_{\alpha\beta}-i_{l}+i_{\alpha\beta}\end{split}$
(1.1)

As it can be seen, inverter dynamics under mentioned modulation strategy match
SG equations. This is the essence of local matching control. An advantage of
proposed local control is in the model-equivalence of SGs and inverters, that
somewhat simplifies the problem of overarching secondary control. The latter
one, in turn, still needs to be designed for network-level objectives, such as
load- and drift-robustness, voltage, frequency control etc.

###  1.4 Modeling and notation

Thanks to local matching, we have a unified framework for local modeling of
inverters and SGs. Now, before tackling the questions stated in Problem
formulation section, appropriate network model has to be defined. A good model
should be simple enough, so that the problem is still tractable, and at the
same time accurate enough, so that the main system aspects are captured.

Throughout the work two main models of the network were used - a reduced and
full one. The adoption of two distinct models is the consequence of a good
model definition above. Indeed, with full model all system aspects are
captured, but analyzing it directly could become a difficult task. Therefore,
as an intermediate step, reduced model - kind of a toy model, which preserves
network, but simplifies node dynamics, could be used to gain some insight on
the system. Matrix representation of a reduced model is depicted on Fig.1.2.

Figure 1.2: Matrix representation of reduced network model

Here we assume that inverter output voltage vector $\xi$ consists of ideal
sinusoidal voltages corresponding to inverters at every node and therefore
there is no need in output $LC$ filter, so that resulting currents are
directly injected into the network, represented as a Laplacian matrix. In SG
notation reduced model neglects stator dynamics, assuming that generated $\xi$
is directly applied to the network. Proposed model, despite being a ”toy”, is
quite general, since, in addition to preserving line dynamics, there are no
constraints on network topology or, in other words, on graph incidence matrix
$\mathbf{E}$.

Matrix representation of a full network model is shown on Fig.1.3. Here output
filter dynamics and local load $G$ are added. In SG notation those correspond
to stator and load dynamics.

Figure 1.3: Matrix representation of full network model

As one can see from Fig.1.3, there are several dynamic components at AC-side
of the inverter/SG. A compact and comprehensive summary of whatever comes at
AC-side at steady state for nominal frequency is given by so called $\Pi$-map:

$\Pi=\begin{bmatrix}\big{(}\mathbf{Z}_{s}+(\mathbf{Y}_{c}+\mathcal{L}_{\mathsf{t}})^{-1}\big{)}^{-1}\\\
(\mathbf{Y}_{c}+\mathcal{L}_{\mathsf{t}})^{-1}\big{(}\mathbf{Z}_{s}+(\mathbf{Y}_{c}+\mathcal{L}_{\mathsf{t}})^{-1}\big{)}^{-1}\\\
\mathbf{Z}_{\mathsf{t}}\mathbf{E}^{\top}(\mathbf{Y}_{c}+\mathcal{L}_{\mathsf{t}})^{-1}\big{(}\mathbf{Z}_{s}+(\mathbf{Y}_{c}+\mathcal{L}_{\mathsf{t}})^{-1}\big{)}^{-1}\end{bmatrix}=\begin{bmatrix}\Pi_{1}\\\
\Pi_{2}\\\ \Pi_{3}\end{bmatrix}$ (1.2)

Provided a fundamental component of SG/inverter output voltage vector $\xi$ at
steady state, one can derive steady state fundamentals of stator currents,
node voltages and transmission line currents throughout the network using this
$\Pi$-map, comprising $\Pi_{1},\Pi_{2},\Pi_{3}$ respectively
($[\hat{i}_{s}^{\top},\hat{v}^{\top},\hat{i}_{\mathsf{t}}^{\top}]^{\top}=\Pi\xi$).

As it was already mentioned, this work, building on previous findings at IfA
lab, continues to develop energy-based framework for decentralized secondary
control in power systems. The notation was inherited from these previous
works, so the last paragraph was taken as a reference to review the
established notation that will be used in this work as well. For full
derivation and description of notation the reader is referred to [10].

Consider the $n$-dimensional angle variable $\theta\in\mathbb{T}^{n}$. Let
$\Phi\in^{2n}$ be its Euclidean-space embedding such that
$\Phi_{i}=\begin{bmatrix}\cos\theta_{i}\\\
\sin\theta_{i}\end{bmatrix}\in^{2}$. Considering $I_{n}$ as the identity in
$n$-dimensional Euclidean space, we denote $j=\begin{bmatrix}0&-1\\\
1&0\end{bmatrix}$. Further define $\boldsymbol{j}=I_{n}\otimes j\in^{2n\times
2n}$,
$\mathbf{R}_{\theta_{i}}=\begin{bmatrix}\Phi_{i}&j\Phi_{i}\end{bmatrix}\in^{2\times
2}$ and
$\mathbf{R}_{\theta}=\text{blockdiag}\left\\{\mathbf{R}_{\theta_{i}}\right\\}\in^{2n\times
2n}$. We will aslo make use of $\mathbf{e}_{1}=\begin{bmatrix}1\\\
0\end{bmatrix}$ and $\mathbf{e}_{2}=j\mathbf{e}_{1}$. Transient error
coordinates are introduced for all real-valued state-space variables, for
example the following equation holds for $\forall t$: $v=\tilde{v}+\hat{v}$,
where $\tilde{v}$ is a transient part and $\hat{v}$ is a steady state function
of $\xi$ through $\Pi_{2}$. $\alpha\beta$-to-$\mathsf{dq}$ transformation is
done using nominal frequency $\omega_{0}$ throughout the work.

###  1.5 Outline

The rest of report is organized as follows: Chapter 2 exploits the
observation, found in [2], to tackle the first question in problem
formulation. Namely, it provides a bijective way of translation between
network-level OPF-prescribed set-point and local controller references used in
model-matching. Chapter 3 builds on obtained knowledge of local references and
introduces the notion of communication graph to deal with circulating power
flow problem within the topology-independence discussion. Inherent robustness
issues of communication-based approach motivate the study of decentralization,
a detailed discussion and substantial results on which are provided in Chapter
4. Finally, the summary of main findings and an outlook are given the
concluding Chapter 5.

##  Chapter 2 Operating point characterization

One of the main duties of central/emergency control level, provided grid state
estimate, is to calculate the network-level OPF set-points. The OPF algorithm
gives a consistent transmission line power flow references $P^{*}$ and $Q^{*}$
that preserve network frequency and bus voltage amplitudes within satisfactory
limits [4]. These power flow references are defined along transmission lines.
The controllers, however, actuate on nodes of the network either in form of
controllable DC current or mechanical torque injections. Therefore, to ensure
these controllers have their references, a translation way from network-level
OPF-prescribed set-points to local SG/inverter references should be
established, especially if one considers further efforts in topology-
independence and decentralization directions.

Both calculation and translation of OPF are difficult tasks, which even have
their dedicated research directions. However the fact that OPF imposition
applies to steady state significantly simplifies the translation task, because
now a simpler controllers/models, that manage to impose an OPF at steady state
and provide equations for bus voltages, can be studied. Taking mentioned
control output bus voltages, one could reverse back into SG/inverter model at
the nodes and derive local references ($I_{r}^{*},\theta^{*}$ for SG. For
inverter these parameters enter $m_{\alpha\beta}$ implicitly (1.1))

One of such useful controllers was derived in [2]. In particular, it managed
to almost globally stabilize the OPF-prescribed set-points on a network with
shunt resistance (load), capacitance and controllable current source. Node
schematic and proposed control law are provided below:

Figure 2.1: Schematics of decentralized control setup [2]

$\begin{split}&\frac{d}{dt}v_{k}=u_{k}(v_{k},i_{o,k})=\omega_{0}jv_{k}+\eta\big{(}K_{k}v_{k}-R(\rho)i_{o,k}\big{)}+\alpha\Phi_{k}(v_{k})v_{k}\
\ \text{where}\\\ &\rho=tan^{-1}(\frac{\omega_{0}l_{jk}}{r_{jk}})\ \
K_{k}=\frac{1}{v_{k}^{*2}}R(\rho)\begin{bmatrix}p_{k}^{*}&q_{k}^{*}\\\
-q_{k}^{*}&p_{k}^{*}\end{bmatrix}\ \
\Phi_{k}(v_{k})=\frac{v_{k}^{*}-||v_{k}||}{v_{k}^{*}}\end{split}$ (2.1)

What else was observed is that under $\frac{R}{L}=const.$ assumption for
transmission lines, system dynamics can be represented current-free or as a
function of node voltages only:

$\begin{split}&\frac{d}{dt}v=\omega_{0}Jv+\eta\big{(}\mathcal{K}-\mathcal{L}\big{)}v+\alpha\Phi(v)v\
\ \text{where}\\\ &\mathcal{K}=diag\big{(}\\{K_{k}\\}_{k=1}^{N}\big{)}\ \
\Phi(v)=diag\big{(}\\{\Phi_{k}(v_{k})\\}_{k=1}^{N}\big{)}\end{split}$ (2.2)

If we slightly rewrite the equation (2.2) we end up with the first vector
field in equation (2.3), which is an explicit operating point characterization
in terms of OPF. Obviously, vector field $f_{1}$ explicitly depends on power
flow set-points $P^{*}$ and $Q^{*}$.

$\begin{split}&\text{Explicit OPF (reference):}\
\dot{v}_{\mathsf{dq}}=f_{1}(P^{*},Q^{*},||v_{\mathsf{dq}}^{*}||,v_{\mathsf{dq}})\\\
&\text{Implicit OPF (this work):}\
\dot{v}_{\mathsf{dq}}=f_{2}(I_{r}^{*},\theta_{\mathsf{dq}}^{*},v_{\mathsf{dq}},...)\end{split}$
(2.3)

Now, in this thesis we have an implicit operating point characterization in
terms of OPF, since our vector field $f_{2}$ depends on local rotor current
and angle references as well as some other parameters, which implicitly carry
information on network-level OPF. Since the whole models in [2] (Fig.2.1) and
this thesis (Figures 1.2-1.3) are very different, corresponding vector fields
$f_{1}$ and $f_{2}$ are also very different, but what’s in common is
transmission line network, modeled in both works as an extended Laplacian
matrix $\mathcal{L}$, or even a $\Pi$-configuration of transmission line in
full model case. Indeed, reference [2] has controllable current sources
attached to every node, which corresponds to the whole part on the left of bus
voltage $v$ in Fig.1.3. The remaining shunt resistance, capacitance and
network, is exactly the same as part of the full model on the right of bus
voltage $v$ and is nothing more that $\Pi$-configuration of transmission
lines. So, using this observation, one could argue that node (or bus) voltage
set-points $v_{\alpha\beta}^{*}$ (or $v_{\mathsf{dq}}^{*}$) that impose OPF on
a network in case of the first vector field at steady state, should also do
the same in case of second vector field at steady state, since these set-
points uniquely determine power flow in the network. Proposed argument can be
used in both reduced and full models.

###  2.1 Reduced model

For a reduced model there is a direct relation between $\xi^{*}$ at steady
state and $I_{r}^{*},\theta^{*}$. Also recall the assumption of $\xi$ being an
ideal sinusoid for a reduced model. So for local-reference-derivation, one
could solve for node voltage set-points $v_{\alpha\beta}^{*}$ (or
$v_{\mathsf{dq}}^{*}$) by solving steady state equation of vector field
$f_{1}$ or, alternatively, by analyzing $ker(\mathcal{K}-\mathcal{L})$,
provided $||v_{\mathsf{dq}}^{*}||$ holds (last part of controllers (2.1) and
(2.2). Now, recall that $\xi_{\mathsf{dq}}^{*}=v_{\mathsf{dq}}^{*}$ for
reduced model, moreover,
$\xi_{\mathsf{dq}}^{*}=\mathbf{R}_{\theta_{\mathsf{dq}}^{*}}(L_{m}\otimes
e_{2})I_{r}^{*}\omega_{0}\boldsymbol{1}$ [10]. From these one can derive
$I_{r}^{*},\theta_{\mathsf{dq}}^{*}$. The steps are summarized below:

$\begin{split}\left.\begin{aligned} &\text{Explicit OPF (reference):}\
\dot{v}_{\mathsf{dq}}=f_{1}(P^{*},Q^{*},||v_{\mathsf{dq}}^{*}||,v_{\mathsf{dq}})\\\
&\text{Implicit OPF (this work):}\
\dot{v}_{\mathsf{dq}}=f_{2}(I_{r}^{*},\theta_{\mathsf{dq}}^{*},v_{\mathsf{dq}},...)\end{aligned}\right\\}&\text{at
steady state}\ \ \dot{v}_{\mathsf{dq}}=0\rightarrow\\\ &\rightarrow
v_{\mathsf{dq}}^{*}\in ker(\mathcal{K}-\mathcal{L})\\\
&\rightarrow\xi_{\mathsf{dq}}^{*}=v_{\mathsf{dq}}^{*}\rightarrow
I_{r}^{*},\theta_{\mathsf{dq}}^{*}\end{split}$ (2.4)

If we plug obtained local references $I_{r}^{*},\theta_{\mathsf{dq}}^{*}$ into
vector field $f_{2}$, we end up imposing very same OPF on our reduced model,
provided that same network topology and parameters are used both in $f_{1}$
and $f_{2}$, including $\frac{R}{L}=const.$ assumption for transmission lines.
Thus, there is a bijective way of translation between network-level OPF-
prescribed set-points and local references for SG/inverter controllers of a
reduced model. The argument is bidirectional, since provided consistent local
references, it is possible (and even much easier) to derive network power
flows.

###  2.2 Full model

For full model, the procedure has an additional step, as steady state node
voltage set-points $v_{\mathsf{dq}}^{*}$ and SG/inverter fundamental of output
voltage vector $\xi_{\mathsf{dq}}^{*}$ are related by $\Pi_{2}$-map
($\Pi_{2}\xi_{\mathsf{dq}}^{*}=v_{\mathsf{dq}}^{*}$). Note that for full model
assumption of noise-free $\xi_{\mathsf{dq}}^{*}$ does not hold, so $\Pi_{2}$
relates only fundamentals of $\xi_{\mathsf{dq}}^{*}$ and
$v_{\mathsf{dq}}^{*}$. However, this is the actual reason of introducing
$\Pi_{2}$, because in inverter notation $L_{s},C$ (SG stator parameters) can
be seen as output filter. Interestingly, it does not even need to be as large,
since we assume simple 3-phase 2-level 6-switch inverter average model in
local matching control. Provided high switching frequency for low-voltage (LV)
range, high corner frequency of output $LC$ filter would suffice. Additional
step for full model is included in (2.5):

$\begin{split}\left.\begin{aligned} &\text{Explicit OPF (reference):}\
\dot{v}_{\mathsf{dq}}=f_{1}(P^{*},Q^{*},||v_{\mathsf{dq}}^{*}||,v_{\mathsf{dq}})\\\
&\text{Implicit OPF (this work):}\
\dot{v}_{\mathsf{dq}}=f_{2}(I_{r}^{*},\theta_{\mathsf{dq}}^{*},v_{\mathsf{dq}},...)\end{aligned}\right\\}&\text{at
steady state}\ \ \dot{v}_{\mathsf{dq}}=0\rightarrow\\\ &\rightarrow
v_{\mathsf{dq}}^{*}\in ker(\mathcal{K}-\mathcal{L})\\\
&\rightarrow\xi_{\mathsf{dq}}^{*}=\Pi_{2}^{-1}v_{\mathsf{dq}}^{*}\rightarrow
I_{r}^{*},\theta_{\mathsf{dq}}^{*}\end{split}$ (2.5)

If we plug obtained local references $I_{r}^{*},\theta_{\mathsf{dq}}^{*}$ into
vector field $f_{2}$, we end up imposing very same OPF on our full model as
well, provided that same network topology and parameters are used both in
$f_{1}$ and $f_{2}$, including $\frac{R}{L}=const.$ assumption for
transmission lines (same conditions). Thus, there is a bijective way of
translation between network-level OPF-prescribed set-points and local
references for SG/inverter controllers of a full model as well. This argument
is also bidirectional, since provided consistent local references, it is
possible (and even much easier) to derive network power flows. Mentioned
bijective relation can be summarized in the following equation:

$P^{*},Q^{*},||v_{\mathsf{dq}}^{*}||\leftrightarrow
I_{r}^{*},\theta_{\mathsf{dq}}^{*}$ (2.6)

Existence of explicit translation (2.6) for both models of this thesis allows
studying decentralization and (to lesser extent) topology-independence,
knowing that local references can impose global OPF on a network.

Note on angles: after estimating desired bus voltage angle at steady state for
vector field $f_{1}$, [2] proposes selecting one node to be reference in terms
of angle. Such node’s desired angle and initial condition are both zero. All
others converge to the difference with respect to reference node angle. This
is because ideally we do not care about exact phase of bus voltages. Power
flow is determined by bus-voltage-differences and not by exact phase of a
particular bus. Reduced and full models of this thesis are perfectly fine with
mentioned angle treatment, since in our case we also care only about
rotor/virtual angle differences.

##  Chapter 3 Topology independence

This chapter is dedicated to the consensus-oriented control of a large-cycled
multi-machine power systems. Reduced model is to be used predominantly, since
additional complication of stator in full model does not affect the control
strategy. An example reduced model of a small, cycled power system is depicted
in Fig.3.1.

(a) Cycled graph model

(b) Single synchronous generator model

Figure 3.1: Reduced power system model with cycled graph of generators

As it can be seen from the figure, cycled graph comprises three nodes and
edges connected in a triangle, which closes the cycle (loop). Graph edges
model three-phase transmission lines in $\alpha\beta$-frame. A synchronous
generator (SG), modeled and depicted in Fig.3.1(b), is connected to every
node. $\xi_{e}$ and $i_{r}^{*}$ denote line-connected SG EMF and reference
rotor current.

When it comes to modeling networks, connecting ideal and fully controllable
voltage/current sources to graph edges is practiced in some literature
[11]-[2]. Similarly, when modeling SG, network is usually modeled as a simple
infinite-bus [12]-[13]. Serving as a bridge, proposed reduced and full models
have this additional advantage. On one hand, they do not oversimplify SG
model, while on the other, preserve network dynamics.

First, the chapter begins with a motivation, introducing the problem of cycled
topologies. Then it proposes the reduced model communication-based control
implementation of identical/shifted rotor angle consensus with special focus
on cycled graphs. Underlying remarks on robustness with respect to
communication and clock drift are included in the end.

###  3.1 Cycled topologies

Several literature sources [8],[3]-[14] report the possibility as well as some
underlying reasons of circulating power flow emergence in cycled topologies
with five or more generators in the loop. Some conclusions are even backed up
by real observations in power systems. It turns out that in large cycled
topologies several local minima of network energy can arise. They correspond
to circulating power flow in either clockwise or counter-clockwise direction.
Interestingly, this power circulation is stable even in dissipative
(resistive) transmission line networks under certain conditions. An example
graph (or network), claimed to have real counterpart with observed circulation
phenomenon, is shown in Fig.3.2:

(a) Nominal power flow

(b) Local power flow

Figure 3.2: Example of two stable load flow solutions [3]

Network topology of Fig.3.2 can have three stable power flow configurations,
corresponding to the sum of rotor angles around the loop being $0,2\pi$ or
$-2\pi$. More generally, the maximum number of possible local minima (stable
power flows) $N$ depends on the number of generators around the loop $n$ as:
$N=2\ \text{Int}[n/4]+1$. Respective rotor angle configurations could sum up
around the loop as follows: $\sum^{loop}\delta_{ij}=\pm 2\pi k\ \text{where}\
k\in\\{-a...0...a\\}\ |\ 2a+1=N$ [3]. From equations above one can see the
reason of at least 5 generators. Indeed, with more nodes it is possible to
have angle loop-sum of an integer multiples of $2\pi$, yet relative angle
difference between neighbors is less than $\pi/2$, thus allowing local
stability without slipping off.

Obviously, circulating power flows are undesired due to excessive power losses
in transmission lines. To deal with this issue, one way to proceed, provided
communication between graph neighbors, would be to introduce the notion of
communication graph. Note that the method to be proposed relies on real-time
angle communication between neighbors. A decentralized and more robust
approach could be possible, however, it requires more thorough
study/understanding of circulating power emergence, which was not main
objective of this work.

###  3.2 Angle consensus in a cycled graph

Identical angle consensus was already reported in [10], which proposes energy-
based control for network flow stabilization in high-order multi-machine power
system. In short, key idea in that work, referred to as an incremental energy
function approach, was to come up with (state- and communicated angle-
dependent) energy function gradient, acting on a machine torque balance
equation as a closed-loop input, which, based on Lyapunov-LaSalle argument,
stabilizes the whole system at identical angle configuration. Inspired by
proposed approach, this section addresses the question of identical angle
consensus in cycled graphs based on reduced model.

Arguments of [10], in general, were independent of graph topology and,
therefore, are applicable to cycled graphs too. Nevertheless, a separate study
of consensus in cycled graphs through a reduced model is still attractive due
to aforementioned issue of cycled topologies.

Now, consider the following scenario: triangular network (Fig.3.1(a)), reduced
model (Fig.1.2) with objective of rotor-angle consensus (all rotor angles are
equal, so no power flow in symmetric network). Reduced model dynamic equations
in $\alpha\beta$-frame look like:

$\begin{split}\dot{\theta}&=\omega\\\
M\dot{\omega}&=-D\omega+\tau_{m}-I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta}^{\top}\mathbf{E}i_{\mathsf{t}}\\\
L_{\mathsf{t}}\dot{i}_{\mathsf{t}}&=-R_{\mathsf{t}}i_{\mathsf{t}}+\mathbf{E}^{\top}\mathbf{R}_{\theta}(L_{m}\otimes
e_{2})I_{r}^{*}\omega\end{split}$ (3.1)

Explanations of various terms in (3.1) are summarized in Table 3.1:

Table 3.1: Summary of model-related variables Variable |  Dimension |  Explanation
---|---|---
$\theta$ |  $\mathbb{T}^{n}$ |  rotor angles of $n$ generators
$\omega$ |  n |  generator angular frequencies
$M,D$ |  n×n |  gen. rotational inertia & damping diag. matrices
$\tau_{m}$ |  n |  mechanical torque inputs
$I_{r}^{*}$ |  n |  rotor current set-points
$L_{m}$ |  n |  mutual coupling inductances
$\mathbf{E}$ |  2n×2m |  expanded incidence matrix $E\otimes I_{2}$
$L_{\mathsf{t}},R_{\mathsf{t}}$ |  2m×2m |  expanded line parameter matrices $L_{t},R_{t}\otimes I_{2}$
$i_{\mathsf{t}}$ |  2m |  transmission line currents in $\alpha\beta$-frame
$I_{r}^{*}(L_{m}\otimes e_{2}^{\top})\mathbf{R}_{\theta}^{\top}\mathbf{E}i_{\mathsf{t}}$ |  n |  electrical torques
$\mathbf{R}_{\theta}(L_{m}\otimes e_{2})I_{r}^{*}\omega$ |  2n |  electro-motive forces in $\alpha\beta$-frame

The same model in $\mathsf{dq}$-frame looks like:

$\begin{split}\dot{\theta}_{\mathsf{dq}}&=\tilde{\omega}\\\
M\dot{\tilde{\omega}}&=-D\tilde{\omega}+\tilde{\tau}_{m}-I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}\mathbf{E}i_{\mathsf{t,dq}}\\\
L_{\mathsf{t}}\dot{i}_{\mathsf{t,dq}}&=-\mathbf{Z}_{\mathsf{t}}i_{\mathsf{t,dq}}+\mathbf{E}^{\top}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}(\tilde{\omega}+\omega_{0}\boldsymbol{1})\end{split}$ (3.2)

Here $\omega=\tilde{\omega}+\omega_{0}\boldsymbol{1}$,
$\tau_{m}=\tilde{\tau}_{m}+D\omega_{0}$, $\tilde{\tau}_{m}$ is a mechanical
torque input in $\mathsf{dq}$-frame, while $D\omega_{0}$ is a steady-state
torque. Rising the question of appropriate $\tilde{\tau}_{m}$, in general, it
is desired to design this closed-loop input such that possible local minima of
energy function in cycled graph case are eliminated. Graph topology affects
rotor dynamics through electrical torque vector, therefore incremental energy
function approach, used in [10], cannot be directly applied. Indeed, mentioned
approach uses line current measurements just to modify electrical torque, or,
in other words, to arrive at the gradient of aforementioned problematic energy
function. Modifying electrical torque in some other way could potentially
solve the problem, but, having the design freedom of our closed-loop input and
angle communication, it seems more wise to eliminate the effect of cycled
graph totally by canceling electrical torque and adding the desired gradient
on top. The gradient of interest should not have local minima problems in its
associated energy function, just as if the graph did not have any cycles at
all.

Bearing in mind the considerations above, one proposes to define a separate
communication graph $\mathbf{B}$ that connects all the nodes of $\mathbf{E}$
graph and has all the desired properties, such as no cycles, appropriate edge
weights, topology etc. Angles are to be communicated based on $\mathbf{B}$
rather than $\mathbf{E}$ graph neighborhood, gradient term of the torque is to
be created as if real connection is described by $\mathbf{B}$ and not by
$\mathbf{E}$ graph. Example modification of the graph in Fig.3.1(a) is
depicted in Fig.3.3.

Figure 3.3: Connection graph $\mathbf{E}$ and communication graph $\mathbf{B}$

So, with the notion of communication graph we kind of ”deceive” our generators
and make them think that they are in an un-cycled graph with incidence matrix
$\mathbf{B}$ rather than in cycled graph $\mathbf{E}$. Summarizing the closed-
loop input properties stated before, one suggests the following candidate:
$\tilde{\tau}_{m}=\tau_{e}-\nabla_{\theta_{\mathsf{dq}}}S_{B}$. Closed-loop
dynamics with such feedback become:

$\begin{split}\dot{\theta}_{\mathsf{dq}}&=\tilde{\omega}\\\
M\dot{\tilde{\omega}}&=-D\tilde{\omega}-\nabla_{\theta_{\mathsf{dq}}}S_{B}\\\
L_{\mathsf{t}}\dot{i}_{\mathsf{t,dq}}&=-\mathbf{Z}_{\mathsf{t}}i_{\mathsf{t,dq}}+\mathbf{E}^{\top}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}(\tilde{\omega}+\omega_{0}\boldsymbol{1})\end{split}$ (3.3)

Line and rotor dynamics are now decoupled (no dependence on $i_{\mathsf{t}}$
of the first two equations of (3.3)) with rotor driving line dynamics. This
approach resembles feedback linearization. Associated energy function can be
defined as
$S_{B}=\frac{1}{2}\Psi_{\mathsf{dq}}^{\top}BKB^{\top}\Psi_{\mathsf{dq}}$,
where $\Psi_{\mathsf{dq}}=\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{1})I_{r}^{*}\boldsymbol{1}$, $BKB^{\top}$ \- new un-cycled communication
graph with edge weight matrix $K$, that has only one stable diagonal on
$\mathbb{T}^{n}$ (property of un-cycled graphs). The gradient can be derived
as

$\nabla_{\theta_{\mathsf{dq}}}S_{B}=I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}BKB^{\top}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{1})I_{r}^{*}\boldsymbol{1}$ (3.4)

The convergence set of rotor dynamics is
$M_{1}=\\{(\theta_{\mathsf{dq}},\omega):\omega=\omega_{0}\boldsymbol{1},\nabla_{\theta_{\mathsf{dq}}}S_{B}=0\\}$.
This set is a main diagonal in $\mathbb{T}^{n}$ (state space of
$\theta_{\mathsf{dq}}$). Line currents, being a simple LTI system under
constant input entering in form of steady state angle
$\theta_{\mathsf{dq}}^{*}$, converge to
$i_{\mathsf{t,dq}}^{*}=\mathbf{Z}_{\mathsf{t}}^{-1}\mathbf{E}^{\top}\mathbf{R}_{\theta_{\mathsf{dq}}^{*}}(L_{m}\otimes
e_{2})I_{r}^{*}\omega_{0}\boldsymbol{1}$, so that entire system state
converges to
$M=M_{1}\times\\{i_{\mathsf{t,dq}}:i_{\mathsf{t,dq}}=i_{\mathsf{t,dq}}^{*}\\}$.
Moreover, in case $L_{m}I_{r}^{*}\in^{n\times n}$ has same elements along its
diagonal, $\mathbf{E}^{\top}\mathbf{R}_{\theta_{\mathsf{dq}}^{*}}(L_{m}\otimes
e_{2})I_{r}^{*}\omega_{0}\boldsymbol{1}=0\Rightarrow
i_{\mathsf{t,dq}}^{*}=0\Rightarrow$ no transmission line currents. For this
controller one requires line current measurements near every node and rotor
angles of neighbors in communication graph. For a simple network of two
generators and one transmission line normalized energy function
$S_{B}(\theta_{\mathsf{dq}})$ can be drawn as in Fig.3.4:

Figure 3.4: Two-node gradient-creating normalized energy function $S_{B}$

For three-node case, depicted in Fig.3.3, normalized energy function can still
be drawn, but in that case all three axes depict rotor angles, while
normalized energy value is conveyed through color. So called ”potential” plot
for both connection $\mathbf{E}$ and communication $\mathbf{B}$ graphs is
shown in Fig.3.5.

(a) 3-node normalized $S_{E}$

(b) 3-node normalized $S_{B}$

Figure 3.5: Three-node gradient-creating normalized energy functions

As it can be seen from Fig.3.5(a) no local non-diagonal minima are predicted
for 3-node cycled graph, as expected. Nevertheless, comparing Fig.3.5(a) and
Fig.3.5(b) interesting observation can be made. In both cases angle consensus
(main diagonal at which all angles are equal) is almost globally
asymptotically stable, but in $S_{B}$ case gradient is not as steep (which can
be seen from less-bright regions around diagonal) as for $S_{E}$ case,
potentially resulting in longer and larger transients.

###  3.3 Shifted angle consensus in a cycled graph

Having discussed identical angle consensus, this section focuses on a more
interesting case of shifted angle consensus, corresponding to non-zero
transmission line currents and, consequently, power transfer between graph
nodes. The general question of interest could be formulated as follows: given
the desired angle difference in $\mathsf{dq}$-frame, how to design the control
input such that mentioned angle configuration is achieved? Note that only
angle differences and not their exact values in $\mathsf{dq}$-frame are of
interest because transmission line currents and associated power transfer
depend on the difference of applied node EMFs. Exact angle values affect only
the phase shift of voltage and current vectors back in $\alpha\beta$-frame, as
was discussed at the end of Chapter 2.

To solve the issue of shifted consensus, desired angle configuration should be
somehow incorporated into the gradient-creating energy function. At least
these three energy-function-candidates can be considered due to non-Euclidean
nature of angle variable:

$\displaystyle\bar{S}_{B}$
$\displaystyle=S_{B}(\theta_{\mathsf{dq}}-\theta_{\mathsf{dq}}^{*})=\frac{1}{2}\Psi^{\top}(\theta_{\mathsf{dq}}-\theta_{\mathsf{dq}}^{*})BKB^{\top}\Psi(\theta_{\mathsf{dq}}-\theta_{\mathsf{dq}}^{*})$
(3.5a) $\displaystyle\tilde{S}_{B}$
$\displaystyle=\frac{1}{2}\big{(}\Psi^{\top}(\theta_{\mathsf{dq}})-\Psi^{\top}(\theta_{\mathsf{dq}}^{*})\big{)}BKB^{\top}\big{(}\Psi(\theta_{\mathsf{dq}})-\Psi(\theta_{\mathsf{dq}}^{*})\big{)}$
(3.5b) $\displaystyle\hat{S}_{B}$
$\displaystyle=S_{B}(\theta_{\mathsf{dq}})-S_{B}(\theta_{\mathsf{dq}}^{*})-\nabla_{\theta_{\mathsf{dq}}}^{\top}S_{B}(\theta_{\mathsf{dq}}^{*})(\theta_{\mathsf{dq}}-\theta_{\mathsf{dq}}^{*})$
(3.5c)

Here
$\Psi(\theta_{\mathsf{dq}}-\theta_{\mathsf{dq}}^{*})=\mathbf{R}_{(\theta_{\mathsf{dq}}-\theta_{\mathsf{dq}}^{*})}(L_{m}\otimes
e_{1})I_{r}^{*}\boldsymbol{1}$ and
$\Psi_{\mathsf{dq}}^{*}=\Psi(\theta_{\mathsf{dq}}^{*})=\mathbf{R}_{\theta_{\mathsf{dq}}^{*}}(L_{m}\otimes
e_{1})I_{r}^{*}\boldsymbol{1}$. Equation (3.5c) is referred to as a Bregman
energy function [15]. The Matlab plots of these energies for various angle
differences are depicted below:

Figure 3.6: Energy-function-candidate plots for various angle differences

As one can see, there are 12 realizations with every row corresponding to
respective energy function in equation (3.5), while every column being an
increasing angle difference: $(0,0),(0,10),(0,30),(0,90)$. The first energy
function nicely stabilizes the shifted diagonal in angle space (xy-plane in
this case). Second energy function also has some angles stabilized, but the
emergence of several valleys implies that invariant set is disjount, meaning
that depending on initial condition, one could end up in wrong valley (angle
configuration) and, consequently, wrong power flow. Third candidate, called
Bregman energy function turns out to be discontinuous in its $0\leftrightarrow
2\pi$ transition, meaning that it is not a good choice for our LaSalle-based
stability proof, which requires continuity [11]. All in all, provided
communication, the best way (among considered ones) to stabilize the shifted
angle consensus in a cycled graph of SG/inverters is to apply the proposed
closed-loop input with the gradient of the first energy function.

Regarding the interpretation of proposed candidates, first of all, remember
that we are talking about virtual energy functions - fictitious virtually
implemented functions that aim to stabilize prescribed configuration in angle
space at steady state. They may not always have real energy counterparts,
nevertheless, by looking at their structure and noticing the structural
similarities with some real energy functions, yet to be utilized in
decentralization chapter, a comprehensive insight can be gained on the meaning
of those candidates.

The first one, $\bar{S}_{B}$ (3.5a), is the ”energy” of angle difference. Its
idea is to minimize the difference between $\theta_{ij}$ and $\theta_{ij}^{*}$
in angle space. Since angle difference directly in angle space and not in the
Euclidean embedding is achievable by multiplication of respective rotation
matrices $\mathbf{R}_{\theta_{\mathsf{dq}}}$ and
$\mathbf{R}_{\theta_{\mathsf{dq}}^{*}}$, the resulting function looks like a
somehow bended manifold in Euclidean space - a result of multiplication by
desired angle rotations $\mathbf{R}_{\theta_{\mathsf{dq}}^{*}}$ and,
consequently, does not seem to have a real energy counterpart. From matrix
perspective, Laplacian of communication graph $\mathbf{B}$ in equation (3.5a)
provides necessary permutation and scaling of sine and cosine terms that,
after application of trigonometric rules, results in components having the
structure such as $sin(\theta_{ij}-\theta_{ij}^{*})$, which is almost globally
stable at desired angle consensus $\theta_{ij}=\theta_{ij}^{*}$.

The second candidate, $\tilde{S}_{B}$ (3.5b), is the ”energy” of flux linkage
difference. Its idea is to minimize the difference between flux linkages
$\Psi(\theta_{\mathsf{dq}})\ \text{and}\ \Psi(\theta_{\mathsf{dq}}^{*})$ in
Euclidean space. From structural perspective transmission line energy stored
in inductances $L_{\mathsf{t}}$ is very similar, which will become evident in
simplified model decentralized control section of next chapter, so
$\tilde{S}_{B}$ can be given the interpretation of transmission line energy.
Moreover, proposed interpretation explains the emergence of disjoint invariant
set in Fig.3.6. It is known that transmission line currents (and therefore
energy) depend on node voltage differences rather than their absolute value.
Undesired valley in angle space occurs due to constant step up in reference
voltages throughout the graph. This effect is kind of circulating current
counterpart in voltage notation, since it adds to node voltages, but not to
line currents. From matrix perspective, valley is explained by non-zero kernel
of incidence matrix. Indeed, the difference of angle embeddings in Euclidean
space could have a constant shift multiplied by vector of all ones
($\boldsymbol{1}$), which lies in kernel of incidence matrix and therefore
will not be observed neither in energy function nor in its gradient.

Regarding Bregman energy function, $\hat{S}_{B}$ (3.5c), its virtue is the
capability of shifting the gradient itself, rather than the references of its
entries in some form. This candidate was not studied extensively after
observing aforementioned discontinuity, thus exploring this option still
remains an interesting open topic for further work.

###  3.4 Remarks on full model, robustness and decentralization

Some brief remarks, observations and further considerations conclude Chapter
3:

$\bullet$ Regarding full model, control strategy does not change. The only
difference is in driven system, which has stator (or output filter in inverter
notation) in addition to line dynamics. Driven system still goes to its own
invariant set in angle-frequency state space, while enlarged driven LTI system
goes to its own steady state under effect of constant angle input, entering as
$\mathbf{R}_{\theta_{\mathsf{dq}}^{*}}$.

$\bullet$ As it was mentioned, gradient creation requires angle communication,
which is not a robust implementation in real time. As a decentralized
alternative, a fully decoupled energy function below could be considered:

$S_{B}=\frac{1}{2}(\Psi_{\mathsf{dq}}-\Psi_{\mathsf{dq}}^{*})^{\top}(L_{m}\otimes
I_{2})^{-1}(\Psi_{\mathsf{dq}}-\Psi_{\mathsf{dq}}^{*})$ (3.6)

where
$\Psi_{\mathsf{dq}}^{*}=\mathbf{R}_{\theta_{\mathsf{dq}}^{*}}(L_{m}\otimes
e_{1})I_{r}^{*}\boldsymbol{1}$. With this as an input, generators do not care
about the network and simply converge to their angle references, stabilizing
only a single point and not the whole diagonal in angle space, potentially
resulting in longer and larger transients.

$\bullet$ Although what was proposed in (3.6) is a decentralized realization,
which eliminates the necessity of communication altogether, it is not robust
at least with respect to controller clock drift. If one allows a tiny drift in
clocks of neighboring controllers, as time evolves the relative angle and,
consequently, the whole network, will drift away from its OPF set-points with
controllers being not even aware of the problem, because what they do after
all is compare self-measured (or -created in case of inverter) angle with
their references, both evolving using same local clock. Of course one way to
solve the communication problem would be to synchronize clocks and continuie
the work in that direction, but ideally it is not a good practice to fully
decouple subsystems. Some kind of feedback is still desired to preclude the
adverse effects of inevitable parameter fluctuations and model inaccuracies.

Robustness issues and coupling arguments above motivate studying ”coupled”
decentralization, which is the culminating chapter of the thesis.

##  Chapter 4 Decentralization

This chapter is dedicated to the decentralized control framework of small-
cycled or un-cycled multi-machine power systems. For now large-cycled graphs
were eliminated from considerations to ensure no sticking into local minimum
corresponding to circulating power flows during transients is possible. The
work is conducted using both reduced and full models interchangeably and the
chapter is structured in chronological order. First, several unsuccessful
trials of designing decentralized controller for full models are shown with
respective issues discussed. Second, after facing difficulties in tackling
full model directly, decentralized controller for reduced model, proposed
recently by authors of [10], is reviewed in search of inspiration. Third, led
by an idea of simplifying the full model to the level of reduced one, from
where a modified controller of previous section could be applied, a simplified
full model decentralized controller and associated simplification conditions
are discussed. It will be evident that there are several major problems
associated with simplification approach, which motivate studying full model
back again. Finally, inspired by last section’s controller as well as
topology-independence remarks, the chapter culminates with a ”coupled”
decentralized controller for full model. Numerical simulations, advantages and
challenges of the proposed controller are discussed in the end. Note that due
to excessive amount of formulas, several obvious steps are omitted for the
sake of readability and simplicity.

###  4.1 Full model approaches

####  4.1.1 Feedback linearization

This approach was already discussed in a concluding remark of topology-
independence chapter. Below are some useful reminder matrix equations as well
as full model dynamics under one-way-coupled (or feedback-linearized) input:
$[\hat{i}_{s,\mathsf{dq}}^{\top},\hat{v}_{\mathsf{dq}}^{\top},\hat{i}_{\mathsf{t,dq}}^{\top}]^{\top}=\boldsymbol{\Pi}\xi_{\mathsf{dq}}$
$[\hat{i}_{s}^{\top},\hat{v}^{\top},\hat{i}_{\mathsf{t}}^{\top}]^{\top}=\boldsymbol{\Pi}\xi$
$\boldsymbol{\Pi_{2}}=\big{[}\mathbf{Z}_{s}(\mathbf{Y}_{c}+\mathcal{L}_{\mathsf{t}})+I\big{]}^{-1}$
$v_{\mathsf{dq}}^{*}=\boldsymbol{\Pi_{2}}\xi^{*}$

$\begin{split}\dot{\theta}&=\omega\\\
M\dot{\omega}&=-D\omega+\tau_{m}-I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta}^{\top}i_{s}\\\
L_{s}\dot{i}_{s}&=-R_{s}i_{s}+\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}\omega-v\\\ C\dot{v}&=i_{s}-\mathbf{E}i_{\mathsf{t}}-Gv\\\
L_{\mathsf{t}}\dot{i}_{\mathsf{t}}&=-R_{\mathsf{t}}i_{\mathsf{t}}+\mathbf{E}^{\top}v\end{split}$
(4.1)

Same equation in $\mathsf{dq}$-frame and error coordinates for driving system
looks like:

$\begin{split}\dot{\theta}_{\mathsf{dq}}&=\tilde{\omega}\\\
M\dot{\tilde{\omega}}&=-D\tilde{\omega}+\tilde{\tau}_{m}-I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}i_{s,\mathsf{dq}}\\\
L_{s}\dot{i}_{s,\mathsf{dq}}&=-\mathbf{Z}_{s}i_{s,\mathsf{dq}}+\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}(\tilde{\omega}+\omega_{0}\boldsymbol{1})-v_{\mathsf{dq}}\\\
C\dot{v}_{\mathsf{dq}}&=i_{s,\mathsf{dq}}-\mathbf{E}i_{\mathsf{t,dq}}-\mathbf{Y}_{c}v_{\mathsf{dq}}\\\
L_{\mathsf{t}}\dot{i}_{\mathsf{t,dq}}&=-\mathbf{Z}_{\mathsf{t}}i_{\mathsf{t,dq}}+\mathbf{E}^{\top}v_{\mathsf{dq}}\end{split}$
(4.2)

Closed-loop input comes from Section 3.2:
$\tilde{\tau}_{m}=\tau_{e}-\nabla_{\theta_{\mathsf{dq}}}S_{B}=I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}i_{s,\mathsf{dq}}-\nabla_{\theta_{\mathsf{dq}}}S_{B}$,
where one could use real $\mathbf{E}$ or communication $\mathbf{B}$ graph.
Remaining arguments are as before: driving system convergence to its own
almost globally asymptotically stable invariant set, driven system, containing
now stator dynamics and still being an LTI system under constant input,
converges to angle-dictated steady state. The problem with this approach,
namely real-time angle communication, still persists.

####  4.1.2 Incremental energy function

Define closed-loop input: $\tilde{\tau}_{m}=I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}(Ki_{s,\mathsf{dq}}+i_{s,\mathsf{dq}}^{*im})$,
where
$K=\mathbf{Re}(\Pi_{1})\cdot\Pi_{1}^{-1}=\Pi^{\top}\begin{bmatrix}R_{s}&&\\\
&G&\\\ &&R_{\mathsf{t}}\end{bmatrix}\Pi\cdot\Pi_{1}^{-1}$ and
$i_{s,\mathsf{dq}}^{*im}=\mathbf{Im}(\Pi_{1})\cdot\xi^{*}=\Pi^{\top}\begin{bmatrix}\boldsymbol{j}\omega_{0}L_{s}&&\\\
&-\boldsymbol{j}\omega_{0}C&\\\
&&\boldsymbol{j}\omega_{0}L_{\mathsf{t}}\end{bmatrix}\Pi\cdot\xi^{*}$. Closed-
loop dynamics look like:

$\begin{split}\dot{\theta}_{\mathsf{dq}}&=\tilde{\omega}\\\
M\dot{\tilde{\omega}}&=-D\tilde{\omega}-I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}\big{(}I-K\big{)}\tilde{i}_{s,\mathsf{dq}}-\underbrace{I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}\big{(}\hat{i}_{s,\mathsf{dq}}^{im}-i_{s,\mathsf{dq}}^{*im}\big{)}}_{\text{$\nabla
S$}}\\\
L_{s}\dot{\tilde{i}}_{s,\mathsf{dq}}&=-\mathbf{Z}_{s}\tilde{i}_{s,\mathsf{dq}}+\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}\tilde{\omega}-\tilde{v}_{\mathsf{dq}}-\omega_{0}\boldsymbol{j}L_{s}\Pi_{1}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}\tilde{\omega}\\\
C\dot{\tilde{v}}_{\mathsf{dq}}&=\tilde{i}_{s,\mathsf{dq}}-\mathbf{E}\tilde{i}_{\mathsf{t,dq}}-\mathbf{Y}_{c}\tilde{v}_{\mathsf{dq}}-\omega_{0}\boldsymbol{j}C\Pi_{2}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}\tilde{\omega}\\\
L_{\mathsf{t}}\dot{\tilde{i}}_{\mathsf{t,dq}}&=-\mathbf{Z}_{\mathsf{t}}\tilde{i}_{\mathsf{t,dq}}+\mathbf{E}^{\top}\tilde{v}_{\mathsf{dq}}-\omega_{0}\boldsymbol{j}L_{\mathsf{t}}\Pi_{3}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}\tilde{\omega}\end{split}$ (4.3)

where
$S=\frac{1}{2}\big{(}\hat{\xi}-\xi^{*}\big{)}^{\top}\Pi^{\top}\begin{bmatrix}L_{s}&&\\\
&-C&\\\ &&L_{\mathsf{t}}\end{bmatrix}\Pi\big{(}\hat{\xi}-\xi^{*}\big{)}$ can
be seen as a Lagrangian function. Consider now the following full-system
energy equation:

$\tilde{H}=\frac{1}{2}\tilde{\omega}^{\top}M\tilde{\omega}+\frac{1}{2}\tilde{i}_{s,\mathsf{dq}}^{\top}L_{s}\tilde{i}_{s,\mathsf{dq}}+\frac{1}{2}\tilde{v}_{\mathsf{dq}}^{\top}C\tilde{v}_{\mathsf{dq}}+\frac{1}{2}\tilde{i}_{\mathsf{t,dq}}^{\top}L_{\mathsf{t}}\tilde{i}_{\mathsf{t,dq}}+S(\theta_{\mathsf{dq}})$
(4.4)

Expression for its Lie derivative looks like:

$\begin{split}&\frac{d}{dt}\tilde{H}=\begin{bmatrix}\tilde{\omega}\\\
\tilde{i}_{s,\mathsf{dq}}\\\ \tilde{v}_{\mathsf{dq}}\\\
\tilde{i}_{\mathsf{t,dq}}\end{bmatrix}^{\top}\underbrace{\begin{bmatrix}-D&&&\\\
\Big{[}\big{(}\omega_{0}\boldsymbol{j}L_{s}\big{)}^{\top}\Pi_{1}+K^{\top}\Big{]}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}&-R_{s}&&\\\
\big{(}\omega_{0}\boldsymbol{j}C\big{)}^{\top}\Pi_{2}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}&&-G&\\\
\big{(}\omega_{0}\boldsymbol{j}L_{\mathsf{t}}\big{)}^{\top}\Pi_{3}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}&&&-R_{\mathsf{t}}\\\
\end{bmatrix}}_{\text{-Q}}\begin{bmatrix}\tilde{\omega}\\\
\tilde{i}_{s,\mathsf{dq}}\\\ \tilde{v}_{\mathsf{dq}}\\\
\tilde{i}_{\mathsf{t,dq}}\end{bmatrix}<0.\end{split}$ (4.5)

From equation (LABEL:eq:two-way-one-lie-derivative) one could derive
explicitly the inequality condition on damping term $D$ such that almost
global asymptotic stability is achieved. Note that whenever angle space
penetrates system state-space (which is the case for both reduced and full
model), almost global asymptotic stability can be achieved at most because the
angle space itself is closed and bounded, meaning that it is inherently
impossible to have a radially-unbounded energy function with level sets
projected onto the state-space. Nonetheless, maximum energy points with zero
gradient (hills of energy plot) are unstable solutions, so there is no need to
worry about those in practice.

The problem with defined closed-loop input is that matrix $K$ is not block-
diagonal, meaning that it is impossible to create such input using local
stator current measurements only. Difficulties in tackling full model directly
motivate reviewing reduced model decentralized control.

###  4.2 Reduced model decentralized control

This section reviews the decentralized controller proposed recently by authors
of [10]. Reduced model open-loop dynamics in $\mathsf{dq}$-frame are the same
as in (3.2):

$\begin{split}\dot{\theta}_{\mathsf{dq}}&=\tilde{\omega}\\\
M\dot{\tilde{\omega}}&=-D\tilde{\omega}+\tilde{\tau}_{m}-I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}\mathbf{E}i_{\mathsf{t,dq}}\\\
L_{\mathsf{t}}\dot{i}_{\mathsf{t,dq}}&=-\mathbf{Z}_{\mathsf{t}}i_{\mathsf{t,dq}}+\mathbf{E}^{\top}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}(\tilde{\omega}+\omega_{0}\boldsymbol{1})\end{split}$ (4.6)

Proposing the following closed-loop input in $\mathsf{dq}$-frame:
$\tilde{\tau}_{m}=I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}\mathbf{E}\mathbf{Z}_{\mathsf{t}}^{-\top}R_{\mathsf{t}}i_{\mathsf{t,dq}}$,
one could arrive at the following closed-loop system dynamics in
$\mathsf{dq}$:

$\begin{split}\dot{\theta}_{\mathsf{dq}}&=\tilde{\omega}\\\
M\dot{\tilde{\omega}}&=-D\tilde{\omega}+I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}\mathbf{E}\mathbf{Z}_{\mathsf{t}}^{-\top}\boldsymbol{j}\omega_{0}L_{\mathsf{t}}\tilde{i}_{\mathsf{t,dq}}\\\
&+I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}\mathbf{E}\mathbf{Z}_{\mathsf{t}}^{-\top}\boldsymbol{j}\omega_{0}L_{\mathsf{t}}\mathbf{Z}_{\mathsf{t}}^{-1}\mathbf{E}^{\top}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}\omega_{0}\boldsymbol{1}=\\\
&=-D\tilde{\omega}+I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}\mathbf{E}\mathbf{Z}_{\mathsf{t}}^{-\top}\boldsymbol{j}\omega_{0}L_{\mathsf{t}}\tilde{i}_{\mathsf{t,dq}}-\nabla_{\theta_{\mathsf{dq}}}S\\\
L_{\mathsf{t}}\dot{\tilde{i}}_{\mathsf{t,dq}}&=-\mathbf{Z}_{\mathsf{t}}\tilde{i}_{\mathsf{t,dq}}+R_{\mathsf{t}}\mathbf{Z}_{\mathsf{t}}^{-1}\mathbf{E}^{\top}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}\tilde{\omega}\end{split}$ (4.7)

Here
$\nabla_{\theta_{\mathsf{dq}}}S=\frac{\partial^{\top}}{\partial\theta_{\mathsf{dq}}}(\frac{1}{2}\hat{i}_{\mathsf{t,dq}}^{\top}L_{\mathsf{t}}\hat{i}_{\mathsf{t,dq}})$.
Now, due to assumption of $\frac{L_{\mathsf{t}}}{R_{\mathsf{t}}}=const$ the
proposed closed-loop input is decentralized. Indeed, thanks to aforementioned
assumption, a constant block-diagonal matrix
$\mathbf{Z}_{\mathsf{t}}^{-\top}\boldsymbol{j}\omega_{0}L_{\mathsf{t}}$ now
commutes with incidence matrix $\mathbf{E}$ in (4.7), which is quite useful,
since the resulting $\mathbf{E}i_{\mathsf{t,dq}}$ is nothing else that just a
vector of current injections into the network, readily available for
measurements locally. The closed-loop dynamics above are almost globally
asymptotically stable, which can be proven by LaSalle argument in the
following manner:

$\tilde{H}=\frac{1}{2}\tilde{\omega}^{\top}M\tilde{\omega}+\frac{1}{2}\tilde{i}_{\mathsf{t,dq}}^{\top}L_{\mathsf{t}}\tilde{i}_{\mathsf{t,dq}}+S(\theta_{\mathsf{dq}})$
(4.8)

Expression for proposed energy function’s Lie derivative looks like:

$\begin{split}&\frac{d}{dt}\tilde{H}=-\tilde{\omega}^{\top}D\tilde{\omega}-\tilde{\omega}^{\top}\nabla_{\theta_{\mathsf{dq}}}S+\tilde{\omega}^{\top}I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}\mathbf{E}\mathbf{Z}_{\mathsf{t}}^{-\top}\boldsymbol{j}\omega_{0}L_{\mathsf{t}}\tilde{i}_{\mathsf{t,dq}}-\tilde{i}_{\mathsf{t,dq}}^{\top}R_{\mathsf{t}}\tilde{i}_{\mathsf{t,dq}}+\\\
&+\tilde{i}_{\mathsf{t,dq}}^{\top}R_{\mathsf{t}}\mathbf{Z}_{\mathsf{t}}^{-1}\mathbf{E}^{\top}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}\tilde{\omega}+(\nabla_{\theta_{\mathsf{dq}}}S)^{\top}\tilde{\omega}=\\\
&=\begin{bmatrix}\tilde{\omega}\\\
\tilde{i}_{\mathsf{t,dq}}\end{bmatrix}^{\top}\underbrace{\begin{bmatrix}-D&I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}\mathbf{E}\mathbf{Z}_{\mathsf{t}}^{-\top}\boldsymbol{j}\omega_{0}L_{\mathsf{t}}\\\
R_{\mathsf{t}}\mathbf{Z}_{\mathsf{t}}^{-1}\mathbf{E}^{\top}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}&-R_{\mathsf{t}}\end{bmatrix}}_{\text{-Q}}\begin{bmatrix}\tilde{\omega}\\\
\tilde{i}_{\mathsf{t,dq}}\end{bmatrix}<0\end{split}$ (4.9)

From here, based on $(Q+Q^{\top})/2\succ 0$ condition, the explicit inequality
for damping $D$ that ensures almost global asymptotic stability and should be
satisfied for all $\theta_{\mathsf{dq}}$, can be derived:

$D\succ I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}\mathbf{E}\frac{R_{\mathsf{t}}^{-1}}{4}\mathbf{E}^{\top}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}$ (4.10)

Proposed controller is of an incremental energy function-type. It stabilizes
the consensus in angle space. To shift the desired references an additional
locally-measurable term should be added to the closed-loop input:
$\tilde{\tau}_{m}=I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}\mathbf{E}\mathbf{Z}_{\mathsf{t}}^{-\top}R_{\mathsf{t}}i_{\mathsf{t,dq}}-I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}\mathbf{E}\mathbf{Z}_{\mathsf{t}}^{-\top}\boldsymbol{j}\omega_{0}L_{\mathsf{t}}\mathbf{Z}_{\mathsf{t}}^{-1}\mathbf{E}^{\top}\mathbf{R}_{\theta_{\mathsf{dq}}^{*}}(L_{m}\otimes
e_{2})I_{r}^{*}\omega_{0}\boldsymbol{1}$. Now the shifted gradient and
associated energy function would look like:
$\nabla_{\theta_{\mathsf{dq}}}S=\frac{\partial^{\top}}{\partial\theta_{\mathsf{dq}}}(\frac{1}{2}(\hat{i}_{\mathsf{t,dq}}-i_{\mathsf{t,dq}}^{*})^{\top}L_{\mathsf{t}}(\hat{i}_{\mathsf{t,dq}}-i_{\mathsf{t,dq}}^{*}))$.
Since this energy function does not enter the equation on full energy
derivative (LABEL:eq:_reduced-model-lie-derivative1), inequality condition on
damping term $D$ remains unchanged.

Note that $S(\theta_{\mathsf{dq}})$ has the structure of second energy
function $\tilde{S}_{B}$ in our topology-independence discussion. From
structural similarity one could argue that virtual energy function
$\tilde{S}_{B}$ back in Chapter 3 had the interpretation of transmission line
inductive energy.

###  4.3 Simplified full model decentralized control

The idea behind this section is to simplify the full model to the level of
reduced one from where a modified version of last section’s controller could
be applied. It turns out that simplification is possible under the following
three conditions, which will become evident in further analysis:

$\bullet$ Static stator assumption or sufficient time scale separation
($\tau_{stator}<<\tau_{line}$)

$\bullet$ $L_{\mathsf{t}}/R_{\mathsf{t}}=const$ for transmission lines as in
Chapter 2

$\bullet$ Block diagonal dominance of $\Pi_{2}$-map (or $\Pi_{1}$)

We first review the full model in both stationary and rotating frames:

$\begin{cases}\begin{split}\dot{\theta}&=\omega\\\
M\dot{\omega}&=-D\omega+\tau_{m}+\tau_{e}\\\
\dot{\lambda}_{r}&=-R_{r}i_{r}+u_{r}\\\ \dot{\lambda}_{s}&=-R_{s}i_{s}+v\\\
C\dot{v}&=-Gv-i_{s}-\mathbf{E}i_{\mathsf{t}}\\\
L_{\mathsf{t}}\dot{i}_{\mathsf{t}}&=-R_{\mathsf{t}}i_{\mathsf{t}}+\mathbf{E}^{\top}v\end{split}\end{cases}\Leftrightarrow\begin{cases}\begin{split}\dot{\theta}_{\mathsf{dq}}&=\tilde{\omega}\\\
M\dot{\tilde{\omega}}&=-D\tilde{\omega}+\tilde{\tau}_{m}+\tau_{e}\\\
\dot{\lambda}_{r}&=-R_{r}i_{r}+u_{r}\\\
\dot{\lambda}_{s,\mathsf{dq}}&=-\mathbf{Z}_{s}i_{s,\mathsf{dq}}+v_{\mathsf{dq}}-\omega_{0}\boldsymbol{j}\Psi_{s,\mathsf{dq}}\\\
C\dot{v}_{\mathsf{dq}}&=-\mathbf{Y}_{c}v_{\mathsf{dq}}-i_{s,\mathsf{dq}}-\mathbf{E}i_{\mathsf{t,dq}}\\\
L_{\mathsf{t}}\dot{i}_{\mathsf{t,dq}}&=-\mathbf{Z}_{\mathsf{t}}i_{\mathsf{t,dq}}+\mathbf{E}^{\top}v_{\mathsf{dq}}\end{split}\end{cases}$
(4.11)

Here we have used the following notation for angle embedding
$\Psi_{s,\mathsf{dq}}=\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes\mathbf{e}_{1})i_{r}$.
Under assumption of static stator, simplified model dynamics in dq-reference
frame can be written as:

$\begin{split}\dot{\theta}_{\mathsf{dq}}&=\tilde{\omega}\\\
M\dot{\tilde{\omega}}&=-D\tilde{\omega}+\tilde{\tau}_{m}-I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}i_{s,\mathsf{dq}}\\\
\mathbf{Z}_{s}i_{s,\mathsf{dq}}&=\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}(\tilde{\omega}+\omega_{0}\boldsymbol{1})-v_{\mathsf{dq}}\\\
\mathbf{Y}_{c}v_{\mathsf{dq}}&=i_{s,\mathsf{dq}}-\mathbf{E}i_{\mathsf{t,dq}}\\\
L_{\mathsf{t}}\dot{i}_{\mathsf{t,dq}}&=-\mathbf{Z}_{\mathsf{t}}i_{\mathsf{t,dq}}+\mathbf{E}^{\top}v_{\mathsf{dq}}\end{split}$
(4.12)

As it can be seen from equation (4.13) static stator assumption eliminates
respective inductance and capacitance dynamics, leaving their steady state
equations only. Expressing $i_{s,\mathsf{dq}}$ and $v_{\mathsf{dq}}$ in terms
of $i_{\mathsf{t,dq}}$, system (4.13) can be reduced to:

$\begin{split}\dot{\theta}_{\mathsf{dq}}&=\tilde{\omega}\\\
M\dot{\tilde{\omega}}&=-D\tilde{\omega}+\tilde{\tau}_{m}-I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}(\mathbf{Z}_{s}+\mathbf{Y}_{c}^{-1})^{-1}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}(\tilde{\omega}+\omega_{0}\boldsymbol{1})\\\
&-I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}(I+\mathbf{Y}_{c}\mathbf{Z}_{s})^{-1}\mathbf{E}i_{\mathsf{t,dq}}\\\
L_{\mathsf{t}}\dot{i}_{\mathsf{t,dq}}&=-\mathbf{Z}_{\mathsf{t}}i_{\mathsf{t,dq}}+\mathbf{E}^{\top}(I+\mathbf{Y}_{c}\mathbf{Z}_{s})^{-1}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}(\tilde{\omega}+\omega_{0}\boldsymbol{1})-\mathbf{E}^{\top}(\mathbf{Z}_{s}^{-1}+\mathbf{Y}_{c})^{-1}\mathbf{E}i_{\mathsf{t,dq}}\end{split}$
(4.13)

Recall that mechanical torque input previously was defined as
$\tau_{m}=D\omega_{0}+\tilde{\tau}_{m}$. In reduced static stator model (4.16)
a stator-dependent damping term allows redefinition of input as follows:
$\tau_{m}=D^{\prime}\omega_{0}+\tilde{\tau}_{m}$. New damping factor looks
like

$D^{\prime}=D+\mathsf{diag}\bigg{[}(L_{m}i_{r}^{*})^{2}(G+\Arrowvert
Y_{c}\Arrowvert^{2}R_{s})\frac{\Arrowvert Y_{c}\Arrowvert^{2}}{\Arrowvert
A\Arrowvert^{2}}\bigg{]}_{k}$ (4.14)

where $\Arrowvert Y_{c}\Arrowvert^{2}=(G^{2}+\omega_{0}^{2}C^{2})$ and
$\Arrowvert A\Arrowvert^{2}=(G+R_{s}\Arrowvert
Y_{c}\Arrowvert^{2})^{2}+\omega_{0}^{2}(C-\Arrowvert
Y_{c}\Arrowvert^{2}L_{s})^{2}$ are scalar values corresponding to every
separate node. Moreover, transmission line impedance matrix
$\mathbf{Z}_{\mathsf{t}}$ can also be redefined as:

$\mathbf{Z}_{\mathsf{t}}^{\prime}=\mathbf{Z}_{\mathsf{t}}+\mathbf{E}^{\top}(\mathbf{Z}_{s}^{-1}+\mathbf{Y}_{c})^{-1}\mathbf{E}$
(4.15)

With these in mind, reduced model simplifies even further:

$\begin{split}\dot{\theta}_{\mathsf{dq}}&=\tilde{\omega}\\\
M\dot{\tilde{\omega}}&=-D^{\prime}\tilde{\omega}+\tilde{\tau}_{m}-I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}(I+\mathbf{Y}_{c}\mathbf{Z}_{s})^{-1}\mathbf{E}i_{\mathsf{t,dq}}\\\
L_{\mathsf{t}}\dot{i}_{\mathsf{t,dq}}&=-\mathbf{Z}_{\mathsf{t}}^{\prime}i_{\mathsf{t,dq}}+\mathbf{E}^{\top}(I+\mathbf{Y}_{c}\mathbf{Z}_{s})^{-1}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}(\tilde{\omega}+\omega_{0}\boldsymbol{1})\end{split}$ (4.16)

As a next step, recall the error-coordinate relation
$\tilde{i}_{\mathsf{t,dq}}=i_{\mathsf{t,dq}}-\hat{i}_{\mathsf{t,dq}}(\theta_{\mathsf{dq}})$
and high-order model steady-state map
$\hat{i}_{\mathsf{t,dq}}=\Pi_{3}\xi=\mathbf{Z}_{\mathsf{t}}^{-1}\mathbf{E}^{\top}\Pi_{2}\xi$.
Using these, one can derive the derivative of steady-state current
$\dot{\hat{i}}_{\mathsf{t,dq}}=\mathbf{Z}_{\mathsf{t}}^{-1}\mathbf{E}^{\top}\Pi_{2}\omega_{0}\boldsymbol{j}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes\mathbf{e}_{2})I_{r}^{*}\tilde{\omega}$
and rewrite (4.16) as:

$\begin{split}\dot{\theta}_{\mathsf{dq}}&=\tilde{\omega}\\\
M\dot{\tilde{\omega}}&=-D^{\prime}\tilde{\omega}+\tilde{\tau}_{m}-I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}(I+\mathbf{Y}_{c}\mathbf{Z}_{s})^{-1}\mathbf{E}i_{\mathsf{t,dq}}\\\
L_{\mathsf{t}}\dot{\tilde{i}}_{\mathsf{t,dq}}&=-\mathbf{Z}_{\mathsf{t}}^{\prime}\tilde{i}_{\mathsf{t,dq}}+\big{(}\mathbf{E}^{\top}(I+\mathbf{Y}_{c}\mathbf{Z}_{s})^{-1}-\omega_{0}\boldsymbol{j}L_{\mathsf{t}}\mathbf{Z}_{\mathsf{t}}^{-1}\mathbf{E}^{\top}\Pi_{2}\big{)}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}\tilde{\omega}\end{split}$ (4.17)

where one used steady state equilibrium condition
$\mathbf{Z}_{\mathsf{t}}^{\prime}\hat{i}_{\mathsf{t,dq}}=\mathbf{E}^{\top}(I+\mathbf{Y}_{c}\mathbf{Z}_{s})^{-1}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}\omega_{0}$. To achieve angle consensus in closed-loop, the
following modified control law is suggested:

$\tilde{\tau}_{m}=I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}(I+\mathbf{Y}_{c}\mathbf{Z}_{s})^{-1}\mathbf{E}i_{\mathsf{t,dq}}\
+\ I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}\Pi_{2}^{\top}\mathbf{E}\mathbf{Z}_{\mathsf{t}}^{-\top}\omega_{0}\boldsymbol{j}L_{\mathsf{t}}i_{\mathsf{t,dq}}$
(4.18)

The first term fully cancels the electrical torque, while the second takes
care of the negative gradient term. Note that the first term can be generated
based on local line current measurements only, thanks to block-diagonal
structure of $(I+\mathbf{Y}_{c}\mathbf{Z}_{s})^{-1}$. The second term’ local
estimation will be addressed later. Proposed control law allows rewriting
(4.17) in closed-loop form as follows:

$\begin{split}\dot{\theta}_{\mathsf{dq}}&=\tilde{\omega}\\\
M\dot{\tilde{\omega}}&=-D^{\prime}\tilde{\omega}+I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}\Pi_{2}^{\top}\mathbf{E}\mathbf{Z}_{\mathsf{t}}^{-\top}\omega_{0}\boldsymbol{j}L_{\mathsf{t}}\tilde{i}_{\mathsf{t,dq}}-\nabla_{\theta_{\mathsf{dq}}}S\\\
L_{\mathsf{t}}\dot{\tilde{i}}_{\mathsf{t,dq}}&=-\mathbf{Z}_{\mathsf{t}}^{\prime}\tilde{i}_{\mathsf{t,dq}}+\big{(}\mathbf{E}^{\top}(I+\mathbf{Y}_{c}\mathbf{Z}_{s})^{-1}-\omega_{0}\boldsymbol{j}L_{\mathsf{t}}\mathbf{Z}_{\mathsf{t}}^{-1}\mathbf{E}^{\top}\Pi_{2}\big{)}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}\tilde{\omega}\end{split}$ (4.19)

where the gradient term is the same as in previous reduced model controller:

$\nabla_{\theta_{\mathsf{dq}}}S=-\bar{I}_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\bar{\mathbf{R}}_{\theta_{\mathsf{dq}}}^{\top}\mathbf{E}\mathbf{Z}_{\mathsf{t}}^{-\top}\omega_{0}\boldsymbol{j}L_{\mathsf{t}}\mathbf{Z}_{\mathsf{t}}^{-1}\mathbf{E}^{\top}\bar{\mathbf{R}}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})\bar{I}_{r}^{*}\omega_{0}\mathbf{1}=\frac{\partial^{\top}}{\partial\theta_{\mathsf{dq}}}(\frac{1}{2}\hat{i}_{\mathsf{t,dq}}^{\top}L_{\mathsf{t}}\hat{i}_{\mathsf{t,dq}})$
(4.20)

A subtle difference between previous section’s controller and present one is
in stator presence. Indeed, if in previous model node electro-motive forces
$\xi=\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}\omega_{0}\mathbf{1}$ influenced transmission lines directly,
now they do it through stator, which has a steady-state map of $\Pi_{2}$:
$\xi\rightarrow v_{\mathsf{dq}}$. Therefore, effective rotation map
$\bar{\mathbf{R}}_{\theta_{\mathsf{dq}}}$ and rotor currents
$\bar{I}_{r}^{*}$, which correspond to $v_{\mathsf{dq}}$, must be used in
gradient. These effective quantities can be interpreted as the amplitudes and
dq-angles of node voltages $v_{\mathsf{dq}}$ throughout the graph.
Corresponding analytical relation is:
$\Pi_{2}\xi=\Pi_{2}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}\omega_{0}\mathbf{1}=\bar{\mathbf{R}}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})\bar{I}_{r}^{*}\omega_{0}\mathbf{1}=v_{\mathsf{dq}}$.

With this, local generation of the second term in closed-loop input (4.18) can
be addressed. Notice
$\hat{i}_{\mathsf{t,dq}}^{\top}/\omega_{0}=\mathbf{1}^{\top}I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}\Pi_{2}^{\top}\mathbf{E}\mathbf{Z}_{\mathsf{t}}^{-\top}=\mathbf{1}^{\top}\bar{I}_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\bar{\mathbf{R}}_{\theta_{\mathsf{dq}}}^{\top}\mathbf{E}\mathbf{Z}_{\mathsf{t}}^{-\top}=\mathbf{1}^{\top}\mathsf{diag}^{\top}(v_{\mathsf{dq}}/\omega_{0})\mathbf{E}\mathbf{Z}_{\mathsf{t}}^{-\top}$
with $\hat{i}_{\mathsf{t,dq}}\in^{2n}$ and
$\mathsf{diag}(v_{\mathsf{dq}}/\omega_{0})\in^{2n\times n}$, meaning that
second term in input can be approximated as:
$\mathsf{diag}^{\top}(v_{\mathsf{dq}}/\omega_{0})\mathbf{E}\mathbf{Z}_{\mathsf{t}}^{-\top}\omega_{0}\boldsymbol{j}L_{\mathsf{t}}i_{\mathsf{t,dq}}$
under diagonal dominance of $\Pi_{2}$. Thus, it also can be generated locally
with node voltage and transmission line current measurements, thanks to
assumption of constant $R_{\mathsf{t}}/L_{\mathsf{t}}$ ratio, which allows
commuting $\mathbf{E}$ and
$\mathbf{Z}_{\mathsf{t}}^{-\top}\omega_{0}\boldsymbol{j}L_{\mathsf{t}}$, and
special structure of $\Pi_{2}$, allowing block-diagonal form of
$\mathsf{diag}(v_{\mathsf{dq}}/\omega_{0})$.

Consider now the following full energy function:

$\tilde{H}=\frac{1}{2}\tilde{\omega}^{\top}M\tilde{\omega}+\frac{1}{2}\tilde{i}_{\mathsf{t,dq}}^{\top}L_{\mathsf{t}}\tilde{i}_{\mathsf{t,dq}}+S(\theta_{\mathsf{dq}})$
(4.21)

Expression for its Lie derivative looks like:

$\begin{split}&\frac{d}{dt}\tilde{H}=-\tilde{\omega}^{\top}D^{\prime}\tilde{\omega}-\tilde{\omega}^{\top}\nabla_{\theta_{\mathsf{dq}}}S+\tilde{\omega}^{\top}I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}\Pi_{2}^{\top}\mathbf{E}\mathbf{Z}_{\mathsf{t}}^{-\top}\boldsymbol{j}\omega_{0}L_{\mathsf{t}}\tilde{i}_{\mathsf{t,dq}}-\\\
&-\tilde{i}_{\mathsf{t,dq}}^{\top}\mathbf{Z}_{\mathsf{t}}^{\prime}\tilde{i}_{\mathsf{t,dq}}+\tilde{i}_{\mathsf{t,dq}}^{\top}\big{(}\mathbf{E}^{\top}(I+\mathbf{Y}_{c}\mathbf{Z}_{s})^{-1}-\omega_{0}\boldsymbol{j}L_{\mathsf{t}}\mathbf{Z}_{\mathsf{t}}^{-1}\mathbf{E}^{\top}\Pi_{2}\big{)}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}\tilde{\omega}+(\nabla_{\theta_{\mathsf{dq}}}S)^{\top}\tilde{\omega}=\\\
&=\begin{bmatrix}\tilde{\omega}\\\
\tilde{i}_{\mathsf{t,dq}}\end{bmatrix}^{\top}\underbrace{\begin{bmatrix}-D^{\prime}&2I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}\Pi_{2}^{\top}\mathbf{E}\mathbf{Z}_{\mathsf{t}}^{-\top}\boldsymbol{j}\omega_{0}L_{\mathsf{t}}\\\
\mathbf{E}^{\top}(I+\mathbf{Y}_{c}\mathbf{Z}_{s})^{-1}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}&-\mathbf{Z}_{\mathsf{t}}^{\prime}\end{bmatrix}}_{\text{-Q}}\begin{bmatrix}\tilde{\omega}\\\
\tilde{i}_{\mathsf{t,dq}}\end{bmatrix}<0\end{split}$ (4.22)

For $\frac{d}{dt}\tilde{H}$ to be negative definite we require that
$A=\frac{1}{2}(Q+Q^{\top})$ be positive definite, which, recalling equation
(4.15), is equivalent to:

$\begin{split}D^{\prime}\succ I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}&\Big{(}\Pi_{2}^{\top}\mathbf{Z}_{\mathsf{t}}^{-\top}\boldsymbol{j}\omega_{0}L_{\mathsf{t}}+\frac{1}{2}(I+\mathbf{Y}_{c}\mathbf{Z}_{s})^{-\top}\Big{)}\mathbf{E}\Big{(}\mathbf{Z}_{\mathsf{t}}+\mathbf{E}^{\top}(\mathbf{Z}_{s}^{-1}+\mathbf{Y}_{c})^{-1}\mathbf{E}\Big{)}^{-1}\mathbf{E}^{\top}\\\
&\Big{(}\frac{1}{2}(I+\mathbf{Y}_{c}\mathbf{Z}_{s})^{-1}+\omega_{0}\boldsymbol{j}^{\top}L_{\mathsf{t}}\mathbf{Z}_{\mathsf{t}}^{-1}\Pi_{2}\Big{)}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}\end{split}$ (4.23)

or, using the relation between $D$ and $D^{\prime}$ (4.14), equivalent to

$\begin{split}D&\succ I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}\Big{(}\Pi_{2}^{\top}\mathbf{Z}_{\mathsf{t}}^{-\top}\boldsymbol{j}\omega_{0}L_{\mathsf{t}}+\frac{1}{2}(I+\mathbf{Y}_{c}\mathbf{Z}_{s})^{-\top}\Big{)}\mathbf{E}\Big{(}\mathbf{Z}_{\mathsf{t}}+\mathbf{E}^{\top}(\mathbf{Z}_{s}^{-1}+\mathbf{Y}_{c})^{-1}\mathbf{E}\Big{)}^{-1}\mathbf{E}^{\top}\\\
&\Big{(}\frac{1}{2}(I+\mathbf{Y}_{c}\mathbf{Z}_{s})^{-1}+\omega_{0}\boldsymbol{j}^{\top}L_{\mathsf{t}}\mathbf{Z}_{\mathsf{t}}^{-1}\Pi_{2}\Big{)}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}-\mathsf{diag}\bigg{[}(L_{m}i_{r}^{*})^{2}(G+\Arrowvert
Y_{c}\Arrowvert^{2}R_{s})\frac{\Arrowvert Y_{c}\Arrowvert^{2}}{\Arrowvert
A\Arrowvert^{2}}\bigg{]}_{k}\end{split}$ (4.24)

Note that in equations (4.23-4.24) $\mathbf{E}$ and
$\mathbf{Z}_{\mathsf{t}}^{-\top}\boldsymbol{j}\omega_{0}L_{\mathsf{t}}$
commute due to $L_{\mathsf{t}}/R_{\mathsf{t}}=\alpha$ assumption. Moreover,
the following relation holds:
$\mathbf{Z}_{\mathsf{t}}^{-\top}\omega_{0}^{2}L_{\mathsf{t}}^{2}R_{\mathsf{t}}^{-1}\mathbf{Z}_{\mathsf{t}}^{-1}=\frac{\omega_{0}^{2}}{\omega_{0}^{2}+\alpha^{2}}R_{\mathsf{t}}^{-1}$.
For shifted angle consensus in closed-loop, the following control law is
suggested:

$\begin{split}\tilde{\tau}_{m}&=I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}(I+\mathbf{Y}_{c}\mathbf{Z}_{s})^{-1}\mathbf{E}i_{\mathsf{t,dq}}+I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}\Pi_{2}^{\top}\mathbf{E}\mathbf{Z}_{\mathsf{t}}^{-\top}\omega_{0}\boldsymbol{j}L_{\mathsf{t}}(i_{\mathsf{t,dq}}-i_{\mathsf{t,dq}}^{*})\\\
&\approx I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}(I+\mathbf{Y}_{c}\mathbf{Z}_{s})^{-1}\mathbf{E}i_{\mathsf{t,dq}}+\mathsf{diag}^{\top}(v_{\mathsf{dq}}/\omega_{0})\mathbf{E}\mathbf{Z}_{\mathsf{t}}^{-\top}\omega_{0}\boldsymbol{j}L_{\mathsf{t}}(i_{\mathsf{t,dq}}-i_{\mathsf{t,dq}}^{*})\end{split}$
(4.25)

Here one adds an imaginary term depending on desired transmission line
currents to complete the gradient of a shifted energy function, just as in
case of reduced model’s controller. Again, thanks to
$L_{\mathsf{t}}/R_{\mathsf{t}}=\alpha$ assumption,
$\mathbf{Z}^{-\top}\omega_{0}\boldsymbol{j}L_{\mathsf{t}}$ is a constant
block-diagonal matrix, which commutes with $\mathbf{E}$, meaning that this
additional input term can also be generated using local current and voltage
information only. Closed-loop dynamics remain same as in (4.19), where

$\begin{split}\nabla_{\theta_{\mathsf{dq}}}S&=-\bar{I}_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\bar{\mathbf{R}}_{\theta_{\mathsf{dq}}}^{\top}\mathbf{E}\mathbf{Z}_{\mathsf{t}}^{-\top}\boldsymbol{j}\omega_{0}L_{\mathsf{t}}\mathbf{Z}_{\mathsf{t}}^{-1}\mathbf{E}^{\top}(\bar{\mathbf{R}}_{\theta_{\mathsf{dq}}}-\bar{\mathbf{R}}_{\theta_{\mathsf{dq}}^{*}})(L_{m}\otimes
e_{2})\bar{I}_{r}^{*}\omega_{0}\mathbf{1}\\\
&=\frac{\partial^{\top}}{\partial\theta_{\mathsf{dq}}}\big{(}\frac{1}{2}(\hat{i}_{\mathsf{t,dq}}-{i}_{\mathsf{t,dq}}^{*})^{\top}L_{\mathsf{t}}(\hat{i}_{\mathsf{t,dq}}-{i}_{\mathsf{t,dq}}^{*})\big{)}\end{split}$
(4.26)

Since gradient-dependent terms in Lie derivative (LABEL:eq:_energy-function-
lie-derivative) cancel, conditions on $D$ and $D^{\prime}$ remain unchanged.

As a miscellaneous remark, to ensure inequality condition on damping term $D$
for all possible angle combinations during transients, a proportional
controller $K_{p}$ could be added in $\mathsf{dq}$-frame to artificially rise
the damping without increasing steady state loss $D\omega_{0}$, acting as a
rotor friction in case of SG and as a DC-link parasitic resistance plus
switching losses in case of inverter.

Several Simulink realizations on a triangular network (Fig.3.3) with
literature-taken specifications show that under sufficiently high $K_{p}$ and
start near steady state one can stabilize the desired OPF set-points with
proposed controller acting on a full model. Time scale separation requires
$\tau_{lines}/\tau_{stator}$ ratio to be about 25 for all possible transient
angle configurations and as small as 8 for some angles. However, there are
many parameters to play with, so these numbers are by no means conclusive.
More precise conclusions require more simulations, but the message to be drawn
from these scarce figures already is that times scale separation does not need
to be large to the extent of static stator approximation in order for
simplified controller to stabilize full model.

From a higher perspective, simplified model controller has several problems.
First, there are three underlying assumptions for it to work, which may not be
practical in reality. Moreover, even if all the assumptions are satisfied,
gradient-creating energy function of such form almost globally stabilizes a
disjoint invariant set in angle space, just as in case of second energy
function in out topology-independence discussion, meaning that depending on
initial conditions one could end up in wrong ”valley”, stabilizing undesired
angle configuration and, consequently, wrong power flow. The same problem also
emerges in reduced model decentralized controller of previous section, which
is not surprising, because both use the same energy function and, in fact, one
is a slight modification of another. These unacceptable deficiencies motivate
returning back to full model altogether. This ”journey” to reduced model
decentralized control, however, should not be seen as useless, since it turned
out to be an intermediate step in full model’s decentralized control design,
addressed in the next section.

###  4.4 Full model decentralized control

In this last section a decentralized controller for a full model is developed.
First, we review the dynamic equations in $\mathsf{dq}$-frame and error
coordinates for rotor dynamics:

$\begin{split}\dot{\theta}_{\mathsf{dq}}&=\tilde{\omega}\\\
M\dot{\tilde{\omega}}&=-D\tilde{\omega}+\tilde{\tau}_{m}-I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}i_{s,\mathsf{dq}}\\\
L_{s}\dot{i}_{s,\mathsf{dq}}&=-\mathbf{Z}_{s}i_{s,\mathsf{dq}}+\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}(\tilde{\omega}+\omega_{0}\boldsymbol{1})-v_{\mathsf{dq}}\\\
C\dot{v}_{\mathsf{dq}}&=i_{s,\mathsf{dq}}-\mathbf{E}i_{\mathsf{t,dq}}-\mathbf{Y}_{c}v_{\mathsf{dq}}\\\
L_{\mathsf{t}}\dot{i}_{\mathsf{t,dq}}&=-\mathbf{Z}_{\mathsf{t}}i_{\mathsf{t,dq}}+\mathbf{E}^{\top}v_{\mathsf{dq}}\end{split}$
(4.27)

Regarding closed-loop input, a good choice of gradient-creating energy
function would preserve coupling, yet avoid non-diagonal matrices as middle
arguments to energy function. Moreover, bijective local measurements (such as
stator currents and node voltages, not transmission line currents) should be
used, because, firstly, these are the only local information sources and
network coupling reflections available, and secondly, bijectivity ensures that
incidence matrix $\mathbf{E}$ will not enter the gradient and create disjoint
invariant sets in angle space (the problem of reduced and simplified-full
model controllers). Here a degree of freedom emerges in the form of choosing
between stator currents and node voltages. From practical perspective, bus
voltages are more important to control and easier to measure/process, so
further work will be continued with those. Another approach would be to create
both current and voltage-related gradients. This option, however, was not
investigated because of possible growth of invariant set in angle space.
Indeed, combined gradient could be zero in two cases: when both its current
and voltage components are zero and when they are counteracting each other in
some directions. This option should not be forgotten though. A suitable
voltage-based candidate fulfilling all of these is given below:

$\begin{split}\tilde{\tau}_{m}&=I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}i_{s,\mathsf{dq}}+I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}\Pi_{2}^{\top}\omega_{0}\boldsymbol{j}C(v_{\mathsf{dq}}-v_{\mathsf{dq}}^{*})\\\
&\approx I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}i_{s,\mathsf{dq}}+\mathsf{diag}^{\top}(v_{\mathsf{dq}}/\omega_{0})\omega_{0}\boldsymbol{j}C(v_{\mathsf{dq}}-v_{\mathsf{dq}}^{*})\end{split}$
(4.28)

Notice that there is no incidence matrix $\mathbf{E}$ in the input, thus no
need for constant matrix commuting, thus no need in
$R_{\mathsf{t}}/L_{\mathsf{t}}=const$. Bearing in mind that now we have
dynamic (not static) stator, only a diagonal dominance assumption of $\Pi_{2}$
persists, which was already used in equation (4.28). Closed-loop now looks
like:

$\begin{split}\dot{\theta}_{\mathsf{dq}}&=\tilde{\omega}\\\
M\dot{\tilde{\omega}}&=-D\tilde{\omega}+I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}\Pi_{2}^{\top}\omega_{0}\boldsymbol{j}C\tilde{v}_{\mathsf{dq}}-\nabla_{\theta_{\mathsf{dq}}}S\\\
L_{s}\dot{\tilde{i}}_{s,\mathsf{dq}}&=-\mathbf{Z}_{s}\tilde{i}_{s,\mathsf{dq}}+\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}\tilde{\omega}-\tilde{v}_{\mathsf{dq}}-\omega_{0}\boldsymbol{j}L_{s}\Pi_{1}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}\tilde{\omega}\\\
C\dot{\tilde{v}}_{\mathsf{dq}}&=\tilde{i}_{s,\mathsf{dq}}-\mathbf{E}\tilde{i}_{\mathsf{t,dq}}-\mathbf{Y}_{c}\tilde{v}_{\mathsf{dq}}-\omega_{0}\boldsymbol{j}C\Pi_{2}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}\tilde{\omega}\\\
L_{\mathsf{t}}\dot{\tilde{i}}_{\mathsf{t,dq}}&=-\mathbf{Z}_{\mathsf{t}}\tilde{i}_{\mathsf{t,dq}}+\mathbf{E}^{\top}\tilde{v}_{\mathsf{dq}}-\omega_{0}\boldsymbol{j}L_{\mathsf{t}}\Pi_{3}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}\tilde{\omega}\end{split}$ (4.29)

where

$\begin{split}\nabla_{\theta_{\mathsf{dq}}}S&=-I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}\Pi_{2}^{\top}\omega_{0}\boldsymbol{j}C\Pi_{2}\big{(}\mathbf{R}_{\theta_{\mathsf{dq}}}-\mathbf{R}_{\theta_{\mathsf{dq}}^{*}}\big{)}(L_{m}\otimes
e_{2})I_{r}^{*}\omega_{0}\boldsymbol{1}\\\
\nabla_{\theta_{\mathsf{dq}}}S&=\frac{\partial^{\top}}{\partial\theta_{\mathsf{dq}}}\big{(}\frac{1}{2}(\hat{v}_{\mathsf{dq}}-{v}_{\mathsf{dq}}^{*})^{\top}C(\hat{v}_{\mathsf{dq}}-{v}_{\mathsf{dq}}^{*})\big{)}\end{split}$
(4.30)

Full energy function can be shown as:

$\begin{split}\tilde{H}&=\frac{1}{2}\tilde{\omega}^{\top}M\tilde{\omega}+\frac{1}{2}\tilde{i}_{s,\mathsf{dq}}^{\top}L_{s}\tilde{i}_{s,\mathsf{dq}}+\frac{1}{2}\tilde{v}_{\mathsf{dq}}^{\top}C\tilde{v}_{\mathsf{dq}}+\frac{1}{2}\tilde{i}_{\mathsf{t,dq}}^{\top}L_{\mathsf{t}}\tilde{i}_{\mathsf{t,dq}}+S(\theta_{\mathsf{dq}})\\\
\tilde{H}&=\frac{1}{2}\tilde{\omega}^{\top}M\tilde{\omega}+\frac{1}{2}\tilde{i}_{s,\mathsf{dq}}^{\top}L_{s}\tilde{i}_{s,\mathsf{dq}}+\frac{1}{2}\tilde{v}_{\mathsf{dq}}^{\top}C\tilde{v}_{\mathsf{dq}}+\frac{1}{2}\tilde{i}_{\mathsf{t,dq}}^{\top}L_{\mathsf{t}}\tilde{i}_{\mathsf{t,dq}}+\frac{1}{2}(\hat{v}_{\mathsf{dq}}-{v}_{\mathsf{dq}}^{*})^{\top}C(\hat{v}_{\mathsf{dq}}-{v}_{\mathsf{dq}}^{*})\end{split}$
(4.31)

Expression for its Lie derivative is:

$\begin{split}\frac{d}{dt}\tilde{H}&=-\tilde{\omega}^{\top}D\tilde{\omega}+\tilde{\omega}^{\top}I_{r}^{*}(L_{m}\otimes
e_{2}^{\top})\mathbf{R}_{\theta_{\mathsf{dq}}}^{\top}\Pi_{2}^{\top}\omega_{0}\boldsymbol{j}C\tilde{v}_{\mathsf{dq}}-\tilde{\omega}^{\top}\nabla_{\theta_{\mathsf{dq}}}S-\tilde{i}_{s,\mathsf{dq}}^{\top}R_{s}\tilde{i}_{s,\mathsf{dq}}-\tilde{i}_{s,\mathsf{dq}}^{\top}\tilde{v}_{\mathsf{dq}}+\\\
&+\tilde{i}_{s,\mathsf{dq}}^{\top}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}\tilde{\omega}-\tilde{i}_{s,\mathsf{dq}}^{\top}\omega_{0}\boldsymbol{j}L_{s}\Pi_{1}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}\tilde{\omega}+\tilde{v}_{\mathsf{dq}}^{\top}\tilde{i}_{s,\mathsf{dq}}-\tilde{v}_{\mathsf{dq}}^{\top}\mathbf{E}\tilde{i}_{\mathsf{t,dq}}-\\\
&-\tilde{v}_{\mathsf{dq}}^{\top}G\tilde{v}_{\mathsf{dq}}-\tilde{v}_{\mathsf{dq}}^{\top}\omega_{0}\boldsymbol{j}C\Pi_{2}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}\tilde{\omega}-\tilde{i}_{\mathsf{t,dq}}^{\top}R_{\mathsf{t}}\tilde{i}_{\mathsf{t,dq}}+\tilde{i}_{\mathsf{t,dq}}^{\top}\mathbf{E}^{\top}\tilde{v}_{\mathsf{dq}}-\\\
&-\tilde{i}_{\mathsf{t,dq}}^{\top}\omega_{0}\boldsymbol{j}L_{\mathsf{t}}\Pi_{3}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}\tilde{\omega}+(\nabla_{\theta_{\mathsf{dq}}}S)^{\top}\tilde{\omega}=\\\
\frac{d}{dt}\tilde{H}&=\begin{bmatrix}\tilde{\omega}\\\
\tilde{i}_{s,\mathsf{dq}}\\\ \tilde{v}_{\mathsf{dq}}\\\
\tilde{i}_{\mathsf{t,dq}}\end{bmatrix}^{\top}\underbrace{\begin{bmatrix}-D&&&\\\
-\big{[}\omega_{0}\boldsymbol{j}L_{s}\Pi_{1}-I\big{]}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}&-R_{s}&&\\\
-2\omega_{0}\boldsymbol{j}C\Pi_{2}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}&&-G&\\\
-\omega_{0}\boldsymbol{j}L_{\mathsf{t}}\Pi_{3}\mathbf{R}_{\theta_{\mathsf{dq}}}(L_{m}\otimes
e_{2})I_{r}^{*}&&&-R_{\mathsf{t}}\\\
\end{bmatrix}}_{\text{-Q}}\begin{bmatrix}\tilde{\omega}\\\
\tilde{i}_{s,\mathsf{dq}}\\\ \tilde{v}_{\mathsf{dq}}\\\
\tilde{i}_{\mathsf{t,dq}}\end{bmatrix}<0.\end{split}$ (4.32)

For $\frac{d}{dt}\tilde{H}$ to be negative definite we require that
$A=\frac{1}{2}(Q+Q^{\top})$ be positive definite, from where explicit
inequality condition on damping term $D$ can be derived. This inequality
should hold for any angle combination during transients, thus additional
control effort in terms of $K_{p}$ could be used as before.

Regarding the interpretation, $S(\theta_{\mathsf{dq}})$ has similar structure
to the forth (fully decoupled) energy function from topology-independence
chapter. Indeed, $S(\theta_{\mathsf{dq}})$ is nothing more than the sum of
capacitor energies throughout the network and block-diagonal structure of
matrix C ensures decentralization. From energy perspective every local
controller cares about its own capacitor only. The problem with fully
decoupled energy function back in Chapter 3 was lack of coupling and
robustness. In this case, although the energy is fully decoupled, coupling
enters not through function structure, but in form of measurements, which
possess implicit information on neighbor states.

As it was already discussed, only diagonal dominance of $\Pi_{2}$ is left as a
precondition for this controller. The next step would be to check whether
slightly violating this condition degrades the controller performance. To
check that, the following scenario was simulated on Simulink: triangular
network (as before), full model, literature (real life) taken specifications
of power system, arbitrary initial conditions and some OPF-prescribed set-
points as a local references (Chapter 2).

First, exact controller, given by the first line in equation (4.28), was
tested. Note that ideally it is not decentralized, but almost global
asymptotic stability proof was shown above for it. Rotor angles and angular
frequency transients are depicted in Fig.4.1. As one can see both frequencies
and angles converge to the prescribed references and to prove that one also
plots the transient difference of node voltages and line currents from their
set-points in $\mathsf{dq}$-frame (Fig.4.2). Second, approximate controller,
given by the second line in equation (4.28), was tested. This controller is
decentralized, but no rigorous stability proof was shown above for it. Rotor
angles and angular frequency transients are depicted in Fig.4.3. Surprisingly,
both frequencies and angles show very similar transients to the ones observed
before. To back that up, one also plots the transient difference of node
voltages and line currents from their set-points in $\mathsf{dq}$-frame
(Fig.4.4).

Regarding the control effort of this decentralized controller, Fig.4.5-4.7
show that major contribution comes from the first term in equation (4.28),
which deals with canceling electrical torque. Thus, following the torque
demand in real time constitutes main control effort, which is not much
different from steady state demand as it turns out from transients, shown in
Fig.4.5-4.7.

(a) Angular frequency transient

(b) Angles transient

Figure 4.1: Angular frequency and angle transients with exact controller

(a) $v_{\mathsf{dq}}-v_{\mathsf{dq}}^{*}$

(b) $i_{\mathsf{t,dq}}-i_{\mathsf{t,dq}}^{*}$

Figure 4.2: Node voltage and line current transients with exact controller

(a) Angular frequency transient

(b) Angles transient

Figure 4.3: Angular frequency and angle transients with approximate controller

(a) $v_{\mathsf{dq}}-v_{\mathsf{dq}}^{*}$

(b) $i_{\mathsf{t,dq}}-i_{\mathsf{t,dq}}^{*}$

Figure 4.4: Node voltage and line current transients with approximate
controller Figure 4.5: Control action measures for generator 1 Figure 4.6:
Control action measures for generator 2 Figure 4.7: Control action measures
for generator 3

Overall, it turns out that slight violations of $\Pi_{2}$’s block-diagonal
dominance do not degrade controller performance, which is reassuring, because
specifications for simulation were taken from real world MV power system.
Those are summarized in Table 4.1. An in-depth theoretical analysis of
$\Pi_{2}$ internal structure is needed to be on safe side, nevertheless, it
can be seen from the first glance that major diagonal dominance violation
would require unrealistic parameter combination of a power system. In fact, to
preserve diagonal dominance, one requires high local load $G$, high
$R_{\mathsf{t}}$ and $L_{\mathsf{t}}$ values as well as low $L_{s},R_{s}$.
These are compatible with inverter model and voltage level, since relatively
low $L_{s},R_{s}$ in combination with moderate $C$ would result in high output
filter corner frequency and, consequently, high switching frequency condition,
which is exactly the assumption of 3-phase 2-level 6-switch inverter average-
switch model within local matching control at low voltage level (LV). Even
with MV parameters from Table 4.1, which do not comply with desired
specification description, controller is able to follow the references fairly
well as shown by simulations.

Table 4.1: Simulated system parameters Variable |  Value |  Explanation
---|---|---
$\omega_{0}$ |  $100\pi$ |  nominal electrical frequency
$M1,M2,M3$ |  $22,10,45\cdot 10^{3}$ |  Rotor moments of inertia ($kg\cdot m^{2}$)
$L_{m1},L_{m2},L_{m3}$ |  $0.04,0.08,0.02$ |  mutual inductances ($H$)
$i_{r1},i_{r2},i_{r3}$ |  $1950,975,3900$ |  reference rotor currents ($A$)
$D_{1},D_{2},D_{3}$ |  $4.0,1.5,8.5\cdot 10^{3}$ |  rotor mechanical damping
$K_{p1},K_{p2},K_{p3}$ |  $5D_{1},5D_{2},5D_{3}$ |  control gain $K_{p}$
$G_{1},G_{2},G_{3}$ |  $0.8,0.4,1.0$ |  nominal local loads ($S$)
$R_{\mathsf{t,1}},R_{\mathsf{t,2}},R_{\mathsf{t,3}}$ |  $0.165,0.166,0.07$ |  line resistances ($\Omega$)
$L_{\mathsf{t,1}},L_{\mathsf{t,2}},L_{\mathsf{t,3}}$ |  $4.7,3.8,2.4\cdot 10^{-3}$ |  line inductances ($H$)
$R_{s1},R_{s2},R_{s3}$ |  $0.166,0.07,0.5$ |  stator resistances ($\Omega$)
$L_{s1},L_{s2},L_{s3}$ |  $0.18,0.10,0.66\cdot 10^{-3}$ |  stator inductances ($H$)
$C_{1},C_{2},C_{3}$ |  $0.01,0.2,4\cdot 10^{-3}$ |  load-side capacitances ($F$)

##  Chapter 5 Conclusion

###  5.1 Summary

To sum up, this thesis presents three major results:

$\bullet$ Bijective translation between network-level OPF-prescribed set-
points and local controller references was proposed. This is a major
decentralization-enabling step, since decentralized controllers need to have a
certificate that converging to local references imposes OPF on the whole
network. Bijective translation provides this certificate under certain
conditions.

$\bullet$ Provided communication, circulating power flow problem can be solved
using the notion of communication graph. There are associated robustness
problems with real-time angle communication and clock drift, which motivate
study of decentralization.

$\bullet$ Motivated by inherent robustness problems associated with
communication, a decentralized controller was proposed for full model applied
on un-cycled or small-cycled network topologies. Interestingly, proposed
controller does not require constant R/L assumption for transmission lines,
stator and DC-side time-scale separation and strict diagonal dominance
conditions on some matrices, as it turns out from simulations with real-life-
taken network parameters.

###  5.2 Outlook

As an outlook, three interesting research directions are indicated:

$\bullet$ First, it would be interesting to investigate the potential of
proposed decentralized controller in dealing with circulating power flows in
large cycled topologies. For that, the emergence of circulating power, which
was not the focus of this thesis, should be studied more thoroughly.

$\bullet$ Second, robustness of the controller with respect to load change
and/or clock drift should be checked. Throughout the work we considered
constant load $G$. Change of parameters in this block-diagonal matrix would
affect the $\Pi$-map and, consequently, the available solution manifold in
angle-space, meaning that further improvements/changes might be necessary for
load-robustness. Regarding clock drift, although decentralized controller
realizes the gradient of a decoupled energy function, coupling is still
preserved through measurements, utilized in gradient-creation. This ”implicit”
coupling gives us a hope that the controller can manage clock drift,
nevertheless, extensive study/simulations need to be done to have a rigorous
proof.

$\bullet$ Third, mentioned controller should be extended to other types of
networks, having only load or generation at some nodes and not both, as it is
now. Various types of ”non-resistive” loads could be considered. Moreover,
inverter average-switch model within local matching control implicitly assumes
high switching frequency, meaning low or lowish-medium voltage levels from
practical perspective. There is nothing wrong with it, as micro-grid type of
systems are most likely to be deployed at residential distribution level,
which operates in stated voltage range. However, for densely-populated cities,
where there is no much place for distributed generation of renewable energy,
it would be useful to step-up voltage level and transmit the power from
renewable energy farms. For that purpose, other suitable inverter topologies,
such as MMC, could be considered. Benefit and applicability of local model-
matching, compatibility with MMC-inherent control objectives, such as voltage-
balancing etc., are naturally arising issues in that case.

## Bibliography

  * [1] C. Arghir, T. Jouini, and F. Dörfler. Grid-forming control for power converters based on matching of synchronous machines. [Online]. Available: https://arxiv.org/pdf/1706.09495.pdf
  * [2] M. Colombino, D. Groß, and F. Dörfler, “Global phase and voltage synchronization for power inverters: A decentralized consensus-inspired approach,” in _2017 IEEE 56th Annual Conference on Decision and Control (CDC)_. IEEE, Dec. 2017, pp. 5690–5695.
  * [3] A. J. Korsak, “On the question of uniqueness of stable load-flow solutions,” in _IEEE Trans. Power Appar. Syst. PAS-91_. IEEE, 1972, p. 1093.
  * [4] H. Bevrani, B. François, and T. Ise, _Microgrid Dynamics and Control_. Wiley.
  * [5] e. a. Ali TAYYEBI, “Grid-forming converters – inevitability, control strategies and challenges in future grids application,” in _CIRED Workshop - Ljubljana_ , June 2018.
  * [6] A. J. V. D. Schaft, _L2-gain and passivity techniques in nonlinear control_. Vol. 2. Springer.
  * [7] M. Maggiore, _Introduction to Nonlinear Control Systems_. Lecture Notes.
  * [8] W. A. Bukhsh, A. Grothey, K. I. M. McKinnon, and P. A. Trodden, “Local solutions of the optimal power flow problem,” in _IEEE Trans. Power Syst., vol. 28_. IEEE, 2013, p. 4780.
  * [9] e. a. B. B. Johnson, “Oscillator-based inverter control for islanded three-phase microgrids,” in _IEEE Journal of Photovoltaics, vol. 4, no. 1_. IEEE, Jan. 2014, pp. 387–395.
  * [10] C. Arghir and F. Dörfler. Energy-based stabilization of network flows in multi-machine power systems. [Online]. Available: https://arxiv.org/abs/1804.09681
  * [11] F. Bullo, J. Cortes, F. Dörfler, and S. Martinez. Lectures on network systems. [Online]. Available: http://motion.me.ucsb.edu/book-Ins/
  * [12] J. W. Chapman and et. al., “Stabilizing a multi-machine power system via decentralized feedback linearizing excitation control,” in _IEEE Transactions on Power Systems 8.3_. IEEE, 1993, pp. 830–839.
  * [13] Chiang and Hsiao-Dong, “Direct methods for stability analysis of electric power systems: theoretical foundation, bcu methodologies, and applications.” John Wiley and Sons, 2011\.
  * [14] R. Delabays, T. Coletta, and P. Jacquod, “Multistability of phase-locking and topological winding numbers in locally coupled kuramoto models on single-loop networks,” in _Journal of Mathematical Physics 57_. American Institute of Physics, March 2016, pp. 032 701–1–21.
  * [15] C. D. Persis and N. Monshizadeh, “Bregman storage functions for microgrid control,” in _IEEE Transactions on Automatic Control_. IEEE, Jan. 2018, pp. 53–68.

See pages - of SemesterProject.pdf See pages - of Originality.pdf