NREL is a national laboratory of the U.S. Department of Energy

Office of Energy Efficiency & Renewable Energy

Operated by the Alliance for Sustainable Energy, LLC

This report is available at no cost from the National Renewable Energy

Laboratory (NREL) at www.nrel.gov/publications.

Contract No. DE-AC36-08GO28308

Conference Paper

NREL/CP-5D00-76367

September 2020

Distributed Minimization of the Power 

Generation Cost in Prosumer-Based 

Distribution Networks

Preprint

Guido Cavraro,1 Andrey Bernstein,1 Ruggero Carli,2 

and Sandro Zampieri2 

1 National Renewable Energy Laboratory

2 University of Padova, Italy 

Presented at the 2020 American Control Conference (ACC)

July 1–3, 2020 
NREL is a national laboratory of the U.S. Department of Energy

Office of Energy Efficiency & Renewable Energy

Operated by the Alliance for Sustainable Energy, LLC

This report is available at no cost from the National Renewable Energy

Laboratory (NREL) at www.nrel.gov/publications.

Contract No. DE-AC36-08GO28308

National Renewable Energy Laboratory

15013 Denver West Parkway

Golden, CO 80401

303-275-3000 • www.nrel.gov

Conference Paper

NREL/CP-5D00-76367

September 2020

Distributed Minimization of the Power 

Generation Cost in Prosumer-Based 

Distribution Networks

Preprint

Guido Cavraro,1 Andrey Bernstein,1 Ruggero Carli,2 

and Sandro Zampieri2 

1 National Renewable Energy Laboratory

2 University of Padova, Italy 

Suggested Citation 

Cavraro, Guido, Andrey Bernstein, Ruggero Carli, and Sandro Zampieri. 2020. Distributed 

Minimization of the Power Generation Cost in Prosumer-Based Distribution Networks: 

Preprint. Golden, CO: National Renewable Energy Laboratory. NREL/CP-5D00-76367. 

https://www.nrel.gov/docs/fy20osti/76367.pdf. 
NOTICE

This work was authored in part by the National Renewable Energy Laboratory, operated by Alliance for Sustainable 

Energy, LLC, for the U.S. Department of Energy (DOE) under Contract No. DE-AC36-08GO28308. Funding 

provided by U.S. Department of Energy and supported in part by the Laboratory Directed Research and 

Development (LDRD) Program at NREL. The views expressed herein do not necessarily represent the views of the 

DOE or the U.S. Government.

This report is available at no cost from the National Renewable 

Energy Laboratory (NREL) at www.nrel.gov/publications. 

U.S. Department of Energy (DOE) reports produced after 1991 

and a growing number of pre-1991 documents are available 

free via www.OSTI.gov. 

Cover Photos by Dennis Schroeder: (clockwise, left to right) NREL 51934, NREL 45897, NREL 42160, NREL 45891, NREL 48097, 

NREL 46526.

NREL prints on paper that contains recycled content. 
Distributed Minimization of the Power Generation Cost in

Prosumer-Based Distribution Networks

Guido Cavraro, Andrey Bernstein, Ruggero Carli, and Sandro Zampieri

Abstract— Traditionally, electrical power was generated in

big power plants. The cost of producing energy was related

to the cost of fuel, e.g., carbon or gas, and by the cost of

maintaining the power plants. With the advent of distributed

energy resources, power can be produced directly at the edge of

the electrical network by a new type of agents: the prosumers.

Prosumers are entities that both consume and generate power,

e.g., by means of photovoltaic panels. The cost of the power

produced by prosumers is no longer related to fuel consumption

since energy coming from distributed generators is essentially

free. Rather, the cost is related to the remuneration that is due

to the prosumers for the services they provide. The proposed

control strategy minimizes the active power generation cost

in the aforementioned scenario. The control scheme requires

that the prosumers measure their voltage and then adjust

the amount of injected power, according to a continuous time

feedback control law that is indeed a projected gradient descent

strategy. Simulations are provided in order to illustrate the

algorithm behavior.

I. INTRODUCTION

The advent of distributed energy resources (DERs), like

wind turbines, photovoltaic panels, or other renewable energy

sources, is going to dramatically change the actual electrical

grid [1] as a massive number of small power generators

will be deployed in the distribution grid. On one hand,

