A power consensus algorithm for DC microgrids ?
C. De Persis a , E.R.A. Weitenberg a , F. Dörfler b

arXiv:1611.04192v2 [math.OC] 24 Feb 2017

a

Engineering and Technology Institute Groningen, University of Groningen, Nijenborgh 4, 9747 AG Groningen, The
Netherlands
b

Automatic Control Laboratory, Swiss Federal Institute of Technology, 8092 Zurich, Switzerland

Abstract
A novel power consensus algorithm for DC microgrids is proposed and analyzed. DC microgrids are networks composed of DC
sources, loads, and interconnecting lines. They are represented by differential-algebraic equations connected over an undirected
weighted graph that models the electrical circuit. A second graph represents the communication network over which the source
nodes exchange information about the instantaneous powers, which is used to adjust the injected current accordingly. This give
rise to a nonlinear consensus-like system of differential-algebraic equations that is analyzed via Lyapunov functions inspired
by the physics of the system. We establish convergence to the set of equilibria consisting of weighted consensus power vectors
as well as preservation of the weighted geometric mean of the source voltages. The results apply to networks with constant
impedance, constant current and constant power loads.
Key words: DC microgrids, Power sharing, Distributed control, Nonlinear consensus, Lyapunov stability analysis

1

Introduction

The proliferation of renewable energy sources and storage devices that are intrinsically operating using the DC
regime is stimulating interest in the design and operation
of DC microgrids, which have the additional desirable
feature of preventing the use of inefficient power conversions at different stages. These DC microgrids might
have to be deployed in areas where an AC microgrid is
already in place, creating what is called a hybrid microgrid [1], for which rigorous analytical studies are still
in their infancy. Furthermore, the envisioned future in
which power generation is far away from the major consumption sites raises the problem of how to transmit
power with low losses, a problem for which High Voltage Direct Current (HVDC) networks perform comparatively better than AC networks. Finally, also mobile
grids on ships, aircrafts, and trains are based on a DC
architecture.
With DC and hybrid microgrids, as well as HVDC net? An abridged version of the paper has been submitted to
the 20th IFAC World Congress, 9-14 July 2017, Toulouse,
France.
Email addresses: c.de.persis@rug.nl (C. De Persis),
e.r.a.weitenberg@rug.nl (E.R.A. Weitenberg),
dorfler@ethz.ch (F. Dörfler).

works, on the rise, we need to develop a deeper systemtheoretic understanding of this interesting class of dynamical networks. In this paper we propose and analyse a control algorithm for a DC microgrid that enforces
power sharing among the different power sources.
1.1

Literature review

The literature on DC microgrids is rapidly growing. We
summarize below the contributions that share a systems
and control-theoretic point of view on these networks.
The work [2] relies on a cooperative control paradigm
for DC microgrids to replace the conventional secondary
control by a voltage and a current regulator. In [3] a
voltage droop controller for DC microgrids inspired by
frequency droop in AC power networks is analyzed,
and a secondary consensus control strategy is added to
prevent voltage drift and achieve optimal current injection. The paper [4] models the DC microgrid via the
Brayton-Moser equations and uses this formalism to
show that with the addition of a decentralized integral
controller voltage regulation to a desired reference value
is achieved. Other schemes achieving desirable power
sharing properties are proposed but no formal analysis
is provided. In [5], a secondary consensus-based control scheme for current sharing and voltage balancing
in DC microgrids is designed in a Plug-and-Play fashion to allow for the addition or removal of generation

27 February 2017

Reactive power sharing algorithms have been first suggested by [17] for network-reduced AC microgrids whose
voltage dynamics show similar features as in DC grids.
In this paper we show that a similar idea can be adopted
also for network preserved DC microgrids. The novelties
of this contribution with respect to [17] are the different
dynamics of the system under study, the explicit consideration of algebraic equations in the model and the use
of Lyapunov arguments to prove the main results.

units. A distributed control method to enforce power
sharing among a cluster of DC microgrids is proposed
in [6]. Other work has focused on the challenges in the
stability analysis of DC microgrids using consensus-like
algorithms due to the interaction between the communication network and the physical one [7]. Finally,
feasibility of the nonlinear algebraic equations in DC
power circuits is studied by [8], [9], and [10].
A closely related research area is that of multi-terminal
HVDC transmission systems. The paper [11] focuses
on cooperative frequency control for these networks. In
[12] distributed controllers that keep the voltages close
to a nominal value and guarantee a fair power sharing
are considered, whereas passivity-based decentralized PI
control for the global asymptotic stabilisation of multiterminal high-voltage is studied in [13]. The paper [14]
studies feasibility and power sharing under decentralized
droop control. We refer to [15, Chapter 4] for an annotated bibliography of HVDC transmission systems.
1.2

1.3

Paper organization

The model of the DC microgrid is introduced in Section
2. The power consensus algorithm is introduced in Section 3. The analysis of the closed-loop system is carried
out in Section 4 for the general case of ZIP loads, and
then specialised to the case of ZI loads, since the latter
permits to obtain stronger results under weaker conditions. The simulations of the algorithm are provided in
Section 5. Conclusions are drawn in Section 6.

Main contribution

This paper focuses on a new control algorithm to stabilize a DC microgrid under different load characteristics while achieving power sharing among the sources.
Our controller is enabled by communicating the instantaneous source power measurements among neighboring
source nodes, averaging these measurements and setting
the voltage at the source terminals accordingly. An additional feature of the algorithm is that a weigthed geometric average of the source voltages is preserved.

1.4

Notation

Given a vector v, the symbol [v] represents the diagonal
matrix whose diagonal entries are the components of v.
The notation col(v1 , v2 , . . . , vn ), with vi scalars, represents the vector [v1 v2 . . . vn ]T . If vi are matrices having
the same number of columns, then col(v1 , v2 , . . . , vn ) denotes the matrix [v1T v2T . . . vnT ]T . The symbol 1n represents the n-dimensional vector of all 1’s, whereas 0m×n
is the m × n matrix of all zeros. When the size of the matrix is clear from the context the index is omitted. The
n × n identity matrix is represented as In . Given a vector v ∈ Rn , the symbol ln(v) denotes the element-wise
logarithm, i.e., the vector [ ln(v1 ) . . . ln(vn ) ]T .

The system dynamics present interesting features. By
averaging the power measurements that the sources communicate amongst each other, the system dynamics becomes an intriguing combination of the physical network (the weighted Laplacian of the electrical circuit
appearing in the power measurements) and the communication network (over which the information about
the power measurements is exchanged). “ZIP” (constant
impedance, constant current and constant power) loads
introduce algebraic equations in the system’s dynamics,
adding additional complexity and nonlinearities.

2

DC resistive microgrid

The DC microgrid is modeled as an undirected connected graph G = (V, E), with V := {1, 2, . . . , n} the set
of nodes (or buses) and E ⊆ V × V the set of edges. The
edges represent the interconnecting lines of the microgrid, which we assume here to be resistive. Associated
to each edge is a weight modeling the conductance (or
reciprocal resistance) 1/rk > 0, with k ∈ E. The set
of nodes is partitioned into the two subsets of ns DC
sources Vs and nl loads Vl , with ns + nl = n.

