Identification of Operation Strategies of Distribution Networks Using a Simulated
Annealing Approach
Jorge Pereira J. Tomé Saraiva Maria Teresa Ponce de Leão
jpereira@inescn.pt jsaraiva@inescn.pt mleao@inescn.pt
FEP & INESC FEUP/DEEC & INESC FEUP/DEEC & INESC
FEP - Faculdade de Economia da Univ. do Porto
FEUP/DEEC - Faculdade de Engenharia da Univ. do Porto
INESC - Instituto de Engenharia de Sistemas e Computadores
Largo Mompilher, 22, 4050 – 392 Porto, PORTUGAL, Fax: +351.2.2084172
Abstract - In this paper we present a model to identify optimal
operation strategies of electric distribution networks
considering that one wants to minimize active power losses.
This objective can be achieved by adequately selecting
transformers taps and sections of capacitor banks that are in
operation. In order to turn the model more realistic the
developed application admits the specification of admissible
voltage ranges for each node and maximum branch currents
for each line and transformer. The integer nature of this
problem virtually turned it impossible to be solved for real
sized networks given its combinatorial nature, its complexity
and thus the involved calculation time. In this paper we
describe the use of a meta-heuristic - Simulated Annealing - to
this problem in order to address the difficulties just referred.
The use of the developed application will be illustrated using a
IEEE test network and a realistic network having 645 nodes
based on a Portuguese distribution system.
Keywords
: Loss Minimization, Volt/VAr Control,
Simulated Annealing, DMS.
I. INTRODUCTION
The operation of distribution networks turned more and
more relevant in recent years given that it is at this level
that more serious problems affecting the quality of service
usually occur and that one aims at operating them with
good economic and efficient indices. From this point of
view, the identification of optimal operation strategies
considering the minimization of active losses is becoming
more and more important together with obtaining good
voltage magnitude profiles. The operation of these
networks is strongly determined by the topology of the
system in operation, by the position of transformer taps and
by the sections of capacitor banks in operation. These
issues are directly related with several aspects associated
with the technical operation conditions of the network and
with the service quality offered to the clients. Among them
one can refer
•
in the first place, the descentralized generation of
reactive power is an important resource leading to the
reduction of reactive branch flows. This enables a
more efficient use of transmission lines since they can
be more intensively used to transmit active power
without violating thermal ratings;
•
secondly, the previous aspect can still be responsible
for the postponement of new investments. This issue
clearly has a new relevance with the advent of
dispersed generation connected to distribution
networks and with the implementation of market
mechanisms in the electric industry;
•
thirdly, the identification of the position of
transformer taps and sections of capacitor banks can
be performed considering a number of constraints not
to be violated. Among them, constraints imposing
admissible ranges for nodal voltage magnitudes and
constraints related to thermal branch flow limits are
surely important ones. The integration of the first type
of constraints can be interpreted as a way to specify
admissible deviations for voltage magnitudes, that is,
maximum values for voltage drops. This way, the
integration of these constraints is an effective way of
contributing to guarantee increased service quality
levels.
The literature on this area describes several models
aiming at optimizing the operation conditions of
distribution networks. For instance, in papers
[
1
]
to
[
5
]
one
can find descriptions of models and algorithms aiming at
identifying optimal reconfiguration strategies -
switching
actions - of distribution networks. In this sense,
reconfiguration is understood as a resource that distribution
companies can use in order to improve both their economic
and technical performances. On the other hand, papers
[
6
]
to
[
9
]
describe algorithms aiming at selecting the most
adequate position of transformer taps and number of
sections of capacitor banks in operation considering that
one wants to minimize active losses. These approaches
assume that the topology of the network in operation is
fixed.
The increasing investments directed to the distribution
sector, namely required by the higher quality standards
impose in the scope of the move to electricity markets,
enable that an higher number of equipments are
telecommanded. This is also related to the implementation
of automation function in Distribution Control Centres and
with the increased number of real time measures available
on those centres. These new capacities together with the
possibility to treat dimensionally larger problems, as the
typically related to distribution networks, are imposing the
upgrade of the traditional SCADA systems to new systems
aiming at supporting and aiding at real time the operation
activities conducted in those centres. The installation of
Distribution Management Systems in distribution networks
means that several concepts and techniques already
available in transmission systems are migrating to the
Paper BPT99–357–17 accepted for presentation at the
IEEE Power Tech ’99 Conference, Budapest, Hungary,
Aug 29 Sept 2, 1999.
distribution sector. The research reported in this paper can
be exactly integrated in this line area since we aim at
developing new applications, in some cases, or adapting
already existing ones so that they can efficiently answer to
the problems emerging in distribution networks.
In this paper we report the application of the Simulated
Annealing technique
[
10
]
to
[
12
]
to the above described
problem considering that the topology of the network is
fixed. This approach has several advantages related to the
possibility of considering the integer nature of several
variables of the problem and also with the integration in the
algorithm of several steps aiming at avoiding local optima.
Apart from this introductory section, Section II includes
a brief description of some theoretical concepts of the
Simulated Annealing technique and Section III presents the
formulation of the problem and several details related with
the practical implementation of that technique. Section IV
includes results for two networks:
•
the IEEE 24 node/36 branch test system. This network
was conveniently adapted in order to integrate 5
transformers having each one 11 taps. In this case, the
Simulated Annealing based solution is compared with
the global optimal solution identified by enumeration
of all possible configurations;
•
a real sized network based on a Portuguese
distribution system. This network has 4 transformers
having each one 21 taps, 1 transformer with 19 taps
and 4 capacitor banks each one with two possible
positions;
Finally, apart from several comments on the results
obtained for those networks, Section V includes some
conclusions.
II. THE SIMULATED ANNEALING PROCESS
A. Introduction
Real sized combinatorial problems are difficult to be
solved for real time purposes. This is the main reason why
the use of meta-heuristic procedures often consists on a fair
compromise between the quality of solutions and the
required computation time.
In the 1980's new contributions to deal with
combinatorial problems started to emerge: genetic
algorithms, neural networks, tabu-search and simulated
annealing. In this scope, the concept of annealing in
combinatorial optimization was introduced by Kirkpatrick,
Gelatt Jr. and Vecchi (1983) and independently by Cerny
(1985).
Simulated annealing (SA) appears like a flexible meta-
heuristic that is an adequate tool to solve a great number of
combinatorial problems. SA consists on a simple
framework that takes advantage of relaxing optimality in a
transitory way in order to escape from local minimum. It is
based on an iterative heuristic improvement together with
a control mechanism that allows solutions to escape from
local optima. The flexibility and simplicity of this
framework turn this meta-heuristic adequate to model
particular and, often, complex constraints, providing
solutions in acceptable computation time.
SA bases on the Metropolis method that allows, with a
certain probability, movements towards worse solutions.
Basically, the Metropolis process consists on two steps. In
the first one the temperature is raised to a state of
maximum energy while in the second step the temperature
is slowly lowered till a minimum energy state, equivalent
to the thermal equilibrium, is reached.
B. Analogy between combinatorial and thermodynamic
process
In order to more clearly explain the SA algorithm, we
will present an analogy between a combinatorial problem
and physical systems integrating a large number of
particles. In this scope, we will use a SA strategy to deal
with the identification of an operation point of a
distribution network in order to conduct an analogy
analysis with the Metropolis method related with a
physical process). This analogy can be stated as follows:
•
the alternative solutions or configurations of the
combinatorial operation distribution problem are
equivalent to the physical system states;
•
the network configuration (alternative solutions)
attributes are equivalent to the energy of different
states;
•
the control parameter is equivalent to the
temperature parameter. This parameter can be set
so that about half the new configurations analysed are
accepted at the start of the process.
The evolution of the solution algorithm is simulated
using probabilistic sampling techniques supported by
successive generation of states. This process begins with
an initial state evaluated by an energy function f
i
. After
generating and analysing a second state - energy f
j
- it is
performed an acceptation test. The acceptance of this new
solution depends on a probability computed with (1).
()
() ()
() ()()
() ()
+
−
ℜ∈
≥
≤
=
c
ifjf ife
ifjf if1
jaccept P
kk
c
j
k
fi
k
f
kk
c
(1)
c = k
B
.T
k
B
Boltzman constant
T Temperature
In order to use this algorithm, some basic decisions will
have to be taken regarding these parameters. In Table I we
present some indications regarding the analogy with the
physical process and in the next subsection we will present
some additional elements to direct the setting of these
parameters.
Table I - Keystroke Decisions For Process Convergence
General (temperature
lowering)
Problem Specifications
T
0
(initial temperature)
L
k
(iteration number)
T
k
(temperature function)
Stop criteria
i
0
(initial solution)
Neighbourhood
generation
Evaluation of solutions.
C. General Decisions
The SA cooling process integrates a number of data
regarding the initial value, the decrement function, the
number of iterations for each temperature level and the
freezing temperature scheme.
The temperature or control parameter (T
0
) should be
large enough to prevent the algorithm from being captured
by local optimum points. On the other hand, the
temperature should not be that high that it could lead to
excessive computation time. In purely planning problems,
the computation time issue may not be so relevant since we
frequently do not need solutions in real time.
The choice of the initial control parameter was
performed so that about half the configurations were
accepted in the beginning of the process. The setting of this
parameter is strongly depends on the nature of the
problem. In any case, expression (2) can be used to obtain
the value of this parameter.
()
5.0lnCT
i0
∆−=
(2)
In this expression:
•
T
0
is the initial temperature parameter
•
∆
C
i
is the cost difference between solution
i
and
solution
i-1
The number of iterations per temperature level can also
be interpreted as the time during which a constant
temperature is maintained. This number is chosen taking
into account experiences from reported in the literature.
This value can be easily obtained from initial temperature
calculations where the possible neighbours must be
explored. This number is highly dependent on the problem
size and the desired accuracy.
Two temperature schedules or cooling schemes can be
used in order to lower the temperature. The first one
corresponds to a stepwise temperature reduction scheme
while the second adopts a continuous temperature
reduction process. Previous experiences on this area
indicate that the quality of results do not benefit from
adopting the continuous reduction scheme. On the other
hand, the computation time is smaller when using a
geometric schedule. The above reasoning justified the
adoption of the first referred process to direct the cooling
scheme, according to (3).
0
i
i
TT
β=
(3)
In this expression, T
i
is the temperature parameter
during period i
and
β
i
is the temperature decreasing
parameter in period
i.
An average value for
β
was also defined. From our
experience we derived a value between 0.85 and 0.92.
Larger values would yield high computation times, and
smaller values, below 0.85, typically lead to poorer quality
solutions.
The freezing point of the problem or the stopping point
of the algorithm is determined so that, for a given
temperature level, the number of worse accepted solutions
is inferior than a pre-specified value.
D. Specific Problem Decisions
The initial solution of the algorithm can be identified
using some auxiliary algorithm. Specifically, in the
problem we are addressing, the initial solution can be
simple corresponding to a situation where all transformer
taps and capacitor banks are at a default position. This
position can, for instance, be derived from the information
available in DMS database or can correspond to the central
tap available in each transformer and to completely
disconnect all capacitor banks.
This solution is evaluated by running a power flow
exercise. From this analysis one gets the values of branch
losses, voltage magnitudes and branch currents that, apart
from other variables, are relevant to the characterization of
this operation point.
Afterwards, it is generated a neighbourhood structure
that will be used to direct the move from the current
solution to the next one. In our particular problem, the
neighbourhood solution simple consists of increasing or
decreasing one step the value of transformer taps or
sections of capacitor banks. The selection of a new
solution is performed using a probabilistic sampling
method that selects one decision variable to be altered -
that is one transformer tap or section of capacitor bank. For
the selected one, the sampling process should also indicate
if the step is to be increased or decreased.
III. FORMULATION OF THE PROBLEM
The minimization problem of active power losses can
be formulated by (4) to (9).
∑
θ−+=
nr
ijji
2
j
2
i
ij
)cos.V.V.2VV.(gz min
(4)
subj.
0)e,t,h(V,
cf
=θ
(5)
maxmin
VVV
≤≤
for each node (6)
max
ijij
II
≤
for each line (7)
{}
,...t,t,tt
3f2f1ff
∈
for each transformer with taps (8)
{}
,...e,e,ee
3c2c1cc
∈
for each cap. bank (9)
In this formulation:
•
h() represents the AC power flow equations;
•
f
t
stands for the available values of transformer taps;
•
c
e
represents the available reactive powers of
capacitor banks;
•
V and
θ
are the voltage magnitudes and phases;
•
g
ij
corresponds to the conductance of branch ij;
•
|I
ij
| is the current magnitude flowing in branch ij.
Traditionally, this problem has been solved assuming
that transformer taps and that the sections of capacitor
banks are represented by continuous variables. Such
formulation consists of a non-linear optimization problem
which can be solved using, for instance, gradient based
methods. The solution obtained this way must be
approximated to the closest discrete solution. This strategy
presents several drawbacks. On one hand the gradient
method can converge to local optima. On the other hand, as
the original problem is an integer one we cannot be sure
that the final integer solution coincides with the global
optimum one.
The application presented in this paper uses Simulated
Annealing to identify a strategy to minimize active power
losses. According to the description presented in section II,
in the next paragraphs we will detail the adopted
algorithm:
i) the initial solution corresponds to the actual
operation position for transformer taps and
capacitor sections, if available, or to default
positions;
ii) the initial solution is considered as actual (
current
x)
and set as present optimal solution identified as x*.
For the current configuration (solution) a power
flow study is run. This power flow study has the
purpose of evaluating global power losses, node
voltage magnitudes and line current flows;
iii) each solution is evaluated by the so-called
evaluation function,
)x(F
current
. In our case,
expression (10) was adopted for F.
∑
+
+
∑
+
+
∑
+
+=
brancheach for
i3
max
i
max
ii
buseach for
i2i
min
i
buseach for
i1
max
i
i
)I(.f
I
I-I
)Vmag().fVmag-(V
)Vmag().fV-(Vmag
esActiveLoss)x(F
(10)
Functions
1
f
,
2
f
and
3
f
are given by (11) to (13).
>
≤
=
max
i
i
max
i
i
i1
VVmag if ,esActiveLoss
VVmag if ,0
)Vmag(f
(11)
<
≥
=
min
i
i
min
i
i
i2
VVmag if ,esActiveLoss
VVmag if ,0
)Vmag(f
(12)
>
≤
=
max
i
i
max
i
i
i3
II if ,esActiveLoss
II if ,0
)I(f
(13)
assign
)x(F
current
to F(x*);
initialise iteration counter n at 1;
initialise iteration worst solution than current one
counter, CountWorse, at 0;
initialise Temperature parameter at 1.0;
iv) Build the neighbourhood of present solution
considering all the possible combinations of
transformer taps and capacitor sections that differ
from current one of one unit. Select a configuration
randomly within the neighbourhood of present
solution (
n
x
);
v) Run a new power flow study for the present
configuration
n
x
and evaluate this solution by
using (10). This leads to
)x(F
n
;
vi) Refresh the information related with the present
solution according to:
if
*)x(F)x(F
n
≤
then
set
n
x
to x* and to
current
x
;
set CountWorse at 0;
else
increment CountWorse;
select a random number p in
[
0.0;1.0
]
;
evaluate p(level s) computing:
×
−
=
re temperatu 00025.0
)F(x)F(x
expp(level s)
n
current
(14)
if
s) level(pp
≤
then
assign
n
x
to
current
x
;
vii) if CountWorse > maximum iteration number
without improvements then go to ix);
viii) if n > maximum iteration number for a
temperature level then
assign temperature = temperature lowering
schedule x temperature;
if temperature < final temperature then
go to ix);
set n at 1;
else increment n by 1;
go to iv);
ix) end.
According to this algorithm each solution is evaluated
using (10) that integrates four terms. The first one
corresponds to the branch active losses as computed by the
power flow run. The remaining three terms correspond to
penalty functions included whenever there are violations of
the specified range of voltage magnitudes or if the
maximum current branch flow are exceeded. It should be
referred that in this expression the values of voltage
magnitudes and branch currents are expressed in p.u.
IV. CASE STUDY
A. IEEE 24 nodes/36 branches network test
In this sub-section we will present the results obtained
for the IEEE 24 nodes/36 branches. The complete data for
this network can be found in reference [13]. This network
is includes by 31 lines and 5 transformers. The tests that
were conducted using this network considered 11 taps in
each transformer. Each tap corresponds to a 0.01 increment
within the interval [0.95,1.05]. Initially, the transformer
taps were considered at the nominal values (position 6
equivalent to 1 p.u.). Loads were set to the values referred
in [13] multiplied by a factor of 1.8. The initial value for
the acceptance function was 139112.49.
The selected Simulated Annealing parameters are:
- the number of iterations for the same temperature
level is 25;
- the initial temperature level was set to 1.0. The
lowering step determining the cooling scheme
(parameter
β
in expression 3) is 95% of previous
temperature;
- the maximum number of iterations without
improvement of the evaluation function is set to 75;
With these parameters the algorithm converged in 312
iterations and the temperature was reduced till 0.54. The
final value for the acceptance function was 137916.45. The
improvement of the evaluation function value regarding
the initial value was 1196.04.
In Figure 1 we present the evolution, along the search
process, of the actual evaluation function value
(continuous line) and for the best solution found till the
current iteration (dotted line).
137800
138200
138600
139000
139400
1 51 101 151 201 251
3
01
F(X)
iteration number
Fig. 1 – Evolution of the evaluation function for the current
solution (continuous line) and for the best-identified solution
(dotted line) for the 24 node network.
B. Network based on a Portuguese distribution system
In this subsection we will present the results obtained
for a Portuguese base distribution network. The network
has 645 nodes, 4 transformers, each one with 21 taps in the
primary level and with a voltage range within [0.85,1.15]
p.u. and a step of 0.015, one transformer with 19 taps in
the secondary level and with a voltage range within
[0.85,1.15] p.u. and a step of 0.0167. The initial position
for the transformers taps was the nominal one (position 11
and 10 respectively corresponding to the value 1.0). We
also considered 4 capacitors with two possible positions
(on or off). In the initial solution capacitors were
considered on. The initial acceptance function value was
825.83.
The Simulated Annealing parameters were set to:
- the number of iterations for the same temperature
level is 45;
- the initial temperature level was set to 1.0. The
lowering step determining the cooling scheme
(parameter
β
in expression 3) is to 95% of previous
temperature;
- the maximum number of iterations without
improvement of the evaluation function is set to
135;
With these parameters the algorithm converged in 256
iterations and the temperature lowered till 0.77. The final
value for the acceptance function was 741.19. This
represents an improvement of the evaluation function of
84.64 regarding the initial value.
7
30
7
55
7
80
8
05
8
30
1
5
1
1
01 151 201
2
51
iteration number
F
(X)
Fig. 2 – Evolution of the evaluation function for the
current
solution (continuous line) and for the best identified
solution (dotted line)
for the 645 node network.
Figure 2 presents the evolution, along the search
process, of the current value of the evaluation function
(continuous line) and for the best identified solution found
till the current iteration (dotted line). In this graphical
representation the two previously referred lines almost
coincide. This is due to the need to adequately scale the
graph in order to fit in 1 column representation. In fact, the
current solution suffered 143 changes along the whole
iterative process while the best-identified solution was only
changed 38 times.
C. Comments
Aiming at evaluating the algorithm performance in
what concerns the efficiency of the search algorithm, we
evaluated the 24 bus system considering all the possible
combinations of the 5 transformer taps. As referred, this
network has 5 transformers each one with 11 taps. This
means that, in the scope of this enumeration process,
161051 power flow studies were run in order to analyse all
possible configurations. As conclusions, we can refer that
the best solution obtained with the enumeration process
coincides with the best one identified using the Simulated
Annealing approach. This aspect is even more remarkable
if one compares the number of power flow runs. The SA
required 312 versus 161051 for the enumeration process.
Therefore, if the computation time of one power flow
exercise is accepted as time unit, the SA approach only