a number of improvements for the electrical system will

be achieved with the integration of the distributed energy

resources in the distribution network, e.g. voltage profile

improvements, reduction of line loss and power generation

cost and other ancillary services [2]. On the other hand,

distribution networks can suffer instabilities or damages

if the (possibly) millions hosted DERs do not coordinate

their power injections. This motivates the recent interest in

developing algorithms that solve optimal power flow (OPF)

problems for distribution networks, while in the past they

were focused only on transmission networks. Utilities have

G. Cavraro and A. Bernstein are with the Power Systems Engineering

Center, National Renewable Energy Laboratory (NREL), Golden, USA.

Email: name.surname@nrel.gov.

R. Carli and S. Zampieri are with the Dept. of Information Engineering, University of Padova, Italy. Email: {carlirug,

zampi}@dei.unipd.it.

This work was authored in part by the National Renewable Energy

Laboratory (NREL), operated by Alliance for Sustainable Energy, LLC,

for the U.S. Department of Energy (DOE) under Contract DE-AC36-

08GO28308 and supported in part by the Laboratory Directed Research

and Development (LDRD) Program at NREL. The views expressed in

the paper do not necessarily represent the views of the DOE or the

U.S. Government. The U.S. Government retains and the publisher, by

accepting the paper for publication, acknowledges that the U.S. Government

retains a nonexclusive, paid-up, irrevocable, worldwide license to publish

or reproduce the published form of this work, or allow others to do so, for

U.S. Government purposes.

to solve optimal power flow (OPF) problems to find a power

system operating point that minimizes a cost function, e.g.,

the generation cost or the line losses, satisfying the power

demand and some operative constraints, like limitation on

the generation capability or voltage magnitude limits.

OPF problems for transmission networks are usually tackled by centralized off-line algorithms which (i) collect all the

necessary field data, (ii) compute the optimal configuration,

and (iii) dispatch the power production to the generators.

Similar off-line solution strategy has been applied to distribution networks too, e.g., algorithms based on the alternating

direction method of multipliers (ADMM) has been developed

[3]. They require a large number of iterations and a high

computational burden to converge, mainly due to the nonlinear relations among powers and voltages which make the

OPF problem non-convex. To overcome these drawbacks, the

OPF problem has been reformulated as a rank-constrained

semidefinite program conveniently convexified by dropping

the rank constraint [4]. Heed that in all the aforementioned

solutions, the OPF solution is applied at the end of the

optimization process, that requires either a global knowledge

of the system state or a number of communication rounds

necessary for the solution numerical computation.

This approach is not practical in distribution networks,

because of the power demand fast variability and in the generators’ generation capabilities, and to the fact that generators

can connect or disconnect, requiring an automatic reconfiguration of the grid control infrastructure (the so called “plug

and play” approach). This explains the recent interest in

online distributed feedback optimization algorithms, where

the feedback is exploited to infer from local measurements,

e.g., voltages or power injections, global information [5]–[7].

The algorithm proposed herein is (i) a continuous time

control algorithm, that exploit as feedback variables the

phasorial voltages and information about the power flowing

through the substation (also termed PCC in the following),

(ii) is inherently an online algorithm, and (iii) can be implemented in a distributed way. The OPF problem we consider

aims at minimizing the global power generation cost. In our

setup we consider two types of cost functions: one associated

with the power produced by conventional power plants

and coming from the external transmission network, and

one associated with the DERs dispersed in the distribution

network. The problem is tackled via a projected gradientbased approach that is shown to be provably convergent to

an approximated optimal solution of the OPF problem. The

mathematical proofs of the resulted presented next can be

found in the paper extended version [8].

This report is available at no cost from the National Renewable Energy Laboratory at www.nrel.gov/publications. 1
II. NOTATION AND MATHEMATICAL PRELIMINARIES

Upper- (lower-) case boldface letters denote matrices (column vectors). Sets are denoted with calligraphic symbols.

Symbol > stands for transposition. Vectors 1 and em represent the all-one vectors and the m-th canonical vector,

respectively. Given a set A, |A| is the number of elements

in A. The projection of a vector x0 onto a convex set X is

denoted as ΠX (x0) = arg minx∈X kx − x0k.

A convex cone K is a set such that, for every x1, x2 ∈ K,