To analyze this system of nonlinear differential-algebraic
equations without going through a linearization of the
dynamics, Lyapunov-based arguments become very convenient. The Lyapunov functions in this case are constructed starting from the power dissipated in the network that is further shaped to take into account the
specifics of the dynamics. The presence of the loads,
which shift the equilibrium of interest, is taken into account by the so-called Bregman function [16]. The level
sets of the Lyapunov functions are used to estimate the
excursion of the state response of these systems and
therefore, combined with the preservation of the geometric average of the source voltages, can be used to obtain
an estimate of the voltage at steady state.

The current-potential relation in a resistive network is given by the identity I = BΓB T V , with
B ∈ Rn×|E| being the incidence matrix of G and
−1
Γ = diag{r1−1 , . . . , rm
} the diagonal matrix of conductances. Considering the partition of the nodes in sources

2

represents an ideal current source that is controlled as a
function of the local voltage Vi and the injected power
Pj = Vj Ij at the neighboring node sources j ∈ Nc,i .

and loads, the relation rewrites as
" #
Is
Il

"
=

Bs ΓBsT Bs ΓBlT

#" #
Vs

Bl ΓBsT Bl ΓBlT

"
=:

Vl

Yss Ysl

#" #
Vs

Yls Yll

,

Vl

The dynamic controllers (3)–(5) are initialised at positive values of the voltage, that is Vi (0) > 0 for all i ∈ Vs .
It will be made evident in later sections that these cons
trollers render the positive orthant Rn>0
positively invariant, thus showing that the positivity of the initial source
voltages yields positivity of these variables for all t ≥ 0.

(1)
where Is = col(I1 , . . . , Ins ), Il = col(Ins +1 , . . . , In ),
Vs = col(V1 , . . . , Vns ), Vl = col(Vns +1 , . . . , Vn ) and
B = col(Bs , Bl ).
Observe that both Yss and Yll are positive definite since
they are principal submatrices of a Laplacian of a connected undirected graph. This allows us to eliminate the
load voltages as Vl = Yll−1 Il − Yll−1 Yls Vs and reduce the
network to the source nodes Vs with balance equations
Is − Yll−1 Yls Il = Yred Vs ,

Remark 1 (Circuit interpretation) The control algorithm has the circuit interpretation given in Fig. 1.
Comparing with [4, (4)], the ideal current source ui can
Ii

(2)
ui

where Yred = Yss − Yll−1 Yls Yls is known as the Kronreduced conductance matrix [18] and −Yll−1 Yls Il is the
mapping of the load current injections to the sources.

Ci (Vi )

Fig. 1. A circuit interpretation of the controller (3).

3

be generated also by a voltage source with value vi in series with a resistance ri provided that vi = ri ui + Vi .
Finally, the dynamic droop controller in [3] corresponds
in our notation to a constant capacitance Ci and current
source ui .

Power consensus controllers

We propose controllers that force the different sources to
share the total power injection in prescribed ratios [17].
For this purpose, a communication network is deployed
to connect the source nodes, through which the controllers exchange information about the instantaneous
injected powers. This communication network is modelled as an undirected unweighted graph (Vc , Ec ), where
Vc = Vs . Associated with the communication graph is
the ns × ns Laplacian matrix Lc = Dc − Ac , where Dc is
the degree matrix and Ac is the adjacency matrix of the
communication graph. Note that the nodes of the communication network (but not necessarily the edges) coincide with the source nodes of the microgrid. For each
node i ∈ Vs , the set Nc,i = {j ∈ Vs : {i, j} ∈ Ec } represents the neighbors connected to node i via the communication graph.

Remark 2 (Digital implementation)The control algorithm above need not be implemented analogically. In
its digital implementation, information regarding the injected powers of neighbouring sources, Pj , j ∈ Nc,i , are
broadcasted and processed along with the current measurement Ii to compute the source voltage value Vi applied at the source terminals.
Multiplying both sides of (3) by Vi2 Dc,i Ci−1 , one arrives
at the closed-loop system
Ci V̇i = −Vi Dci Ci−1 Pi + Vi

−1 2
Ci (Vi ) = Vi−2 Dci
Ci , i ∈ Vs
(4)
can be interpreted as a nonlinear capacitance, Ci > 0 is
a positive parameter of suitable units such that Ci (Vi )
actually has the units of a capacitance, Ii is the injected
current at node i ∈ Vs as defined in (1), and the term

X

Cj−1 Pj ,

i ∈ Vs

(Cj−1 Pj

− Ci−1 Pi ),

i ∈ Vs ,

(6)

that is, the voltage at the source terminal is updated according to a weighted power consensus algorithm scaled
by the voltage. Provided that Vi 6= 0 (a property that will
be established in the next sections), equation (6) shows
that at steady state the algorithm achieves proportional
power sharing according to the Ci ratios, namely

(3)

where

−1
ui = Vi−1 Dci
Ci

X
j∈Nc,i

Controllers. The proposed controllers are of the form
i ∈ Vs ,

Cj−1 Pj

j∈Nc,i

= Vi

Ci (Vi )V̇i = −Ii + ui ,

X

Pj
Pi
=
,
Cj
Ci

∀i, j ∈ Vs .

(7)

A detailed characterisation of the steady-state power signals is given in the next section (Lemma 1).

(5)

j∈Nc,i

3

col(P1 , . . . , Pns ) given by

For interpretation purposes, we write (6) as
X
d
(Cj−1 Pj − Ci−1 Pi ),
Ci ln(Vi ) =
dt

Ps = [Vs ]Is = [Vs ](Yss Vs + Ysl Vl )

i ∈ Vs .

j∈Nc,i

are source power injections and

In a classic power system analysis [19], the term Ci ln(Vi )
is the natural energy representation of a power source of
constant value Ci . The interpretation of the closed loop
(6) is then that the voltage at this constant power source
is adapted according to a power consensus algorithm.

Il (Vl ) = Il∗ − Yl∗ Vl + [Vl ]−1 Pl∗

Remark 4 (Nonlinear consensus algorithms) To
compare the algorithm (6) with related nonlinear consensus algorithms proposed in the literature, we neglect
the algebraic constraints and the differentiation between
sources and loads. This allows us to rewrite (6) as

Ci V̇i = −Ii + pi
(Cj−1 Vj pj − Ci−1 Vi pi ), i ∈ Vs
(8)
where pi is a control variable in units of currents. Any
steady state of this controller would guarantee for all i ∈
Vs that V̇i = 0, and pi = Ii is the steady-state current
injection, and the vector of power injections Cs−1 [Vs ]p has
all identical entries (power sharing). Numerical results
(see Section 5) show that (6) and (8) perform similarly.
Indeed, in the limit Di = 0, near steady-state, and for
nearly unit voltages (in per unit system), the closed-loops
(8) and (6) have similar dynamics. In the rest of the paper
we focus on the analysis of (6).
P