ax1 ∈ K, a ≥ 0, and ax1 + (1 − a)x2 ∈ K, 0 ≤ a ≤ 1.

Given a convex cone K ⊂ RN , its polar cone K∗is the set

K∗ = {y : y

>x ≤ 0, x ∈ K}.

Lemma 1 (Lemma 1, [7]): Consider a convex cone K ⊂

RN and a vector a ∈ RN . Let b = ΠK(a). Then, it holds

a − b ∈ K∗.

The tangent cone of a non-empty convex set X at x ∈ X is

the convex cone

T

x

X = cl{d ∈ RN : ∃ > 0, x + d ∈ X }.

The normal cone of a non-empty convex set X at x ∈ X is

the convex cone

N

x

X = {d ∈ RN : d

>(y − x) ≤ 0, ∀y ∈ X }.

The tangent and the normal cone of X at x are such that

N x

X = (T

x

X )

∗

[9].

III. CYBER-PHYSICAL MODEL OF A DISTRIBUTION GRID

A smart power distribution network can be described as a

cyber-physical system in which

• The physical layer comprises the distribution infrastructure: lines, loads, microgenerators, and the substation;

• The cyber layer consists of intelligent agents provided

with actuation, communication, sensing, and computational capabilities.

The physical layer can be modeled using an undirected

graph G = (V, E); the set V, with cardinality |V| := n + 1,

collects nodes that are associated with electrical buses. The

substation, also termed point of common coupling (PCC),

is labeled as 0. As a consequence, we can write V =

{0, 1, . . . , n}. The set E, E ⊆ V × V, with cardinality

|E| := m, collects the edges that are associated with the

electrical lines. Usually, distribution networks feature a radial

topology; thus, m = n. Define the following quantities:

• u ∈ C

n, with uv being the voltage at bus v;

• i ∈ C

n, where iv is the current injected at node v;

• s = p+iq ∈ C

n, where sv, pv and qv are the complex,

the active and the reactive power injected at node v.

Active and reactive powers will take positive values, i.e.,

pv, qv ≥ 0 when they are injected in the grid by bus v.

Conversely, pv, qv ≤ 0 means that bus v is absorbing

power from the grid.

• z = r + ix ∈ C

m, where z`, r` and x` are the

impedance, the resistance and the reactance of line `.

The system state is described by the following system of

equations:



i0

i



= Y



u0

u



(1a)

u0 = UN (1b)

uv

 ̄iv = pv + iqv v 6= 0 . (1c)

Equation (1a), where Y is the bus admittance matrix,

provides the relation between voltages and currents. Equation

(1b) holds since the substation is modeled as an ideal voltage

generator imposing the nominal voltage u0 = UN . Finally,

equation (1c) comes from the fact that all the nodes, except

the PCC, are modeled as constant power buses. It can be

shown that there exists a unique matrix Z ∈ C

n×n, called the

Green matrix, that allows us to write voltages as a function

of currents via:

u = Zi + 1UN . (2)

Notably, the Z is symmetric and positive semidefinite (see

[10]). Furthermore, it has been shown that the non-linear

relation between power injections and voltages satisfies [10]

u =

(R + iX)

UN

(p − iq) + 1UN + o



1

UN



where the little-o notation means that limUN→∞

o(f(UN ))

f(UN ) =

0. As a consequence, the voltage-powers relation can be

approximated with the linear function

u =

(R + iX)

UN

(p − iq) + 1Un (3)

where the matrices R and X are the real and the imaginary

part of Z, respectively. Power loss equals the sum of the

powers injected by the substation and all the prosumers; thus,

p0 can be written as

p0(p) = −

Xn

v=1

pv + `(p) = −1

>p + `(p). (4)

where `(p) denotes the active power loss. The dependence of

power loss from q is dropped, since in this paper the reactive

power injections q are considered as fixed, non-controllable

variables. The linearization (3) can be used to find convex

approximations of the active power network loss

`(p) = p

> R

U

2

N

p + q

> X

U

2

N

q. (5)

and of the power flowing through the substation

p0(p) = −1

>p + p> R

U

2

N

p + q

> X

U

2

N

q. (6)

We assume that every bus in the power network is a prosumer

[11], except the PCC. A prosumer is a bus able to generate

the power p

g