(11)

are the load currents. The interconnected closed-loop
DC microgrid is then entirely described by equations
(9), (10), (11). An example of a simple closed-loop DC
microgrid with two sources and one constant impedance
load is given in Figure 2.

Remark 3 (Alternative power sharing control) A
possibly more simplistic and obvious power sharing controller inspired by the current-sharing controller in [3] is
based on a distributed averaging integral control given by

Di ṗi = Ii − pi +

(10)

j∈Nc,i

C V̇ = −[V ]Lc C −1 [V ]BΓB T V.
The weighted power mean consensus algorithms of
[20,21], on the other hand, can be written as [W ]V̇ =
[V ]1−r BΓB T V , where W is vector of weights satisfying
1T W = 0 and r ∈ R. In the special case r = 0, we get
[W ]V̇ = [V ]BΓB T V,
which is known to converge to the consensus value
wn
V1w1 . . . VP
n . The analysis
Qn is based on the Lyapunov
n
function i=1 wi Vi − i=1 Viwi .
The nonlinear power consensus algorithm is different in
that it uses another layer of averaging in addition to the
averaging induced by the physical network. This, and the
algebraic constraints, requires a different analysis based
on physically inspired Lyapunov functions.

Loads. Depending on the particular load models, the
term Il in (1) takes different expression and will henceforth be denoted as Il (Vl ) to stress the functional dependence on the load voltages. Prototypical load models
that are of interest include the following:
l
(i) constant current loads: Il (Vl ) = Il∗ ∈ Rn<0
,
∗
(ii) constant impedance: Il (Vl ) = −Yl Vl , with Yl∗ > 0
a diagonal matrix of load conductances, and Vl =
col(Vns +1 , . . . , Vns +nl ), and
l
.
(iii) constant power: Il (Vl ) = [Vl ]−1 Pl∗ , with Pl∗ ∈ Rn<0

4

Power consensus algorithm with ZIP loads

In this section we analyze the closed-loop system (9),
(10), (11). We start by studying its equilibria, namely
the set of points V ∈ Rn>0 that satisfy (10), (11), and

To refer to the three load cases above, we will use the indices “I”, “Z” and “P” respectively. The analysis of this
paper will focus on the more general case of a parallel
combination of the three loads, thus on the case of “ZIP”
loads. Moreover, additional and stronger statements results on the “ZI” case will be reported.

"

0

#

−Il (Vl )

"
=−

[Vs ]Lc Cs−1 Ps
Bl ΓB T V .

#
(12)

Bearing in mind (1), (6), and vectorizing the expressions
to avoid cluttered formulas, the closed-loop system is
"

Cs V̇s

#

−Il (Vl )

"
=−

[Vs ]Lc Cs−1 Ps
Bl ΓB T V

4.1

#
,

Steady-state characterization

(9)
In the following, we show that the equilibria are fully
characterized by power balance equations at the sources
and current balance equations at the loads, respectively.

where V = col(Vs , Vl ), Cs = diag(C1 , . . . , Cns ), Ps =

4

r1

u1

r2

I1

−I2

r3

C1 (V1 )

C2 (V2 )

u2

I3
Fig. 2. Circuit considered in Example 1.

Next, we left-multiply (14) by [ 1Tns 1Tnl ] to obtain

Lemma 1 (System equilibria) The equilibria of the
system (9), (10), (11) are equivalently characterized by

1Tns [Vs ]−1 Cs 1ns p∗s + 1Tnl Il (Vl ) = 0.

EZIP = {V ∈ Rn>0 : IZIP (V ) = 0 , PZIP (V ) = 0},

The latter equation can be solved for p∗s as in (13). From
Il (Vl ) = Bl ΓB T V , we obtain (see (2)) IZIP (V ) = 0 or

where IZIP (V ) = 0 is the current balance at the loads

Vl = −Yll−1 Yls Vs + Yll−1 Il (Vl ),

IZIP (V ) = Il (Vl ) − Yll Vl − Yls Vs ,

which replaced in the first equation of (14) returns

PZIP (V ) = 0 depicts the power balance at the sources

Yss Vs + Ysl (−Yll−1 Yls Vs + Yll−1 Il (Vl )) = [Vs ]−1 Cs 1ns p∗s .

PZIP (V ) = [Vs ]Yred Vs+ [Vs ]Ysl Yll−1 Il (Vl ) − Ps ,
|{z}
| {z } |
{z
}
network
dissipation

By rearranging the terms, we arrive at

source
injections

load
demands

(15)

Yred Vs + Ysl Yll−1 Il (Vl ) − [Vs ]−1 Cs 1ns p∗s = 0,
Yred is the Kron-reduced conductance matrix, Yll−1 Ysl Il (Vl )
is the mapping of the ZIP loads Il (Vl ) to the source buses
in the Kron-reduced network as in (2), and Ps is vector of
power injections by the sources written for V ∈ EZIP as
Ps = −Cs 1

1T Il (Vl )
=: Cs 1p∗s .
1T [Vs ]−1 Cs 1

which can be reformulated as PZIP (V ) = 0 after leftmultiplying by [Vs ] and bearing in mind (13). The latter
and (15) show that V ∈ EZIP .
Conversely, let V ∈ EZIP . Then the equation Il (Vl ) =
Bl ΓB T V in (12) is trivially satisfied. From PZIP (V ) =
0, and Il (Vl ) = Bl ΓB T V written as (15), and going backwards through the passages above, we arrive at

(13)

Yss Vs + Ysl Vl = [Vs ]−1 Cs 1ns p∗s ,
Observe that the steady-state injections (13) achieve indeed power sharing, and the asymptotic power value p∗s
to which the source power injections converge (in a proportional fashion according to the coefficients Ci , i ∈ Vs )
is the total current demand divided by the weighted sum
of the steady-state source voltages. The latter values and
those of the load voltages are interestingly entangled by
the power balance at the sources PZIP (V ) = 0 and the
current balance equations at the loads IZIP (V ) = 0.

or equivalently at [Vs ]Bs ΓB T V = Cs 1ns p∗s . Hence, the
power vector Ps = [Vs ]Bs ΓB T V satisfies Lc Cs−1 Ps = 0,
that is, the first equation in (12). Hence, V ∈ EZIP
implies that the equilibrium equations (12) are met. 
We make the standing assumption that equilibria exist:
Assumption 2 EZIP 6= ∅.
Remark 5 (Existence of the equilibria EZIP ) The
analytical investigation of the existence of the equilibria EZIP is deferred to a future research. This is a topic
of interest on its own and similar problems have been
dealt with in recent work about the solvability of reactive power flow equations [22,8,9,23]. For instance, the
problem in [23] boils down to the solution of quadratic
algebraic equations of the form [Vl ]Yll Vl − [Vl ]Yll Vl∗ +
Ql = 0, where Ql is the vector of constant power load
demands and Vl∗ is the so called vector of open cir-