v ≥ 0 (v is a producer) but with a power demand

p

d

v

, qv ≤ 0 too. That is, a prosumer v is both a producer and a

consumer of energy. Define, pv:= p

d

v+p

g

v

, i.e., pvis the sum

between the power demand and the generation capability of

prosumer v. We have two cases:

• When pv ≥ 0, the generation capability of v exceeds the

power demand and v behaves like a generator. Its power

output pv is controllable and satisfies 0 ≤ pv ≤ pv;

• When pv ≤ 0, the power demand exceeds the generation

capability of v and v behaves like an uncontrollable load

requiring pv = pv

.

This report is available at no cost from the National Renewable Energy Laboratory at www.nrel.gov/publications. 2
Fig. 1. Power cost function associated with a prosumer.

Each prosumer (referred to also as agent) is provided with

sensing capability in the form of a phasor measurement

unit (PMU). Finally, we assume that the power grid has a

Network Supervisor (NS) able to measure the amount of

power flowing into the distribution grid through the PCC.

The NS acts as an intermediary between the prosumers and

the distribution system operator, receiving updates on the

generation cost of p0, and broadcasting to all the prosumers

information that will be used in the control algorithm.

IV. PROBLEM FORMULATION

Next, a distributed feedback control algorithm will be devised to minimize the power generation cost. The generation

cost is often chosen to be a polynomial, e.g., a quadratic

function of the type [4], [12], [13]