Proof. Let V be an equilibrium of (9), (10), (11), that
is let V ∈ Rn>0 satisfy (12). From the first equation, 0 =
[Vs ]Lc Cs−1 Ps , it immediately follows that Ps = Cs 1ns p∗s
for some scalar p∗s . We rewrite the current balances as
"
#
[Vs ]−1 Cs 1ns p∗s
Il (Vl )

=

"
#
Bs ΓB T V
Bl ΓB T V

.

(14)

5

tions PZIP (V ) = 0, are in this case

cuit voltages (again constant). Although similarities between these equations and the equations PZIP (Vs ) =
0= [Vs ]Yred Vs + [Vs ]Ysl Yll−1 Il (Vl ) + Ps could be useful to
investigate the nature of the set EZIP , the non-quadratic
nature of PZIP (Vs ) = 0, as well as the presence of the additional equations Yll−1 Il (Vl ) − Vl = Yll−1 Yls Vs pose additional challenges. Extra insights could come from the
convex relaxation of the DC power flow equations in the
context of optimal DC power flow dispatch [10].

γ1
γ1 γ2
V1 (V1 − V2 ) −
V1 Il (Vl )+
γ1 + γ2
γ1 + γ2
Il (Vl )

γ2
γ1 γ2
V1 (V2 − V1 ) −
V2 Il (Vl )+
γ1 + γ2
γ1 + γ2
Il (Vl )

V14 − r2 Il (Vl )V13 + cr1 Il (Vl )V1 − c2 = 0
V24 − r1 Il (Vl )V23 + cr2 Il (Vl )V2 − c2 = 0,

where ri = γi−1 , i = 1, 2 (the resistance of the transmission line i connecting the source i to the load).
We want to study the solutions of these equations as functions of Il (Vl ). Then these can be regarded as two independent quartic functions for which an analytic, although
involved, expressions of the solutions exist according to
the Ferrari-Cardano’s formula. These expressions simplify if one takes r1 = r2 . Then
√there is a unique positive
solution given by V1 = V2 = c, independent of Il (Vl ).
The value of Vl is obtained from the algebraic equation
0 = Vl Bl ΓB T V −Il (Vl ), solving

1T [Vl ]−1 Pl∗
+ 1T Pl∗ ≥ 0 .
1T [Vs ]−1 Cs 1

This inequality can be reformulated as
X ai
X bi
≥
,
Vi
Vi

i∈Vl

(18)

(16)

implying that the amount of supplied power has to make
up for load demands and resistive losses. In the special case of constant power loads, Il (Vl ) = [Vl ]−1 Pl∗ , we
obtain the total (or average) power inequality 1T Ps +
1T Pl∗ ≥ 0. Equivalently, after using (13), we arrive at
−1T Cs 1

V1 V2
= 0.
V1 + V2

We study solutions to the algebraic equations on the curve
V1 V2 =: c. The reason for this choice will become clear in
Subsection 4.3. On such a curve, the equations simplify
as

Remark 6 (Equilibrium power balance and voltage inequalities) To gain further insights into the equilibrium set EZIP , recall that the vector of power injections
is P = col(P1 , Pl ) = [V ]BΓB T V , where Pl = [Vl ]Il (Vl ).
Thus, we have the inherent power balance
1T Ps + 1T Pl = V T BΓB T V ≥ 0

V1 V2
=0
V1 + V2

0 = Vl (−γ1 V1 − γ2 V2 + (γ1 + γ2 )Vl ) − Il (Vl )
√
= 2γVl2 − 2γ cVl − Il∗ + Yl∗ Vl − [Vl ]−1 Pl∗
√
= 2γVl3 + (Yl∗ − 2γ c)Vl2 − Il∗ Vl − Pl∗ .

(17)

i∈Vs

P
P
∗
∗
with ai = Pl,i
/ i∈Vl Pl,i
and bi = Ci / i∈Vs Ci , which
relates a convex combination of the reciprocals of the voltages at the loads, with a convex combination of the reciprocals of the voltages at the sources, and represents another relation between Vs , Vl in addition to those in (16).
The average voltage inequality (17) implies that the reciprocal of the harmonic average source voltage must be
larger than the reciprocal of the harmonic average load
voltage so that power can flow from sources to loads.

(19)

In the absence of loads, we have three real√roots: a double
root at Vl = 0 and a single root at Vl = c = V1 = V2 .
Since the roots of a polynomial are continuous in the
parameters, at most the double root at Vl = 0 can turn to
a complex-conjugate
root for small loading, and the real
√
root near c persists.

4.2

A Lyapunov function and hidden gradient form

In some special cases reviewed in the examples below, an
explicit characterization of the equilibria can be given.
We pursue a Lyapunov-based analysis of the stability of
the closed-loop system (9), (10), (11). Inspired by the
Lyapunov analysis of the reactive power consensus algorithm in [16], we consider the total power dissipated
through the network resistors, 12 V T BΓB T V , as the first
natural Lyapunov candidate for our analysis, to which
we add the power dissipated through the impedance

Example 1 Consider the case of two sources (ns = 2)
and one load (nl = 1) as in Figure 2, in which the constant impedance load is replaced by a ZIP load. The equa-

6

Hence, the Bregman function (22) satisfies

loads, to obtain the power losses at passive devices as
1
J(V ) = V T
2

"
#!
0 0

T

BΓB +

0 Yl∗

V.

"
#
Bs ΓB T (V − V )
∂M
∂M
∂M
=
+
=
−
∂V
∂V
∂V V =V
Bl ΓB T (V − V )
# "
#
"
([Vs ]−1 − [V s ]−1 )P s
0
−
.
([Vl ]−1 − [V l ]−1 )Pl∗
Yl∗ (Vl − V l )

(20)

Let V ∈ EZIP , and define P s = [V s ]Bs ΓB T V the source
power injection corresponding to the equilibrium source
voltage V (see (13)). To cope with the asymmetry in the
dynamics of the sources and loads we add to J the terms

Bearing in mind the equilibrium condition at the loads
Bl ΓB T V = Il (V l ) =Il∗ − Yl∗ V l + [V l ]−1 Pl∗ ,

T

H(V ) = −P s ln(Vs ),
and

and replacing it in the second line of the identity above
describing ∂M/∂V , we obtain

−Pl∗ T ln(Vl ),

K(V ) =
which is the way classical power systems transient stability analysis absorbs constant power injections [19] into
a so-called energy function defined here as

∂M
= Bl ΓB T V + Yl∗ Vl − [Vl ]−1 Pl∗ − Il∗
∂Vl
= Bl ΓB T V − Il (Vl ),

M (V ) := J(V ) + H(V ) + K(V )
#
"
0 0
T
1 T
T
)V − P s ln(Vs ) − Pl∗ T ln(Vl ).
= 2 V (BΓB +
0 Yl∗
(21)
The natural “energy function” (21) has its minimum at
the trivial zero voltage level. To center the function M
with respect to a non-trivial equilibrium V ∈ EZIP , we
use the following Bregman function [16]
M(V ) = M (V ) − M (V ) −

∂M
∂V

which equals precisely the second equation in (23).
Analogously, for the first line ∂M/∂Vs , we write
∂M
= Bs ΓB T (V − V ) − ([Vs ]−1 − [V s ]−1 )P s
∂Vs
= [Vs ]−1 Ps − [V s ]−1 P s − ([Vs ]−1 − [V s ]−1 )P s
= [Vs ]−1 (Ps − P s ),
(25)
where to write the second equality we have used the
identities Ps = [Vs ]Bs ΓB T V and P s = [V s ]Bs ΓB T V .

T

(V − V ). (22)
V =V

The next result shows a (perhaps surprising) gradient
relation between the dynamics of system (9), (10), (11)
and the Bregman function (22) above:

Now note that
Ps = [Vs ]

Lemma 3 (Gradient dynamics) The following holds
"

Lc Cs−1 Ps
T

#

Bl ΓB V − Il (Vl )

"
=

Lc [Vs ]Cs−1 0
0

I nl

#

and, multiplying both sides by Lc Cs−1 , we obtain

∂M(V )
∂V

∂M
+ Lc Cs−1 P s
∂Vs
∂M
= Lc Cs−1 [Vs ]
,
∂Vs

Lc Cs−1 Ps = Lc Cs−1 [Vs ]

(23)
for all V ∈ Rn>0 . Hence the system (9), (10), (11) can be
rewritten as a weighted gradient flow
"

Cs V̇s
0

#

"
= −

[Vs ]Lc [Vs ]Cs−1 0
0

I nl

#

∂M(V )
.
∂V

∂M
+ Ps
∂Vs

having exploited that V ∈ EZIP implies P s = Cs 1p∗s .
∂M
The identity Lc Cs−1 Ps = Lc Cs−1 [Vs ]
is the first
∂Vs
equation in (23).

(24)

In view of the dynamics (9), (10), (11), one immediately
realizes that
"
# "
#
Lc Cs−1 Ps
−[Vs ]−1 Cs V̇s
=
,
Bl ΓB T V − Il (Vl )
0

Proof. The gradient of the function M (V ) writes as
"
# "
# "
#
0
[Vs ]−1 P s
0
∂M
T
= BΓB V +
−
−
.
∂V
Yl∗ Vl
0
[Vl ]−1 Pl∗

7

showing the identity (24) which concludes the proof.
4.3

contained in the positive orthant. Without loss of generality this compact sublevel set can be taken so that all
the solutions to (9), (10), (11) that originate here locally
exist.



Convergence of solutions

The particular form of the dynamics (9), (10), (11) elucidated in Lemma 3 permits a straightforward analysis
of the convergence properties of the solutions.

The algebraic equations (11) written as in Lemma 1 are
0 =IZIP (V ) = Il (Vl ) − Yll Vl − Yls Vs
=Il∗ −Yl∗ Vl + [Vl ]−1 Pl∗ − Yls Vs − Yll Vl .

Theorem 4 (Main result) Assume that there exists
V ∈ EZIP such that

To study local solvability of these equations, we analyze
Yll + Yl∗ + [V l ]−2 [Pl∗ ] − Yls (Yss + [V s ]−2 [P s ])−1 Ysl > 0
(26)
Then there exists a compact sublevel set ΛZIP of M contained in Rn>0 such that any solution to (9), (10), (11)
that originates from initial conditions V (0) belonging to
ΛZIP exists, always remain in ΛZIP with strictly positive voltages for all times, and asymptotically converges
to the set EZIP ∩ ΛZIP ∩ VZIP , where VZIP specifies the
preserved weighted geometric mean of the source voltages


∂IZIP
= − Yll + Yl∗ + [Vl ]−2 [Pl∗ ] .
∂Vl
In view of (26), nonsingularity of ∂IZIP /∂Vl and therefore regularity of the algebraic condition holds in a neighborhood of V ∈ ΛZIP from the implicit function theorem [24]. The sublevel set ΛZIP can be taken sufficiently
small such that it is contained in the neighborhood of
regularity for the algebraic equations, thus showing the
claim that solutions starting from ΛZIP locally exist in
time, see [25, Theorem 1] and [26, Lemma 2.3].

VZIP := {(Vs , Vl ) ∈ ΛZIP : IZIP (V ) = 0 ,
C

C

V1C1 · . . . · Vnsns = V1C1 (0) · . . . · Vnsns (0)}.

When computed along these solutions, M(V (t)) satisfies
Remark 7 (Interpretation of the main condition)
The main condition (26) guarantees regularity of the algebraic equations and stability of the solutions. Its role is
revealed when converting the constant power loads and the
asymptotically constant power injections at the sources
to the equivalent impedances [V l ]−2 [Pl∗ ] and [V s ]−2 [P s ].
In this case, the equivalent conductance matrix in the
steady-state current-balance equations (1) read as
"
Yeq =

Yss
Yls

Ṁ(V (t)) =

T

V̇s (t) +
V =V (t)

∂M
∂Vl

T

V̇l (t).
V =V (t)

Notice that, by the algebraic constraint (23),
∂M
∂Vl

#
# "
Ysl
[Vs ]−2 [P s ]
0
. (27)
+
0
[Vl ]−2 [Pl∗ ] + Yl∗
Yll

= Bl ΓB T V (t) − Il (Vl (t)) = 0
V =V (t)

for all t for which a solution exists. Hence, we arrive at

By a Schur complement argument, observe that Yeq is a
well-defined (i.e., positive definite) conductance matrix
if and only if the main condition (26) holds.

Ṁ(V (t)) =
=−

Proof. Existence and boundedness of solutions. Observe
first that
#
"
# "
0
0 0
[Vs ]−2 [P s ]
∂2M
T
= BΓB +
+
,
∂V 2
0 Y∗
0
[Vl ]−2 [P ∗ ]
l

∂M
∂Vs

∂M
∂Vs

∂M
∂Vs

T

V̇s (t)
V =V (t)

T

Cs−1 [Vs ]Lc [Vs ]Cs−1
V =V (t)

∂M
∂Vs

≤ 0,
V =V (t)

where the second equality holds because of (24). The
inequality above shows that M(V (t)) is a non-increasing
function of time. By the compactness of the sublevel set
around V , the solutions are bounded, exist and belong
to ΛZIP for all times. Thus, among others the voltages
stay positive for all times.

l

(28)
Let V > 0 be an equilibrium of the system, i.e., V ∈
l
, and V > 0, the steady-state
EZIP . Since Il∗ , Pl∗ ∈ Rn<0
s
power injection at the sources satisfies P s ∈ Rn>0
by
−2
(13). Hence, [Vs ] [P s ] > 0 is positive definite. Then
the Bregman function M has an isolated minimum at
the equilibrium V , in view of (26), (28) and a standard
Schur complement argument. Then there exists a compact sublevel set ΛZIP of M around the equilibrium V

Convergence. Exploiting the regularity of the algebraic
equation, the DAE system can be reduced to an ODE
system and then the standard LaSalle invariance principle for ODE can be used to infer convergence, see also
[26]. We argue as follows. Any solution (Vs , Vl ) to the
DAE system (9), (10), (11) originating in ΛZIP is such

8

gence set EZIP ∩ ΛZIP ∩ VZIP , V1 V2 = V1 (0)V2 (0) =: c
for all t ≥ 0. Hence, as discussed in Example 1, the
expression of the (real and positive) solution to the equations (18) takes on
p a particularly simple form, namely
√
V1 = V2 = c = V1 (0)V2 (0). It follows that any point
on EZIP ∩ VZIP is such that each source voltage is the geometric mean of the initial voltage sources. Accordingly,
the load voltage must satisfy (19).

that its component Vs is a solution to the system of ODE
V̇s = −Cs−1 [Vs ]Lc [Vs ]Cs−1 (Yss Vs + Ysl δ(Vs )),

(29)

where the map Vl = δ(Vs ) denotes the solution of the
algebraic equation IZIP (V ) = 0 in ΛZIP . Define
N (Vs ) := M(Vs , δ(Vs ))

(30)

and observe that
Ṅ (Vs (t)) =

T

∂M
∂Vs

V̇s (t) +
Vs = Vs (t)
Vl = δ(Vs (t))

∂δ
∂Vs
T

∂M
=
∂Vs

∂M
∂Vl

Remark 8 (Capacitors at the loads) If loads are interconnected to the network via capacitors, the load equations are modified as

T

·
Vs = Vs (t)
Vl = δ(Vs (t))

Cl V̇l = −Il (Vl ) + Bl ΓB T V.

V̇s (t)
Notice that the equilibria of the system remain the same.
Bearing in mind (23), the load dynamics read as

Vs =Vs (t)

V̇s (t),
Vs = Vs (t)
Vl = δ(Vs (t))

Cl V̇l = −

since
∂M
∂Vl

∂M
.
∂Vl

It follows that
= Yls Vs (t) + Yll δ(Vs (t)) − Il (δ(Vs (t)))
Vs = Vs (t)
Vl = δ(Vs (t))

Ṁ = −
= Yls Vs (t) + Yll Vl (t) − Il (Vl (t)) = 0

and one can infer convergence to the set EZIP ∩ ΛZIP ∩
VZIP similarly as for the differential-algebraic model.

where the second equality holds because Vl (t) = δ(Vs (t))
on ΛZIP and the third equality because of the algebraic
equation in (9), (10), (11). It then follows that

Remark 9 (Constant voltage buses) Similarly as in
[3, Remark 3.3], one can consider also consider voltagecontrolled buses. For example, consider the scenario of
all load buses having constant (not necessarily identical)
voltages V l (see [3] for a discussion on this load condition). More precisely, a controller adjust the current injection Il depending on Vs to maintain the value of the
voltage at the constant level V l so that system (9) reads
as

Ṅ (Vs ) = (Ps − P s )T [Vs ]−1 V̇s
= −(Ps − P s )T Cs−1 Lc Cs−1 Ps
=

−PsT Cs−1 Lc Cs−1 Ps

∂M T −1
∂M ∂M T −1 ∂M
Cs [Vs ]Lc [Vs ]Cs−1
Cl
−
,
∂Vs
∂Vs
∂Vl
∂Vl

(31)

≤ 0,

where the first equality descends from (25), the second
from (9), and the third from (13).

Cs V̇s = −[Vs ]Lc Cs−1 [Vs ](Yss Vs + Ysl V l )

(32a)

Since Vs is bounded, then the standard La Salle invariance principle for ODEs yields convergence of Vs to the
largest invariant set where Lc Cs−1 Ps = 0. Moreover,
since the solutions evolve in ΛZIP , since they satisfy the
algebraic equations, and since Lc Cs−1 Ps = 0, we have
from Lemma 1 that at steady state (Vs , Vl ) ∈ EZIP .
Since (Vs , Vl ) is a solution to (9), (10), (11) that remains
in ΛZIP , convergence to the set EZIP ∩ ΛZIP is inferred.
Moreover, the quantity V1 · . . . · Vns is conserved, namely
V1 (t) · . . . · Vns (t) = V1 (0) · . . . · Vns (0) for all t. In fact, by
d
(29), Cs dt
lnVs = −Lc [Vs ]Cs−1 (Yss Vs + Ysl δ(Vs )), and
d T
therefore dt 1 Cs lnVs = 0. The thesis then follows. 

The only relevant equations for stability of (32) are the
ordinary differential equations (32a) driven by the constant term V l . We study their stability using a similar Lyapunov argument as before. Since V l is now constant, we consider a simplified version of the function
M, namely M̃(Vs ) = J˜(Vs ) + H(Vs ), where

Example 2 Consider again the case of two sources
(ns = 2) and one load (nl = 1) connected in a “T” configuration, as in Example 1. If Cs = I2 C, for some positive
real number C, the result above shows that on the conver-

∂J
are the (shifted) network losses so that ∂V
= Yss (Vs −
s
V s ) = (Yss Vs + Ysl V l ) − (Yss V s + Ysl V l ) = [Vs ]−1 Ps −
T
[V s ]−1 P s . Together with H(Vs ) = −P s ln(Vs ) +

−Il = Yls Vs + Yll V l .

(32b)

1
J˜(Vs ) = (Vs − V s )T Yss (Vs − V s )
2
˜

9

T

T

∂ M̃
∂Vs =
M̃
−1
.
−[Vs ]Lc Cs [Vs ] ∂∂V
s

P s ln(V s ) + P s [V s ]−1 (Vs − V s ), we obtain that
−1

regularity of the algebraic equation IZI (V ) = 0. Additionally, the limit set in case of ZI loads is EZI ∩ΛZI ∩VZI ,
where the set of equilibria EZI is characterized in Lemma
5, ΛZI is a sublevel set associated with the Lyapunov
function M with Pl∗ = 0, and the set VZI is defined as

[Vs ] (Ps − P s ) and thus Cs V̇s =
The convergence analysis of the solutions of the system
(32) is now analogous to the proof of Theorem 4.
4.4

The case of ZI loads

VZI := {(Vs , Vl ) ∈ ΛZI :
C

C

V1C1 · . . . · Vnsns = V1C1 (0) · . . . · Vnsns (0),

In the case of ZI loads the previous results can be
strengthened. First, the set of equilibria can be characterized by two systems of equations, one depending on
the source voltages only and the other one allowing for a
straightforward calculation of the load voltages once the
source voltages are determined. Second, the convergence
result can be established without any extra condition
on the equivalent conductance matrix in (27). Finally,
the convergence is to a point rather than to a set.

Vl = (Yll + Yl∗ )−1 (Il∗ − Yls Vs )}.
Finally, a stronger convergence result can be established,
namely any trajectory converges to a point depending on
the initial condition. This can be formalized as follows:

Theorem 6 (Point convergence) The solutions to
(9), (10), (11) with Pl∗ = 0 which originate from any
initial condition V (0) belonging to a sublevel set ΛZI
of the shifted Lyapunov function M (22) with Pl∗ = 0
contained in Rn>0 always remain in ΛZI and converge
to an asymptotically stable equilibrium belonging to
EZI ∩ ΛZI ∩ VZI .

The first result we present concerns the set of equilibria,
which follows by adapting the proof of Lemma 1.
Lemma 5 (Equilibria for ZI loads) The set of equilibria of system (9), (10), (11) with Il (Vl ) = Il∗ − Yl∗ Vl is
EZI = {V ∈ Rn>0 : PZI (Vs ) = 0,
Vl = (Yll + Yl∗ )−1 (Il∗ − Yls Vs )},
where PZI (Vs ) depicts the power balance at the sources

Proof. First of all we observe that the proof of Theorem
4 holds for the case of ZI loads (it suffices to set Pl∗ = 0
and Il (Vl ) = Il∗ − Yl∗ Vl throughout the proof). As an
additional feature of ZI loads (to be exploited below) we
can explicitly construct δ(Vs ) = (Yll +Yl∗ )−1 (Il∗ −Yls Vs ).

PZI (V ) = [Vs ]Ŷred Vs+ [Vs ]Ysl (Yll + Yl∗ )−1 Il∗ − Ps ,
|{z}
| {z } |
{z
}
network
dissipation

load
demands

source
injections

Ŷred = Yss − Ysl (Yll + Yl∗ )−1 Yls is the Kron-reduced
conductance matrix that also absorbed the constant
impedance loads, and Ps is vector of power injections by
the sources written for V ∈ EZI as Cs 1p∗s , with
p∗s = −1T

From the proof of Theorem 4 (specialized to the case of
ZI loads), it is known that any solution Vs of the ODE
(29) is bounded. By Birckhoff’s Lemma ([27, Lemma
3.1]) the positive limit set Ω(Vs ) associated with a solution Vs (t) is non-empty, compact, and invariant. Moreover, it is contained in EZI ∩ ΛZI ∩ VZI . We would like
to prove that Ω(Vs ) is a singleton. To this end, and similarly to [16] we appeal to [28, Proposition 4.7], which
states that if the positive limit set Ω(Vs ) of a trajectory contains a Lyapunov stable equilibrium V s , then
Ω(Vs ) = {V s }. To see this first notice that V s being in
Ω(Vs ) and hence in EZI ∩ ΛZI ∩ VZI , it is indeed an equilibrium of the system. Thus, following (30), one can construct a shifted function N (Vs ) associated to V s . The
explicit expression of N (Vs ) is given by

Il∗ − Yl∗ (Yll + Yl∗ )−1 (Il∗ − Yls Vs )
.
1T [Vs ]−1 Cs 1

We remark that in the ZI case the equations PZI (Vs ) = 0
depend on the source voltages only, and once a solution
to it is determined, the corresponding voltages at the
loads are obtained as Vl = (Yll + Yl∗ )−1 (Il∗ − Yls Vs )
thereby explicitly solving previous IZI (V ) = 0.
Our second result concerns the convergence of the dynamics. In the case of ZI loads, convergence can be established without the definiteness condition on the equivalent conductance matrix Yeq in (27). Indeed, for Pl∗ = 0,
the condition (26) is automatically satisfied. Before, this
condition was needed to certify strict convexity of the
shifted Lyapunov function M (see (28)) as well as the

T

T

T

N (Vs ) = −P s ln(Vs )+P s ln(V s )+P s [V s ]−1 (Vs −V s )
"
#T "
#"
#
Vs − V s
Yss
Ysl
Vs − V s
1
+
.
2 δ(Vs ) − δ(V s )
Yls Yll + Y ∗ δ(Vs ) − δ(V s )
l

10

The gradient of N (Vs ) is given by
∂N
= −[Vs ]−1 P s + [V s ]−1 P s +
∂Vs
T

∂δ
(Vs − V s )+
Yss + Ysl
∂Vs

T
∗ ∂δ
Yls + (Yll + Yl )
(δ(Vs ) − δ(V s )) .
∂Vs
∂δ
Since ∂V
= −(Yll + Yl∗ )−1 Yls , the last summand above
s
vanishes. With the shorthand Ŷred = Yss − Ysl (Yll +
Yl∗ )−1 Yls , the gradient simplifies as

∂N
= −[Vs ]−1 P s + [V s ]−1 P s + Ŷred (Vs − V s )
∂Vs
Fig. 3. Voltage plots of the simulation

∂N
Note that the gradient ∂V
vanishes if Vs = V s and N
s
has a strict local minimum at V s since

∂2N
= Ŷred + [Vs ]−2 P s .
∂Vs2
By (31), Ṅ ≤ 0, and these two properties (properness and the nonnegative time derivative) show that V s
is a Lyapunov stable equilibrium. Therefore, Ω(Vs ) =
{V s }, and the solution Vs (t) converges to an equilibrium
point. Because Vs (t) is the Vs component of the solution to the DAE, and since Vl satisfies Vl = δ(Vs ) =
(Yll + Yl∗ )−1 (Il∗ − Yls Vs ) we also see that the solution
(Vs (t), Vl (t)) of the DAE (9), (10), (11) converges to a
point in EZI ∩ ΛZI ∩ VZI . Since this equilibrium point
is Lyapunov stable by (28) (with Pl∗ = 0) and (31), the
limit point is also asymptotically stable.

Fig. 4. Power plots of the simulation

5

Simulations

In this section, we present simulation results comparing the proposed control strategy to an averaging-based
control method. We use an example network obtained
from [4]. The network topology is sketched in Fig. 5,
and the physical parameters are given in Table 1. It can
be checked that condition (26) is satisfied for this network. As in the reference experiment, there are seven
constant power loads, five of which are initially turned
off and are turned on gradually between 9.5 and 10.5 ms.
This means that there is a gradual increase of the total power load from 70 W to 245 W. We simulate both
the proposed control strategy (5), and the distributed
averaging integral controller (8) for a comparison. The
power measured at the source nodes is shown for both
control strategies in Fig. 4. As predicted by the analysis,
at steady state proportional power sharing is achieved
by the power sources in conformity with (7). We also observe that the two controllers perform similarly, only a

5

4

6

1

2

7

8

9

3

10

Fig. 5. The node network used for the simulations. Sources
are depicted as circles, loads as rectangles. Solid lines denote
the interconnecting lines, while dashed blue lines represent
the communication graph used by the controllers.

slight overshoot for the integral controller at the power
source 2 can be observed. The voltage evolution both at
the sources and at the loads is depicted in Fig. 3.
6

Conclusions

We have proposed controllers for DC microgrids that
average power measurement at the sources. The results
apply to network preserved model (systems of DAE) of

11

Parameter

Value

Transmission line weights Γi

6 × 10−1 Ω
√
4 × 10−2 kgm/s
√
8 × 10−2 kgm/s

Capacitance weight Ci , i = 1, 3
i=2
Nominal voltage V

∗

[6] S. Moayedi and A. Davoudi. Distributed tertiary control
of DC microgrid clusters. IEEE Transactions on Power
Electronics, 31(2):1717–1733, Feb 2016.
[7] L. Meng, T. Dragicevic, J. Roldán-Pérez, J. C. Vasquez,
and J. M. Guerrero. Modeling and sensitivity study of
consensus algorithm-based distributed hierarchical control
for DC microgrids. IEEE Transactions on Smart Grid,
7(3):1504–1515, May 2016.
[8] N. Barabanov, R. Ortega, R. Grino, and B. Polyak. On
existence and stability of equilibria of linear time-invariant
systems with constant power loads. IEEE Transactions on
Circuits and Systems I: Regular Papers, 63(1):114–121, Jan
2016.
[9] J. W. Simpson-Porco, F. Dörfler, and F. Bullo. On resistive
networks of constant-power devices. IEEE Transactions on
Circuits and Systems II: Express Briefs, 62(8):811–815, Aug
2015.
[10] Javad Lavei, Anders Rantzer, and Stephen Low. Power flow
optimization using positive quadratic programming. IFAC
Proceedings Volumes, 44(1):10481 – 10486, 2011. 18th IFAC
World Congress.
[11] Alain Sarlette, Jing Dai, Yannick Phulpin, and Damien Ernst.
Cooperative frequency control with a multi-terminal highvoltage DC network. Automatica, 48(12):3128 – 3134, 2012.
[12] Martin Andreasson, Mohammad Nazari, Dimos V.
Dimarogonas, Henrik Sandberg, Karl H. Johansson, and
Mehrdad Ghandhari. Distributed voltage and current control
of multi-terminal high-voltage direct current transmission
systems. IFAC Proceedings Volumes, 47(3):11910 – 11916,
2014. 19th IFAC World Congress.
[13] Daniele Zonetti, Romeo Ortega, and Abdelkrim Benchaib.
Modeling and control of HVDC transmission systems from
theory to practice and back. Control Engineering Practice,
45:133 – 146, 2015.
[14] Daniele Zonetti, Romeo Ortega, and Johannes Schiffer. A
tool for stability and power sharing analysis of a generalized
class of droop controllers for high-voltage direct-current
transmission systems. CoRR, abs/1609.03149, 2016.
[15] Daniele Zonetti. Energy-based modelling and control of
electric power systems with guaranteed stability properties.
Theses, Université Paris-Saclay, April 2016.
[16] C. De Persis and N. Monshizadeh. A modular design of
incremental Lyapunov functions for microgrid control with
power sharing. In 2016 European Control Conference (ECC),
pages 1501–1506, July 2016.
[17] J. Schiffer, T. Seel, J. Raisch, and T. Sezi. Voltage stability
and reactive power sharing in inverter-based microgrids
with consensus-based distributed voltage control. IEEE
Transactions on Control Systems Technology, 24(1):96–109,
Jan 2016.
[18] F. Dörfler and F. Bullo. Kron reduction of graphs with
applications to electrical networks. 60(1):150–163, January
2013.
[19] H.-D. Chiang. Direct Methods for Stability Analysis of
Electric Power Systems. John Wiley & Sons, 2011.
[20] D. Bauso, L. Giarre, and R. Pesenti. Non-linear protocols for
optimal distributed consensus in networks of dynamic agents.
Systems & Control Letters, 55(11):918 – 928, 2006.
[21] Jorge Cortes. Distributed algorithms for reaching consensus
on general functions. Automatica, 44(3):726 – 737, 2008.
[22] S. Bolognani and S. Zampieri. On the existence and
linear approximation of the power flow solution in power
distribution networks. IEEE Transactions on Power Systems,
31(1):163–172, Jan 2016.

48 V

Integral controller weights Di

1 × 10−4

Load values −Pl∗

35 W

Table 1
Simulation parameter values.

the microgrid in the presence of ZIP loads. Capacitors
at the terminals of the grid that model either Π-models
of lines or power converter components can be included
by means of passivity-based analysis.
Many interesting new research directions can be taken.
The first one is to consider more complex scenarios such
as the inclusion of dynamical (inductive) lines and loads.
Another one is the extensions of the controllers to network preserved AC microgrids. Moreover, although the
preservation of the geometric mean of the voltages allows for an estimate of the voltage excursion, no active
voltage regulation is present in the proposed scheme. An
addition of voltage controllers to the power consensus algorithm is an interesting open problem. The power consensus algorithms lead to a new set of power flow equations, whose solvability still needs to be investigated,
e.g., starting from recent advances concerning power flow
feasibility and approximations; see [22,8,23] and references therein. The distributed averaging integral controller (8) discussed in Remark 3 enjoys the nice feature
of not requiring power measurements and could be an
enthralling algorithm to investigate further. Finally, the
power consensus algorithms preserves the weighted geometric mean of the voltages and is thus a compelling
application for nonlinear consensus schemes [20,21]. We
believe this connection deserves a deeper investigation.
References
[1] P. C. Loh, D. Li, Y. K. Chai, and F. Blaabjerg. Autonomous
operation of hybrid microgrid with AC and DC subgrids.
IEEE Transactions on Power Electronics, 28(5):2214–2223,
May 2013.
[2] V. Nasirian, S. Moayedi, A. Davoudi, and F. L. Lewis.
Distributed cooperative control of DC microgrids. IEEE
Transactions on Power Electronics, 30(4):2288–2303, 2015.
[3] Jinxin Zhao and Florian Dörfler. Distributed control and
optimization in DC microgrids. Automatica, 61:18–26, 2015.
[4] Julia A Belk, Wardah Inam, David J Perreault, and
Konstantin Turitsyn. Stability and control of ad hoc DC
microgrids. arXiv preprint arXiv:1603.05289, 2016.
[5] Michele Tucci, Lexuan Meng, Josep M Guerrero, and
Giancarlo Ferrari-Trecate. Consensus algorithms and plugand-play control for current sharing in DC microgrids. arXiv
preprint arXiv:1603.03624, 2016.

12

[23] J. W. Simpson-Porco, F. Dörfler, and F. Bullo. Voltage
collapse in complex power grids. Nature Communications,
7(10790), 2016.

[26] J. Schiffer and F. Dörfler. On stability of a distributed
averaging PI frequency and active power controlled
differential-algebraic power system model. In 2016 European
Control Conference (ECC), pages 1487–1492, July 2016.

[24] R. Abraham, J. E. Marsden, and T. S. Ratiu. Manifolds,
Tensor Analysis, and Applications, volume 75 of Appl. Math.
Sci. Springer, 2 edition, 1988.

[27] H.K. Khalil. Nonlinear Systems. Prentice Hall, 1996.

[25] David J Hill and Iven MY Mareels. Stability theory
for differential/algebraic systems with application to power
systems. IEEE Transactions on Circuits and Systems,
37(11):1416–1423, 1990.

[28] W.M. Haddad and V.S. Chellaboina. Nonlinear dynamical
systems and control: a Lyapunov-based approach. Princeton
University Press, 2008.

13