fv(pv) = (

αv,2p

2

v + αv,1pv + αv,0 pv ≥ 0

0 pv < 0

. (7)

This is a common choice in the literature that well describe

the case in which power is produced in power plants;

however, it does not capture the case in which power is

produced by a DER, like a photovoltaic panel. In the latter

case, energy is produced at zero cost; rather, the power cost is

the remuneration due to the generator owners for the power

injection. In summary,

• When a prosumer is behaving like a generator, the

remuneration is proportional to the quantity of energy

injected into the grid, with a proportionality constant

that depends on the contractual agreement with the

utility;

• When a prosumer is behaving like a load, it will not

receive a remuneration; instead, the prosumer has to pay

the power consumed, whose price is proportional to the

power absorbed. As a result, the utility must sustain a

negative price, i.e. it earns money.

The model we use for the remuneration (see Figure 1) comes

straightforwardly from the previous reasoning, and it is

ˆfv(pv) = (

cv,G pv pv ≥ 0

cv,L pv pv < 0

∀v ∈ V (8)

In this framework there are mainly two possible scenarios.

The first one is related to the “prosumers point of view”,

in which the prosumers compete in order to maximize the

profit they earn from selling power to the grid. Overall, this

boils down to a game among the agents. A first treatment of

this scenario can be found in [14].

The second scenario is instead related to the “utility point

of view”, where the total cost accounts for the production

cost of the energy injected by the PCC (that comes from big

generation plants such as nuclear or hydroelectrical plants)

and for the remuneration to be paid to the owners of DRES.

In this framework, the goal of the utility is to minimize the

total cost while satisfying some operative constraints.

In this paper we focus on the second scenario. To have

a lighter notation, we assume that all the agents are paid in

the same way, that is,

cv,G = cG ∀v ∈ V\{0} (9)

and we model f0(p0) as

f0(p0) = (

c0p0 p0 ≥ 0

0 p0 < 0

(10)

However, our results can be easily generalized to the case

in which (9) does not hold and the power generation cost at

the PCC is modeled with classic cost functions, e.g., (7). We

make also the following assumption.

Assumption 2: The reward cG is positive, that is cG > 0.

Formally, the problem we are interested into can be stated

as the following optimization problem

min

p

Xn

v=1

fv(pv) + f0(p0) (11a)

s. t. p0(p) = −1

>p + p> R

U

2

N

p + q

> X

U

2

N

q (11b)

p ∈ B (11c)

where

• equation (11b) models active power conservation;

• B is the feasible set, that is

B = {p : p

v

≤ pv ≤ pv}.

If prosumer v behaves as a load, then p

v

= pv < 0.

If prosumer v behaves as a generator, then p

v

= 0 and

0 < pv < ∞. Thus, the set B is bounded.

• the cost functions fv(pv)’s are defined as

fv(pv) := cGpv. (12)

The function fv(pv) is different from the piece-wise

actual cost ˆfv(pv): fv(pv) is linear. Nonetheless, the

minimizer of (11) is the same as if the ˆfv(pv)’s would

have been used, even though the minimum is different.

In fact, fv(pv) = ˆfv(pv) for pv ≥ 0, i.e., where pv is an

actual optimization variable; while fv(pv) 6= ˆfv(pv) for

pv < 0, where constraint (11c) forces pv = p

v

= pv.

Define f(p) = Pn

v=1 fv(pv) + f0(p0). Thanks to (10) and

(12), we can write

f(p) = (

f

+(p) if p ∈ S+ ∪ S0

f

−(p) if p ∈ S− ∪ S0

(13)

This report is available at no cost from the National Renewable Energy Laboratory at www.nrel.gov/publications. 3
where

f

+(p) = (cG − c0)1T p + c0



p

> R

U

2

N

p + q

> X

U

2

N

q



f

−(p) = cG1T p

S

+ = {p : p0(p) > 0} S− = {p : p0(p) < 0}

S

0 = {p : p0(p) = 0}.

S

0

is a variety separating S

+ from S−. The cost f(p) is

a continuous function: f

+(p) and f−(p) are continuous in

S

+ ∪ S0

and in S

− ∪ S0

, respectively; moreover, f

+(p) =

f

−(p) for all p ∈ S0

. Furthermore, f(p) is differentiable in

S

+ and in S−, but not in S0

. The solution of problem 11 is

characterized by the following proposition, proved in [6].

Proposition 3: Consider problem (11), let Assumption 2

hold and denote with p

∗

an optimal solution of (11). Then,

every optimal configuration p

∗

is such that p0(p

∗

) ≥ 0, i.e.,

it belongs to S

+ ∪ S0

.

V. A CONTINUOUS-TIME ALGORITHM

Problem (11) could be solved by a continuous-time projected gradient descent of the cost function. However, remind

that Proposition 3 ensures that the optimum configuration lies

in S

+∪S0

. To speed up the rate of convergence, a reasonable

strategy is to steer the prosumers power injections p in the

region S

+ ∪ S0

as fast as possible when p ∈ S−, instead of

executing the mere gradient descent on f. Hence, we propose

the following strategy:

• if p0(p) < 0, agents aims at increasing p0 following

the rule

p ̇(t) = ΠT

p(t)

B

(∇p0(p(t))) (14)

• if p0(p) ≥ 0, agents perform the cost function gradient

descent.

p ̇(t) = ΠT

p(t)

B

(−∇f

+(p(t))) (15)

The overall control rule can be written as

p ̇(t) = ΠT

p(t)

B

(φ(p(t))) (16)

where

φ(p) = (

−∇f

+(p) = (c0 − cG)1 −

c0

U2

N

Rp p ∈ S+

∇p0(p) = −1 +

1

U2

N

Rp p ∈ S−.

(17)

Heed that the vector field (17) is continuous on the sets

S

+ and S−, but it is discontinuous and not defined on

S

0

, since f is not differentiable in S

0

. A widely accepted

way to extend (17) on S

0

is given by the Filippov convex

method, which considers the following set-valued extension,

or convexification, of φ(p) [15]. Precisely, let φ

0

(p) be

the extension of ∇f(p) on S

0

. According to the Filippov

method,

φ

0

(p) ∈ co{∇p0(p), −∇f

+(p)} (18)

where co{∇p0(p), −∇f

+(p)} denotes the smallest convex

set that contains ∇p0(p) and −∇f

+(p). Hence, for p ∈ S0

,

we can write

φ

0

(p) = −α(p)∇f

+(p) + (1 − α(p))∇p0(p) (19)

where α(p)) is a suitable non-negative variable, 0 ≤ α(p) ≤

1, which can take different values leading to:

• Attractive sliding mode. When both −∇f

+(p) and

∇p0(p) are pointing towards S

0

, or

−∇f

+(p)>∇p0(p) < 0 (20)

an attractive sliding mode occurs. When (20) is satisfied, p will move along S

0 giving rise to a sliding

motion. During the sliding motion, α(p) is such that

φ

0

(p) is tangent to S

0

, or φ

0

(p)

>∇p0(p) = 0. It can

be easily shown that

α(p) = k∇p0(p)k

2

∇p0(p)>(∇f +(p) + ∇p0(p)) (21)

and that α(p) > 0.

• Transversal intersection. When ∇f

+(p) and ∇p0(p)

have the same orientation, that is

−∇f

+(p)>∇po(p) > 0. (22)

there is no sliding motion and the trajectory of p(·) will

leave the surface S

0

. The parameter α(p) is chosen as

α(p) = 1.

Based on the previous considerations, the gradient field (17)

can be extended as

φ(p) =







−∇f

+(p) p ∈ S+

∇p0(p) p ∈ S−

φ

0

(p) p ∈ S0

(23)

Beside the discontinuity of φ(p) in S

0

, the dynamic (16)

may be discontinuous also on the boundaries of the feasible

set B, where the gradient is projected on the cone of the

feasible directions. For these reasons, the standard notion of

solution for ordinary differential equations does not apply. In

this paper, the trajectory p : [0, t] → R

n is a Caratheodory

solution on the interval [0, t] if it is absolutely continuous

and satisfies (16) almost everywhere in [0, t]. The system

considered in this paper guarantees the uniqueness of the

Caratheodory solution, see [16].

The next Proposition, which is the main result of this

paper, shows that the dynamic (16) is asymptotically stable

and that steers the system towards the solution of (11).

Proposition 4: Consider problem (11) and the dynamic

(16). Let Assumption 2 hold. Then, (11) has a unique

solution p

∗

and the trajectory t → p(t) is asymptotically

stable and converges to the unique optimum p

∗

.

VI. A DISCRETE-TIME FEEDBACK CONTROL

IMPLEMENTATION

The control strategy (16) cannot be actually implemented,

since inverters or power electronic devices in general perform control actions on discrete times. For this reason, we

introduce the discrete-time version (16):

p(t + 1) = (

ΠB



p(t) − +∇f

+(p(t)



if p0(p) > 0

ΠB



p(t) + −∇p0(p(t)



if p0(p) < 0

(24)

This report is available at no cost from the National Renewable Energy Laboratory at www.nrel.gov/publications. 4
where +, − are a suitable positive constant. The continuous

time dynamic (16) can be seen as the limit of the discrete

time dynamic (24) when + and − tend to zero.

To avoid the computational burden of a centralized algorithm, we seek for a distributed implementation of (24), i.e.,

we aim at designing a control scheme that can be performed

by each prosumer having just a local and partial (thus not

global) knowledge of the grid state.

Note that the gradients ∇f

+(p(t)) and ∇p0(p(y)) cannot

be computed by each prosumer using only local information.

Rather, global information, i.e., the whole power injection

vector p and the grid structure R, is needed, in contrast

with our desire for designing a distributed algorithm. However, an estimate of Rp can actually be computed locally.

Assume that all the power lines in the grid have the same

inductance/resistance ratio, but possibly different impedance

magnitude, i.e., z` = e

iθ|z`| for any ` in E and for a fixed

θ. In this case, equation (3) becomes

u = e

iθ|Z|(p − iq) + 1UN (25)

where the matrix |Z| collects the absolute values of the

entries of Z. This assumption is satisfied when the inductance/resistance ratio of the power lines of the grid is relatively homogeneous, which is reasonable in many practical

cases, e.g., see the IEEE standard testbeds [17] and [18, Table

I]. Under (25), we can estimate Rp as

Rp = cos θ<(e

−iθ(u − 1UN )). (26)

Equation 26 can be used finally to obtain the approximations

∇f

+(p) = (cG − c0)1 +

c0 cos θ

U

2

N

<(e

−iθ(u − 1UN ))1

∇p0(p) = −1 +

cos θ

U

2

N

<(e

−iθ(u − 1UN ))1

Note that every prosumer can compute its component of

∇f

+(p) knowing its own voltage phasor, the prices cG and

c0, and the parameter θ, since

[∇f

+(p)]v = (cG − c0) + c0 cos θ

U

2

N

<(e

−iθ(uv − UN ))

(27)

[∇p0(p)]v = −1 +

cos θ

U

2

N

<(e

−iθ(uv − UN )). (28)

By plugging (27) and (28) into equation (24), we obtain

the control rule performed by the prosumers. Trivially, the

projection onto the feasible set can be done locally. To

perform (24), besides local information, prosumers need to

know the current power flowing the substation, the substation

voltage, and the current price c0. The NS is in charge

of metering p0 and UN and broadcasting their value and

c0 to the prosumers. Based on the above description, it

is clear what is the feedback scheme that underlies the

procedure we propose: during each iteration agents and the

NS sense the grid, communicate, compute the power setpoint and then actuate it. The control strategy, formally

described in Algorithm I, require only an estimate of the

resistance/reactance ratio to be executed. Hence, it can be

considered as a model free optimization algorithm.

ALGORITHM I

Let  be a suitable positive constant. At the (t + 1)-th cycle.

The NS performs the following actions:

1. Measure the power p0(p(t)) and the voltage UN .

2. Broadcast p0(p(t)), UN and c0.

Prosumer v performs the following actions

1. Measure its own voltage uv(t)).

2. Receive the quantities p0(p(t)), UN and c0.

3. Computes its gradient component through (27) or (28).

4. Computes, using equation (24), and applies the active

power setpoint.

1

2

3

4

5

6

7 8

Fig. 2. Schematic representation of the IEEE 37-bus test feeder. Red

nodes represent prosumers behaving as generators, as their power generation

capability exceed their power requirement. Black nodes instead represent

prosumer behaving as loads, whose power demand is bigger than their

generation capabilities.

VII. SIMULATIONS

The algorithm has been tested on a single-phase equivalent

of the IEEE 37-bus test feeder, reported in Figure 2. The load

buses are a mixture of constant-current, constant-impedance

and constant-power loads, with a total active power power

demand of almost 2 MW and reactive power demand of

almost 1 MVAR [17]. The line impedances differ from edge

to edge, however, the inductance/resistance ratio exhibits a

smaller variation, ranging from ∠z` = 0.47 to ∠z` = 0.59.

Eight prosumers, able to generate more power than they

require, behave as generators with controllable power output.

Their maximum active power capabilities have been set to

values that go from 85 to 490 kW. The algorithm presented

in Section VI was simulated on a nonlinear exact power

flow solver [19]; on the same solver, the true solution of

the optimization p

∗ was computed. The step-size  was

chosen small to make the trajectory similar to the continuous

time version. To study the convergence of the feedback

optimization algorithm, introduce the variable d = kp−p

∗k.

Firstly, we simulate the case in which cG < c0, i.e.,

the case in which it is cheapest to import power from

the external grid. Precisely, we normalized c0 = 1 and

This report is available at no cost from the National Renewable Energy Laboratory at www.nrel.gov/publications. 5
Fig. 3. An algorithm numerical simulation with cG = 1.5 and c0 = 1.

Fig. 4. An algorithm numerical simulation with cG = 0.5 and c0 = 1.

set cG = 1.5. In this case, p

∗ ∈ S+. The results are

reported in Fig. 3. Secondly, the case where c0 < cG,

i.e., when buying energy from the distributed generators is

more convenient, was inspected. Fig. 4 shows the results for

c0 = 1, cG = 0.5. In this case, p

∗ ∈ S0

and a sliding motion

appeared. Some comments are in order. First, the cost of

the equilibria configurations is slightly higher than the cost

of the true minimum, i.e., the equilibrium is suboptimum,

see Fig. 3 and Fig. 4, left panels. This is because of the

approximations used to distribute the algorithm. Second, in

the case c0 = 1, cG = 0.5, the reaching phase, i.e., the part

of the trajectories in which p reaches the surface S

0

, lasts

about 400 iterations. Once p starts sliding on the surface S

0

,

the trajectory of p slows down, see Fig. 3, right panel. Third,

when c0 = 1, cG = 0.5, and p is on the sliding surface, the

trajectory experienced the typical sliding mode chattering,

which appears because the control strategy is implemented

in discrete time rather than in continuous time, and it is

evident on the trajectory of p0, see Fig. 5.

VIII. CONCLUSION

This paper presents and analyzes a feedback control algorithm for the minimization of the generation cost in a

distribution network. The main features of our approach are

the following. Firstly, it considers the case in which the cost

of the power is not related to actual generation cost, but rather

to economical agreements between utilities and prosumers.

Secondly, it can be referred to as a “quasi local” approach.

Each prosumer is able to implement the control rule by

using only local information and some common knowledge

that it is broadcasted by the NS. Finally, simulations on

a standard testbed prove the effectiveness of the proposed

Fig. 5. Trajectory of p0(p) in the same numerical simulation depicted in

Figure 4.

control scheme.

REFERENCES

[1] A. Ipakchi and F. Albuyeh, “Grid of the future,” Power and Energy

Magazine, IEEE, vol. 7, no. 2, pp. 52–62, 2009.

[2] F. Katiraei and M. Iravani, “Power management strategies for a

microgrid with multiple distributed generation units,” Power Systems,

IEEE Transactions on, vol. 21, no. 4, pp. 1821–1831, 2006.

[3] T. Erseghe, “Distributed optimal power flow using admm,” IEEE

Trans. on Power Systems, vol. 29, no. 5, pp. 2370–2380, Sep. 2014.

[4] J. Lavaei and S. H. Low, “Zero duality gap in optimal power flow

problem,” IEEE Transactions on Power Systems, vol. 27, no. 1, pp.

92–107, Feb 2012.

[5] G. Cavraro and R. Carli, “Local and distributed voltage control

algorithms in distribution networks,” IEEE Transactions on Power

Systems, vol. 33, no. 2, pp. 1420–1430, March 2018.

[6] G. Cavraro, R. Carli, and S. Zampieri, “A distributed control algorithm

for the minimization of the power generation cost in smart micro-grid,”

in Conference on Decision and Control (CDC14), 2014.

[7] C. Chang, M. Colombino, J. Corte, and E. Dall’Anese, “Saddle-flow  ́

dynamics for distributed feedback-based optimization,” IEEE Control

Systems Letters, vol. 3, no. 4, pp. 948–953, Oct 2019.

[8] G. Cavraro, A. Bernstein, R. Carli, and S. Zampieri, “Distributed minimization of the power generation cost in prosumer-based distribution

networks,” 2020.

[9] R. T. Rockafellar and R. J.-B. Wets, Variational analysis. Springer

Science & Business Media, 2009, vol. 317.

[10] S. Bolognani and S. Zampieri, “A Distributed Control Strategy for Reactive Power Compensation in Smart Microgrids,” IEEE Transactions

on Automatic Control, vol. 58, no. 11, pp. 2818–2833, 2013.

[11] S. Grijalva and M. U. Tariq, “Prosumer-based smart grid architecture

enables a flat, sustainable electricity industry,” in Innovative Smart

Grid Technologies (ISGT), 2011 IEEE PES. IEEE, 2011, pp. 1–6.

[12] E. Devane and I. Lestas, “Stability and convergence of distributed

algorithms for the opf problem,” in 52nd IEEE Conference on Decision

and Control, 2013.

[13] E. Mallada and A. Tang, “Dynamics-aware optimal power flow,” in

52nd IEEE Conference on Decision and Control, 2013.

[14] G. Cavraro and L. Badia, “A game theory framework for active

power injection management with voltage boundary in smart grids,”

in European Control Conference ECC 2013, 2013.

[15] L. Dieci and L. Lopez, “Sliding motion in filippov differential systems:

theoretical results and a computational approach,” SIAM Journal on

Numerical Analysis, vol. 47, no. 3, pp. 2023–2051, 2009.

[16] J. Cortes, “Discontinuous dynamical systems,” IEEE Control Systems

Magazine, vol. 28, no. 3, pp. 36–73, June 2008.

[17] W. H. Kersting, “Radial distribution test feeders,” in Power Engineering Society Winter Meeting, 2001. IEEE, vol. 2. IEEE, 2001, pp.

908–912.

[18] G. Cavraro, V. Kekatos, and S. Veeramachaneni, “Voltage analytics for

power distribution network topology verification,” IEEE Trans. Smart

Grid, vol. PP, no. 99, pp. 1–1, 2017.

[19] R. D. Zimmerman, C. E. Murillo-Sanchez, and R. J. Thomas, “MAT-  ́

POWER: steady-state operations, planning and analysis tools for

power systems research and education,” IEEE Transactions on Power

Systems, vol. 26, no. 1, pp. 12–19, Feb. 2011.

This report is available at no cost from the National Renewable Energy Laboratory at www.nrel.gov/publications. 6