arXiv:1601.04150v1 [] 16 Jan 2016

Selection of Most Effective Control Variables for
Solving Optimal Power Flow Using Sensitivity
Analysis in Particle Swarm Algorithm
Mohamed Abuella

Constantine J. Hatziadoniu

Departmtent of Electrical and Computer Engineering
University of North Carolina at Charlotte
Charlotte, USA
Email: mabuella@uncc.edu

Department of Electrical and Computer Engineering
Southern Illinois University
Carbondale, USA
Email: hatz@siu.edu

Abstract—Solving the optimal power flow problem is one of the
main objectives in electrical power systems analysis and design.
The modern optimization algorithms such as the evolutionary
algorithms are also adopted to solve this problem, especially when
the intermittency nature of generation resources are included,
as in wind and solar energy resources, where the models are
stochastic and non-linear. This paper uses the particle swarm
optimization algorithm for solving the optimal power flow for
IEEE-30 bus system. In addition to selection of the most effective
control variables based on sensitivity analysis to alleviate the
violations and return the system back to its normal state.
This adopted strategy would decrease the optimal power flow
calculation burden by particle swarm optimization algorithm,
especially with large systems.
Keywords— Optimal power flow; Particle swarm; Sensitivity
analysis.

I. I NTRODUCTION
Optimal power flow (OPF) is a major tool that has been
extensively researched since it has introduced. It deals with
the minimum cost of power production in electrical power systems with certain constraints [1], [2]. The studies about OPF
methods can be traced back to the 1960s when Carpentier and
Siroux firstly discussed the OPF problem, and then Domme
and Tinney presented a simplified derivative algorithm which
is the first practicable algorithm for OPF problem [3]. But in
this algorithm the vibration phenomenon is appearing closing
to the optimum. Since then, various kinds of mathematical
programming approaches, based on linear and nonlinear programming were proposed in succession, including Newton
method, quadratic programming, and interior-point method
[3].
The aforementioned methods utilize first or second derivative in essence. By this way, it is apt to fall into local optima.
Furthermore, there is a difficulty of applying derivative-based
optimization techniques to solve OPF for systems that including variable generations. Therefore, various non-classical optimization methods have emerged to cope with some of these
shortcomings [4]. The main modern optimization techniques
are genetic algorithm (GA), evolutionary programming (EP),
artificial neural network (ANN), simulated annealing (SA), ant

colony optimization (ACO), and particle swarm optimization
(PSO). They have been successfully applied to wide range
of optimization problems in which global solutions are more
preferred than local ones [5], [6].
Kennedy and Eberhart introduced PSO as a new heuristic
method [7]. Different PSO applications in power systems
are covered in [6]. Abido introduced PSO to solve the OPF
problem [8]. References [9], [10] have the mechanism of
implementing PSO for solving OPF.
The sensitivity analysis of OPF is applied to deduce the
sensitivity matrices of the voltage and the current states. The
fundamentals of sensitivity analysis of OPF can be found
with details in [11], [12], and [13]. The sensitivity analysis
is applied in power systems for voltage stability studies as
a fast indication of the voltage collapse and a safe voltage
margin [14], [15], [16], and [17].
The rest of the paper is organized as follows. Section II gives
a description of the optimal power flow problem. Section III
introduces the mathematical derivation of the sensitivity analysis. Section IV particle swarm optimization (PSO) algorithm
is introduced with its parameters and its mechanism is also
explained. In Section V implementation of PSO and the
sensitivity analysis to find the OPF solution for several study
cases are presented and the discussion of the results is included
as well. Finally, Section VI gives the conclusions.
II. O PTIMAL P OWER F LOW
The optimal power flow is an optimization problem to find
the optimal allocation of output power among the available
generators with given constraints without violations of the
operation and the security constraints of the system. The optimal allocation depends on various factors, such as operating
cost, system security (or risk) and CO2 emissions, in general
they are called cost factors. The objective function of the
optimization problem in this paper is to minimize the operating
cost of real power generation, while keeping the electrical
system in normal and secure operation condition.
The objective function and the constraints are mostly nonlinear and many methods and algorithms have been developed

on the basis of cost factors; generation source type, conventional or renewable; uncertainty treatment (i.e. deterministic or
stochastic). For instance, Lagrangian relaxation, direct search
method, evolution programming, particle swarm optimization,
genetic algorithms, and simulated annealing are some of the
solution methods for the optimal power flow problem [18].
In mathematical terms, the optimal power flow (OPF) is an
optimization problem set up to minimize an objective function
subject to equality and inequality constraints.
The equality constraints are the power flow equations, while
the inequality constraints are the limits on control variables
and the operating limits of power system dependent (state)
variables. The control variables include the generator real
power, the generator bus voltages, the transformer tap settings,
and the reactive power of switchable VAR sources. On the
other hand, the dependent or state variables include the load
bus voltages, the generator reactive powers, and the line flows.
Generally, the OPF problem is a large-scale highly constrained
nonlinear non-convex optimization problem.
Mathematically, the OPF problem can be formulated as
follow:
M in J(x,u)

u is the vector of independent (control) variables consisting
of generator voltages VG , generator real power outputs PG
except the slack bus PG1 , transformer tap settings T, and shunt
VAR compensations QC .
Hence, u can be expressed as:
uT = [VG1 ...VGN G , PG2 ...PGN G , T1 ...TN T , QC1 ...QCN C ]
(7)
Where indices NT and NC are the number of the regulating
transformers and shunt compensators respectively.
• Generating constraints:
Generator voltages, real power outputs, and reactive
power outputs are restricted by their lower and upper
limits as follows:

•

h(x,u) ≤ 0

(3)

Where J is the objective function to be minimized, it could
be the cost of real power of thermal units, wind-powered
units, or mix of them. g is the equality constraints represent
the power flow equations. h is the inequality constraints that
represent the operating limits of the system.
Where g(x,u) = 0 are the balanced power flow equations as
following:
Pi − Vi

n
X

Vj Yij cos(δi − δj − θij ) = 0

Qi − Vi

Vj Yij sin(δi − δj − θij ) = 0

,
≤ PGi ≤ PGmax
PGmin
i
i

i = 1, ..., N G

(9)

(5)

j=1

Where Pi is the specified real power at bus i, and it equals to
the difference between the generation and demand real power
(PGi -PDi ) at bus i, and the same for Qi . Yij is the admittance
between buses i and j. Vi is the voltage magnitude of bus i
and δi is the phase angle of the voltage at bus i.
In equations (1-3), x is the vector of dependent (state)
variables consisting of slack bus power PG1 , load bus voltages
VL , generator reactive power output QG , and transmission line
ratings (loadings) Sline . Here x can be expressed as:
xT = [PG1 , VL1 ...VLN L , QG1 ...QGN G , Sline1 ...Slinenl ] (6)
Where the indices NL, NG, and nl are number of load
buses, number of generators, and number of transmission lines
respectively.

i = 1, ..., N G

Timin ≤ Ti ≤ Timax ,
•

max
Qmin
Ci ≤ QCi ≤ QCi ,
•

(10)

i = 1, ..., N T

(11)

Shunt VAR constraints:
Shunt VAR compensations are restricted by their limits
as follows:
i = 1, ..., N C

(12)

Security constraints:
These include the constraint of voltages at load buses and
the transmission line loadings as follows:
,
≤ VLi ≤ VLmax
VLmin
i
i
max
,
Slinei ≤ Sline
i

(4)

j=1

n
X

(8)

Transformer constraints:
Transformer tap settings are bounded as follows:

Subject to:
(2)

i = 1, ..., N G

max
Qmin
Gi ≤ QGi ≤ QGi ,

(1)

g(x,u) = 0

,
≤ VGi ≤ VGmax
VGmin
i
i

i = 1, ..., N L

(13)

i = 1, ..., nl

(14)

It is worth to mention that the control variables are selfconstrained. The hard inequalities of PG1 , VL , QG , and Sline
can be incorporated in the objective function as quadratic
penalty terms (penalty functions). Therefore, the objective
function in equation (1) can be augmented as follows:
)2 + λV
Jaug = J + λP (PG1 − PGlim
1

NL
X

)2
(VLi − VLlim
i

i=1

+λQ

NG
X
i=1

2
(QGi − Qlim
Gi ) + λS

nl
X

max 2
)
(Slinei − Sline
i

i=1

(15)

Where λP , λV , λQ , and λS are penalty factors and xlim is
the limit value of the dependent variable x given as:
 max
x
;
(x > xmax )
xlim =
(16)
min
x
;
(x < xmin )

III. S ENSITIVITY A NALYSIS FOR O PTIMAL P OWER F LOW
Earlier research on the application of sensitivity analysis in
power system belongs to Peschon et al [11]. They introduced
two methods. First one can be applicable to normal power
flow problems for small changes in the variables such as active
generation, and the second method considers the minimization
of objective function satisfying some constraints such as power
flow equation. Similaresearch was carried out by Gribik et al
[12].

where p1 , p2 ,...,p2N are the elements of p.
When changes are small, solution for ∆x will be:

A. Mathematical Formulation

If p variables are not changed then (24) can be rewritten as:

The method of calculating the sensitivities of voltages and
currents are determined simultaneously, which are further used
to determine the changes in power flows [13]. Considering the
generalized equations of the form:
g(x,u,p) = 0

(17)

where g is 2N dimensional vector, and N is number of buses.
The variables mentioned in equation (17) can be categorized
as:
(x) are dependent (state) variables, these are the controlled
variables and they are unknown. x is a 2N dimensional vector.
(u) are independent control variables, these are the operating
variables or imposed variables of the system. u is an M
dimensional vector.
(p) are parameter variables, these are uncontrollable variables
and are normally specified in the power flow problem such as
the admittance and the loads.
Depending on the variables to be determined, the variables
in the power flow problem can be selected as x, u, and p. One
might be interested in controlling M variables out of the 2N
variables.
If x0 , u0 , and p0 are the initial state vectors, rewriting equation
(17) as:
g(x0 , u0 , p0 ) = 0
(18)
The changes ∆x corresponding to small changes ∆u and ∆p,
will satisfy the new equations:
g(x0 + ∆x, u0 + ∆u, p0 + ∆p) = 0

(19)

Expanding (19) by Taylor’s series and neglecting higher order
terms,
g(x0 +∆x, u0 +∆u, p0 +∆p) = g(x0 , u0 , p0 )+gx ∆x+gu ∆u+gp ∆p
(20)

where, gx , gu and gp are the partial derivatives of g with
respect to x, u and p respectively and are given by:
gx =

∂(g1 , g2 , ..., g2N )
∂(x1 , x2 , ..., x2N )
∂(g1 , g2 , ..., g2N )
∂(u1 , u2 , ..., uM )

(22)

∂(g1 , g2 , ..., g2N )
∂(p1 , p2 , ..., p2N )

Su = −gx −1 gu

(25)

Sp = −gx −1 gp

(26)

∆x = Su ∆u

(27)

The set of dependent and independent variables can be chosen as the system requirements and the problem formulation.
Some of the parameters of a type may belong to the set of
dependent whereas remaining parameters of same type may
belong to the set of independent variables. for instance, as bus
voltages they might be considered as independent variables
when they are at generator buses while they are considered
dependent at load buses.
B. Determination of Voltage Sensitivities at Buses
Power flow equations are comprising of 6 variables namely
P, Q, V, δ ,Y and θ. All the variables can be assumed to
be obtained or specified at the base condition. The variables
Y and θ are normally specified and are constant. The other
variables are not always constant and they are either specified
or determined, depending upon the type of buses. The variables
for which changes are specified are grouped as independent
variables and the variables which are determined against these
changes are grouped as dependent variables [13].
For the slack bus, V and δ are specified and P and Q are
subjected to change. For generator bus, P and V are specified
and Q and δ are subjected to change. For load buses, P and Q
are specified and V and δ are changed. Now consider the power
system of N buses and B branches. Power flow equations can
be described by (4, 5). There are 2N set of equations and a
set of 2N variables can be selected as state variables (x) and
remaining as control variables (u).
Consider that only M control variables are changed and for
these changes, it is desired to obtain the changes in the real
and reactive power at slack buses, reactive power and angles at
generator buses and voltages and angles at load buses. Then,
the power flow equations can be written as following:
(28)

Let
Psl , PG , PL ∈ Pi
Qsl , QG , QL ∈ Qi
Vsl , VG , VL ∈ Vi
Grouping the variables of (28) as:

where u1 , u2 ,...,uM are the elements of u.
gp =

(24)

where Su and Sp are the sensitivities of x with respect to u
and p respectively and are obtained as:

g(Vi , Vj , δi , δj , Pi , Qi , Yij , θij ) = 0
(21)

where x1 , x2 ,...,x2N are the elements of x.
gu =

∆x = Su ∆u + Sp ∆p

(23)

x = [Psl , Qsl , QG , δG , VL , δL ]

(29)

u = [Vsl , δsl , PG , VG , PL , QL ]

(30)

p = Yij , θij

(31)

From (27), the changes in dependent variables can be obtained
[∆Psl , ∆Qsl , ∆QG , ∆δG , ∆VL , ∆δL ] =
S[∆Vsl , ∆δsl , ∆PG , ∆VG , ∆PL , ∆QL ]

(32)

where S is the sensitivity matrix of order 2Nx2N and can be
obtained as given by (25). For slack bus and generator buses
following substitution can be made in (32):
∆Vsl = ∆VG = ∆δsl = 0

(33)

After determining the changes in the load bus voltages, load
bus angles and generator bus angles from (32) and with the
substitutions from (33) all the bus voltages and angles can be
arranged as:
[∆V, ∆δ] = [∆Vsl , ∆VG , ∆VL , ∆δsl , ∆δG , ∆δL ]

(34)

C. Determination of Current Sensitivities in the Lines
It is well known that the changes in voltage angles and
voltage magnitudes will cause changes in branch currents [13].
These currents in complex form can be expressed as:
Iij = Yij [Vi (cosδi + jsinδi ) − Vj (cosδj + jsinδj )]

(35)

Where Yij = |Yij |∠θij and Iij ∈ B, since B is the number of
branches. Equation (35) can be written in the form:
gij (Iij , Yij , θij , Vi , Vj , δi , δj ) = 0

(36)

Grouping the variables of (36) as
x = Iij
u = Vi , Vj , δi , δj (i.e. V and δ variables at all buses)
p = |Yij |, θij .
Sensitivities of Iij for the changes in Vi , Vj , δi , δj can be
obtained from (27) as:
∆Iij = R[∆V, ∆δ]

(37)

where R is sensitivity matrix obtained by (25) which is given
as:
−1
R = −gijx
giju
(38)
With gijx is Jacobian of gij with respect to x (i.e. Iij ).
While giju is Jacobian of gij with respect to u (i.e.
Vi , Vj , δi , δj ) Substituting from (34) and (38), the equation
(37) can be rewritten as:

A. Standard Algorithm
PSO, as an optimization tool, provides a swarm-based
search procedure in which particles change their positions with
time. In a PSO system, particles fly around in a multidimensional search space. During flight, each particle adjusts its
position according to its own experience, and the experience of
neighboring particles, making use of the best position encountered by itself and its neighbors. When improved positions are
being discovered these will then come to guide the movements
of the swarm. The process is repeated and by doing so it
is hoped, but not guaranteed, that a satisfactory solution will
eventually be discovered [10].
The following is the conventional terminology of the parameters in PSO: Let x and v denote a particle coordinates
(position) and its corresponding flight speed (velocity) in
a search space, respectively. Therefore, the ith particle is
represented as xi = [xi1 , xi2 , ...., xim ]. Since m is the last
dimension or coordinate of the position of the the ith particle
in the search space and so that d = 1, 2, ..., m.
The best previous position of the ith particle is recorded and
represented as,
pbesti = [pbesti1 , pbesti2 , ...., pbestim ].
The position of the best particle among all the particles in the
group is represented by the gbest. In a particular dimension d
there is a group best position which is gbestd.
The velocity for the ith particle is represented as,
vi = [vi1 , vi2 , ...., vid ]. The modified velocity and position
of each particle can be calculated by using the following
formulas:
k+1
k
vid
= w∗vid
+c1 ∗U ∗(pbestkid −xkid )+c2 ∗U ∗(gbestkd −xkid )
(40)
k+1
k
xk+1
=
x
+
v
(41)
id
id
id

i = 1, 2, ...., n;
d = 1, 2, ..., m
Where
n
number of particles in a group;
m
number of members in a particle;
k
pointer of iterations (generations);
w
inertia weight factor;
c1 , c2 acceleration factors;
U
uniform random number in the range [0,1];
k
xkid ,vid
the position and velocity of the ith particle in the dth
dimension at iteration k;
The search mechanism of the PSO using the modified
velocity and position of individual based on (40) and (42)
is illustrated in Fig. 1.

[∆Iij ]B×1 = [R]B×2N [∆Vsl , ∆VG , ∆VL ∆δsl , ∆δG , ∆δL ]2N×1
(39)

Where ∆Vsl = ∆VG = ∆δsl = 0.
IV. PARTICLE S WARM O PTIMIZATION A LGORITHM
The original PSO suggested by Kennedy and Eberhart is
based on the analogy of swarm of bird and school of fish
[7]. The algorithm was simplified and it was observed to be
performing a solution to an optimization problem.

Fig. 1.

PSO search mechanism

In the above procedures, the velocity should between
vdmin ≤ vid ≤ vdmax If vdmax is too high, particles might
fly past good solutions. If vdmax is too small, particles may
not explore sufficiently beyond local solutions. In many experiences with PSO, it is often set at 10 - 20% of the dynamic
range of the variable on each dimension [10].
The constants c1 and c2 represent the weighting of the
stochastic acceleration terms that pull each particle toward the
pbest and gbest positions. Low values allow particles to move
far from the target regions before being dragged back. On the
other hand, high values result in sudden movement toward, or
past, target regions. Hence, the acceleration constants c1 and
c2 are often set to be 2 according to empirical experience [10].
Suitable selection of inertia weight w in (40) provides
a balance between global and local explorations, to find a
sufficiently optimal solution. As originally developed w, often
decreases linearly from about 0.9 to 0.4 during the run. In
general, the inertia weight is set according to the following
equation:
w = wmax −

(wmax − wmin )
× iter
itermax

(42)

Where itermax is the maximum number of iterations (generations), and iter is the current number of iterations.
Start

Initialize particles with
random position and velocity
vectors

For each particle’s position
(x) evaluate objective
function (J)

If J(x) better than
J(pbest)
then pbest = x

Set best of pbest as gbest

Update
vi = ω vi + c1 U (pbesti-xi) + c2U (gbest-xi)
x i = x i + vi

o

Is Iter.=maxiter?
Is adequate objective
function reached?
Yes
Stop: giving gbest, optimal
solution

M in J(x,u)

(43)

g(x,u) = 0

(44)

h(x,u) ≤ 0

(45)

u∈U

(46)

x = [PG1 , VL , QG , Sline ]

(47)

u = [PG , VG , T, QC ]

(48)

subject to:

Where:

The equality constraint in (44) are the nonlinear power flow
equations as in (4 and 5)
The inequality constraints (45) are the functional operating
constraints, such as transmission line limits, load bus voltage
magnitude limits, generator reactive capabilities, and slack
bus active power output limits. Constraints (46) define the
feasibility region of the control variables of the problem such
as the active power output limits of the generators (except the
generator at the slack bus), generation bus voltage magnitude
limits, transformer-tap setting limits, and bus shunt admittance
limits.
Each particle in PSO is a vector containing the control
variables u, suggesting a possible solution to the OPF problem.
Then the position of the ith particle xi can be represented
as xi = ui = (ui1 , ui2 , ..., uim ), where m is the number of
dimensions and it is also represented the number of control
max
variables, d ∈ [1, m], uid ∈ [umin
id , uid ].
min
max
uid , and uid are the lower and upper bounds of uid . The
particles are moving in an m dimensional space.
For consistency’s sake, the general definition of the swarm
particle is which used in the rest of the paper, as in equations
(40 and 42). Therefore, the ith particle is represented as
xi = [xi1 , xi2 , ...., xim ] instead of ui .
Each particle attempts to minimize the following objective
function:
#
" NS
NG
X
X
µi ∗ hi (x,u)
(49)
Fi (PGi ) + λ
Jaug =
µi =

Fig. 2.

PSO algorithm flowchart

B. Implementation of PSO for OPF Problem
The primary objective of the OPF problem is to determine
the optimal combination of power outputs of all generating
units so as to meet the required load demand at minimum
operating cost while satisfying system equality and inequality
constraints.
The formulation of OPF for applying PSO is done by
separating the problem variables to state variables, x, and
control variables, u, as it was described in equations (1, 2,
and 3), they also mentioned here as follows:

i=1

i=1

Since:



1;
0;

hi (x,u) > 0
hi (x,u) ≤ 0

(50)

Here the objective function becomes unconstrained or augmented objective function by using the classical penalty functions principle. All inequality constraints in equation (45)
replaced by penalty terms. While the power balance equations
(44), which are the equality constraints, is solved for each
particle and in every iteration by Newton-Raphson power flow
algorithm, therefore no need to use a penalty function for this
equality constraint in equation (49).
Jarg is the penalized objective function and Fi (PGi ) is
the cost of the real power from the generator PGi while
λ is the penalty factor for operating constraints. µi is an
indicator of occurring any violations and work outside the

feasibility region of the solution. It has only two values as
in equation (50), it is either 1 when a violation of the limits
occurred in the corresponding constraint or 0 when there is no
violation. The penalty terms that have been used are quadratic
penalty functions as those in equation (15). Whereas N G is
the number or generators while N S is the number of the state
variables to be bounded within their limits. The penalty factor
λ is used to penalize the cost proportional to the amount of
constraint violations, the suitable value of the penalty factor
is chosen after some runs of the algorithm [19]. According to
the equations (40 and 42) in every iteration each particle of
the swarm updates its position coordinates (dimensions) until
the termination condition of the algorithm is met.
V. S TUDY C ASES

AND

S IMULATION R ESULTS

For analysis and investigation aim, the PSO is applied to find
the optimal power flow for the standard IEEE 30-bus system,
besides the implementation of the sensitivity analysis to find
the most effective control variables for solving the OPF by
PSO but in this case within a reduced space of dimensions. The
simulation is performed by MATLAB software, the flowcharts
of the main parts of simulation codes can be found in [20] .
A. The Data of The System
The system is shown in Fig. 3 and the data of the buses,
lines, and generators are given in Appendix (A). It consists
of six conventional thermal generators at buses 1, 2, 5, 8, 11,
and 13, and 41 branches, four of them are transformers with
off-nominal tap ratios in branches 6-9, 6-10, 4-12, and 28-27.
In addition, the buses 10, 12, 15, 17, 20, 21, 23, 24, and 29 are
equipped with shunt VAR compensators. The limits of control
variables are indicated in table (I). For the other operating
(state) variables such as voltages at load buses, the limits are
[0.95-1.1]. The limits of reactive power of generators QG and
the transmission lines ratings are both indicated in generator
data and line data tables respectively in Appendix (A).
TABLE I
C ONTROL VARIABLES AND THEIR LIMITS
The control variable

PG1
(MW)

PG2
(MW)

PG3
(MW)

PG4
(MW)

PG5
(MW)

PG6
(MW)

VGs

The upper limit
The lower limit

200
50

80
20

50
15

35
10

30
10

40
12

0.95
1.1

(pu)

Ts

Qshs
(MVAR)

0.9
1.1

0
5

Ts: transformer tap ratio; VGs: voltage at generator bus; Qsh: VAR compensation.

B. The PSO Algorithm and Its Parameters for Solving OPF
The skeleton of PSO algorithm is taken form these references [21], [22] after it has been modified for solving OPF.
Initially several runs are carried out besides to the helpful
information from [8], [9], [10] to select the suitable parameters
for PSO algorithm. The inertia weight is decreasing linearly
with iterations from its initial value at 0.9 to its ultimate value
at 0.4, while the acceleration factors C1 and C2 are equal to 2
and the number of particles is 10. The termination condition
occurs when the 5 significant digits after the decimal point of
the optimal solution have not changed for last 50 iterations,
the algorithm will then consider this as an optimal solution. In

Fig. 3.

Single-line diagram of IEEE 30-bus test system [8]

addition, the maximum number of iterations after which the
algorithm also terminates is 500.
C. The Objective Function
The objective function is to minimize the operating cost of
the system as in equation (49) which is used in PSO algorithm
and its constraints are represented by (50).
The penalty functions are quadratic penalty functions as
those in equation (15) and its constraints as (16).
D. Study of Base Case
The running power flow of the initial operating point, which
is the base case of loading and it is given in bus data table
in Appendix (A), yielding violations in the lower limit of
voltage at load buses 18, 19, 20, 21, 22, 23, 24, 25, 26,
27, 29, 30. There is no violation for bus 28 although it is
also far from generators, that because this bus is fed by two
branches and one of them directly from generator at bus 8.
Furthermore, there is a rating violation of the transmission
line which connects the buses 1 and 2.
Let’s first apply PSO algorithm with all of control variables
(i.e. PGs , VGs , Ts , and Qshs ) for solving the optimal power
flow (OPF) of this base case of loading (283.4M W ) . For
sure, in this case all the violations can be easily alleviated
and the voltage at load buses and the transmission lines rating
within their limits as they are shown in Fig. 4. The cost of the
real output power of generators is minimized to 798.43 $/hr
as it is illustrated in table (II).
By adjusting the voltages at both ends of the transmission line
within their limits, the rating constraint of the transmission line
1-2 is alleviated as well for this branch.

TABLE II
G ENERATORS OUTPUTS OF BASE CASE (283.4 M W )

combination is done in PSO for adjusting the violations
in operating constraints of the power system by using less
number of dimensions. When control variables increase the

PG1
(MW)

PG2
(MW)

PG3
(MW)

PG4
(MW)

PG5
(MW)

PG6
(MW)

Losses
(MW)

Cost
($/hr)

176.94

48.71

21.27

21.09

11.83

12.00

8.4382

798.43

TABLE III
R ANKED CONTROL VARIABLES FOR VIOLATED VOLTAGES OF LOAD BUSES
AT B ASE C ASE

1.15

Upper limit

1.1

P.U.

1.05

1

Lower limit 0.95

0.9

VL3

VL4

VL6

VL7

VL9 VL10 VL12 VL14 VL15 VL16 VL17 VL18 VL19 VL20 VL21 VL22 VL23 VL24 VL25 VL26 VL27 VL28 VL29 VL30

Load Bus Voltages

Fig. 4.



Voltage levels at load buses for base case of IEEE 30-bus system

VL19
V8
V2
V1
V5
V11
T9
T27
V13
T12
Qsh23
Qsh29
T10
Qsh24
Qsh12
Qsh15
Qsh20
Qsh17
Qsh21
Qsh10

VL20
T27
V5
V8
V2
V1
V11
T9
V13
T12
T10
Qsh29
Qsh24
Qsh23
Qsh17
Qsh21
Qsh10
Qsh12
Qsh20
Qsh15

VL21
V2
V8
V1
V5
T27
V11
T9
V13
T12
Qsh29
Qsh23
T10
Qsh15
Qsh12
Qsh20
Qsh17
Qsh24
Qsh10
Qsh21

VL22
T27
V1
V2
V5
V8
V13
T12
V11
T9
Qsh29
Qsh23
T10
Qsh15
Qsh12
Qsh24
Qsh21
Qsh20
Qsh17
Qsh10

VL23
V1
V2
V8
V5
T27
Qsh29
V11
T9
V13
T12
T10
Qsh20
Qsh21
Qsh17
Qsh12
Qsh10
Qsh24
Qsh15
Qsh23

VL24
V8
V5
V2
V1
V11
T9
T27
T10
V13
T12
Qsh15
Qsh17
Qsh20
Qsh29
Qsh12
Qsh23
Qsh10
Qsh21
Qsh24

VL25
V1
V5
V2
V13
T12
V8
V11
T9
T10
Qsh20
Qsh15
T27
Qsh12
Qsh17
Qsh23
Qsh29
Qsh24
Qsh21
Qsh10

VL26
V8
V2
V5
V1
V13
T12
T10
V11
T9
Qsh15
Qsh23
T27
Qsh12
Qsh17
Qsh20
Qsh24
Qsh21
Qsh29
Qsh10

VL27
T27
V13
T12
Qsh23
Qsh29
Qsh24
V1
V2
V8
V5
T10
Qsh15
V11
T9
Qsh12
Qsh20
Qsh17
Qsh21
Qsh10

VL29
V8
V1
V2
V5
V13
T12
Qsh23
V11
T9
T10
Qsh15
T27
Qsh24
Qsh12
Qsh17
Qsh20
Qsh29
Qsh21
Qsh10

VL30
V8
V1
V2
V5
V13
T12
V11
T9
T10
Qsh15
T27
Qsh29
Qsh12
Qsh23
Qsh24
Qsh20
Qsh17
Qsh21
Qsh10



E. Implementation of Sensitivity Analysis With PSO
Firstly, the sensitivity analysis of optimal power flow is
conducted for the base case to find the sensitivity matrices of
the voltage and the current and to select the most effective
control variables to adjust the violations in the OPF at the
base case within a less-dimensions space of PSO algorithm.
The mathematical approaches of calculating the sensitivity
matrices of voltage Su and current R are discussed in (III-B
and III-C), the resulting rank of state and control variables,
in addition to the visualized sparsity pattern of elements
in matrices Su and R are shown in Fig. 5a and Fig. 5b
respectively.
The sensitivity analysis of voltages

Ěϭϰ
Ěϭϱ
Ěϭϲ
Ěϭϳ
Ěϭϴ
Ěϭϵ
ĚϮϬ
ĚϮϭ
ĚϮϮ
ĚϮϯ
ĚϮϰ
ĚϮϱ
ĚϮϲ
ĚϮϳ
ĚϮϴ
ĚϮϵ
ĚϯϬ
sϯ
sϰ





The sensitivity analysis of Currents

Su=-Jx-1Ju

∆x 

WŐϭ
YŐϭ
YŐϮ
YŐϱ
YŐϴ
YŐϭϭ
YŐϭϯ
ĚϮ
Ěϯ
Ěϰ
Ěϱ
Ěϲ
Ěϳ
Ěϴ
Ěϵ
ĚϭϬ
Ěϭϭ
ĚϭϮ
Ěϭϯ



sϲ
sϳ
sϭϰ
sϭϱ
sϭϲ
sϭϳ
sϭϴ
sϭϵ
sϮϬ
sϮϭ
sϮϮ
sϮϯ
sϮϰ
sϮϱ
sϮϲ
sϮϴ
sϮϵ
sϯϬ

∆u
WŐϮ
WŐϱ
WŐϴ
WŐϭϭ
WŐϭϯ
sŐϭ
sŐϮ
sŐϱ
sŐϴ
sŐϭϭ
sŐϭϯ
dϵ
dϭϬ
dϭϮ
dϮϳ
YƐŚϭϬ
YƐŚϭϮ
YƐŚϭϱ
YƐŚϭϳ
YƐŚϮϬ
YƐŚϮϭ
YƐŚϮϯ
YƐŚϮϰ
YƐŚϮϵ

0

10

20

30

40

50

0

5

10
15
nz = 1344

20

25

R=-Jx-1Ju

∆x 
I16_17
I15_18
I18_19
I19_20
I10_20
I10_17
I10_21
I10_22
I21_22
I15_23
I22_24
I23_24
I24_25
I25_26
I25_27
I28_27
I27_29
I27_30
I29_30
I8_28
I6_28

I1_2
I1_3
I2_4
I3_4
I2_5
I2_6
I4_6
I5_7
I6_7
I6_8
I6_9
I6_10
I9_11
I9_10
I4_12
I12_13
I12_14
I12_15
I12_16
I14_15

∆u

0
5
10
15
20
25
30
35
40
0

10

20

30
nz = 164

40

50

60

V1
.
.
.
.
.
.
.
V30
d1
.
.
.
.
.
.
.
.
d30



Where d stands for phase-angle of bus voltage (į).


(a) Voltage sensitivity matrix Su

Where I1-2 stands for the power flow in the branch (1-2).


(b) Current sensitivity matrix R

Fig. 5. Sensitivity matrices and their state ∆x and control ∆u variables
at base case 283.4M W

dimensions of the particles in PSO algorithm also increase and
that can cause a complexity for manipulating the dimensions of
the particles [9]. Thereby, finding and using the most effective
control variables to adjust and correct the violations can
decrease the dimensions in PSO and enhance its performance.

F. PSO Solution for Base Case Using Most Effective Control
Variables
Now let’s using only the most effective control variables
to adjust the violations in initial operating point of IEEE
30-bus system. Several combinations of dominant control
variables (as they are ranked in table (III)) can be chosen,
some of them are sufficient combination to bring back the
violations in the voltage at load buses within their limits for
base case of 283.4 M W loading level. The PSO’s results
with the combinations of most effective control variables are
listed in table (IV). As it is shown in table (IV), except
TABLE IV
PSO RESULT OF COMBINATIONS OF MOST EFFECTIVE CONTROL
VARIABLES AT BASE CASE

Table (III) presents the most dominant control variables for
every violated bus voltage are ranked from top to down at base
case.
Notice that the real power of generators Pgs are not listed
in Table (III) because they are also used in the OPF solution
search by PSO algorithm to find the minimum cost, although
Pgs are effectively contributing in adjusting the violations
since the constraints are augmented with the main objective
function of OPF as in (49).
This order of control variables are achieved by sensitivity
analysis between the voltage sensitivity matrix Su and its
corresponding column vectors of state and control variables.
Once the control variables that have most effect on state
variables are determined and ranked from the most powerful
to the less. Next to that, the selection of these control variable

Case
283.4MW
283.4MW
283.4MW
283.4MW
283.4MW
283.4MW
283.4MW

Violation
VL18…30, Sline1
VL18…30, Sline1
VL18…30, Sline1
VL18…30, Sline1
VL18…30, Sline1
VL18…30, Sline1
VL18…30, Sline1

*

Control variables
(**Pgs ), (Vg1, Vg5, Vg8)
(Pgs ), (Vg1, Vg2, Vg8)
(Pgs ), (Vg1, Vg2, Vg5, Vg8, Vg13)
(Pgs ), (Vg1, Vg2, Vg8), (T27)
(Pgs ), (Vg1, Vg2, Vg5, Vg8)
(Pgs ), (Vgs)
All (i.e. Pgs, Vgs, Ts, Qshs)

Still exist violations
VL26, VL30, Sline10
VL30, Sline10
VL30
-----------------

Cost
818.92
805.51
800.17
802.41
800.22
799.86
798.43

* VL18…30 stands for violations at buses 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 29, 30.
** Pgs stands for all real power of generators except the first generator at slack bus.

the case of using all control variables, it’s obvious that the
control variables combination of Pgs and Vgs is the most
effective with minimized cost 799.86 $/hr. While for other
less number of control variables, the combination of most
effective control variables of Pgs , Vg1 , Vg2 , Vg5 and Vg8 and
the other of Pgs , Vg1 , Vg2 , Vg8 and T27 succeed to adjust all
violations, but the former combination produces lower cost
800.22 $/hr.

G. PSO Solution for Different Loading Levels
For IEEE 30-bus test system with other cases of loading,
higher and lower than the base case, PSO is used with only
the most effective control variables to adjust the violations, if
they exist. The results are shown in table (V).
ϮϱDt
dŚĞŐĞŶĞƌĂƚŽƌŽƵƚƉƵƚƐƐĞƚĂƚƚŚĞŝƌůŽǁĞƌůŝŵŝƚƐсϭϭϳDt
TABLE V
PSO RESULT FOR SEVERAL LOADING CASES OF IEEE 30- BUS TEST

SYSTEM
Total loading

Violation

Control variables

Cost

125MW

----

(*Pgs )

309.060

150MW

VL30

(Pgs ), Vg1, Vg2, Vg8

374.577

200MW

VL26, VL29, VL30

(Pgs ), Vg1, Vg2, Vg5, Vg8

517.827

250MW

----

(Pgs )

681.602

275MW

Sline1

(Pgs ), Vg1, Vg2

769.947

300MW

Sline1

(Pgs ), Vg1, Vg2

861.65

325MW

VL26, VL30, Sline1

(Pgs ), Vg1, Vg2, Vg5, Vg8, Vg13

956.922

VL26, VL29, VL30, Sline1

(Pgs ), Vg1, Vg2, Vg5, Vg8, Vg11,
Vg13, T9, T10, T12, Qsh15

1058.605

350MW

effort is needed for PSO algorithm and enhancing its performance, especially for large systems including many stochastic
generation resources.
The following could be included for further work: PSO
algorithm needs some work on selecting proper parameters
and it also needs more accurate mathematical description
for its convergence. PSO can be applied in wind power
bid marketing between electric power operators. In addition
to operating cost, the environment effects and security or
risk of wind power penetration can be included by using
multi-objective models. Using singular value decomposition
and pseudo-inverse techniques could be considered for further
study to find the effective control variables.
A PPENDIX
A. The Data for IEEE 30-Bus Test System
The data for IEEE 30-bus test system as following [23]:

* Pgs stands for Pg2, Pg5, Pg8, Pg11, and Pg13.



Note: Every OPF solution of the previous loading case
is considered
an initial operating point for the following
Ϯϴϯ͘ϰDt asŽƐƚсϴϬϬ
loading case. for example, if the OPF of the loading case 125
M
W has been solved then the OPF of the successive loading
EŽŐƵĂƌĂŶƚĞĞŽĨǁŚŝĐŚƐĞůĞĐƚŝŽŶƚŚĂƚǁŽƵůĚǁŽƌŬĂƐŝŶĐĂƐĞϯϬϬDtƚŚĞĨŽůůŽǁŝŶŐĐŽŵďŝŶĂƚŝŽŶŚĂǀĞ
ďĞĞŶƐĞůĞĐƚĞĚ;
Vg1, W
Vg2, Vg5,
Vg8Ϳ΀ďƵƚƚŚĞǀŝŽůĂƚŝŽŶƐǁĞƌĞŶŽƚĐŽƌƌĞĐƚĞĚ΁;
Vg2, Vg5, Vg8,
Vg13) Žƌ
case
150 M
considers
the solution of theVg1,
previous
loading
Vg1, Vg2, Vg5, Vg8,Vg13, T12ͿŽƌ;Vg1, Vg2, Vg5, Vg8, Vg11, Vg13, T12)
case
125 M W as its initial point and so on.
VzŽƵŬŶŽǁǁŚŝĐŚĞůĞŵĞŶƚŚĂƐƚŚĞŵŽƐƚĞĨĨĞĐƚďƵƚƵĐĂŶŶŽƚŬŶŽǁǁŚŝĐŚĐŽŵďŝŶĂƚŝŽŶƚŚĂƚǁŝůůǁŽƌŬ͊
L26 , VL29 , and VL30 are the weakest buses in the system that
are susceptible to violations more than other buses. While the
transmission line that connect bus 1 and bus 2 is the weakest
transmission line and it suffers from violation of its rating
for several loading cases. In the last case 350 M W a variety
of control variables are needed to keep the system in secure
operation, but they are still less than using all control variables.

Bus No.

Bus Code*

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30

1
2
0
0
2
0
0
2
0
0
2
0
2
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0

DATA OF

TABLE VI
IEEE 30- BUS TEST SYSTEM

Voltage

Angle

PL

QL

PG#

Q G#

Q G_low

Q G_high

Q sh

(pu)

(pu)

(MW)

(MVAR)

(MW)

(MVAR)

(MVAR)

(MVAR)

(MVAR)

1.05
1.04
1
1
1.01
1
1
1.01
1
1
1.05
1
1.05
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1

0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0

0
21.7
2.4
7.6
94.2
0
22.8
30
0
5.8
0
11.2
0
6.2
8.2
3.5
9
3.2
9.5
2.2
17.5
0
3.2
8.7
0
3.5
0
0
2.4
10.6

0
12.7
1.2
1.6
19
0
10.9
30
0
2
0
7.5
0
1.6
2.5
1.8
5.8
0.9
3.4
0.7
11.2
0
1.6
6.7
0
2.3
0
0
0.9
1.9

0
80
0
0
50
0
0
20
0
0
20
0
20
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0

0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0

-20
-20
0
0
-15
0
0
-15
0
0
-10
0
-15
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0

250
100
0
0
80
0
0
60
0
0
50
0
60
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0

0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0

* Bus code: 1 for slack bus; 2 for PV bus; 0 for Load bus.
#These generators' Data are an initial guess.

376.2
376

961

375.8

960

375.6

959

375.4

Cost

Cost

B US

375.2

TABLE VII
G ENERATORS DATA OF IEEE 30- BUS TEST SYSTEM

958
957

375
956
374.8
955

374.6
374.4

Gen. No.

954
0

100

200

300

400

500

20

40

Iterations

80

100

120

140

160

Iterations

(a) 150 M W
Fig. 6.

60

(b) 325 M W

Cost vs. iterations of two loading cases

1
2
5
8
11
13

a

b

($/MW^2.Hr)

($/MW.Hr)

0.00375
0.0175
0.0625
0.00834
0.025
0.025

2
1.75
1
3.25
3
3

c
0
0
0
0
0
0

PG_low

PG_high

(MW)

(MW)

50
20
15
10
10
12

200
80
50
35
30
40

VI. C ONCLUSION
The implementation of PSO algorithm to solve the OPF
problem is useful and worth of investigation. Moreover, PSO
algorithm is easy to apply and simple since it has fewer
number of parameters to deal with, comparing to other modern
optimization algorithms. In addition, PSO algorithm is appropriate for solving the optimal power flow for systems that
include variable generation resources. Using most effective
control variables by applying sensitivity analysis reducing the
space dimensions of PSO and hence improving the computing

R EFERENCES
[1] A. J. Wood and B. F. Wollenberg, Power Generation, Operation and
Control. New York: Wiley, 2nd ed., 1996.
[2] X.-F. Wang, Y.-H. Song, and M. Irving, Modern Power Systems Analysis.
Springer, 2008.
[3] C.-R. Wang, H.-J. Yuan, Z.-Q. Huang, J.-W. Zhang, and C.-J. Sun,
“A modified particle swarm optimization algorithm and its application
in optimal power flow problem,” in Proc. Int Machine Learning and
Cybernetics Conf, vol. 5, pp. 2885–2889, 2005.

TABLE VIII
L INE DATA OF IEEE 30- BUS TEST SYSTEM
Branch o.

From Bus

To Bus

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41

1
1
2
3
2
2
4
5
6
6
6
6
9
9
4
12
12
12
12
14
16
15
18
19
10
10
10
10
21
15
22
23
24
25
25
28
27
27
29
8
6

2
3
4
4
5
6
6
7
7
8
9
10
11
10
12
13
14
15
16
15
17
18
19
20
20
17
21
22
22
23
24
24
25
26
27
27
29
30
30
28
28

B

Transformer

Ratings

(pu)

(pu)

(pu)

(turnns ratio)

(MVA)

0.0192
0.0452
0.057
0.0132
0.0472
0.0581
0.0119
0.046
0.0267
0.012
0
0
0
0
0
0
0.1231
0.0662
0.0945
0.221
0.0824
0.1073
0.0639
0.034
0.0936
0.0324
0.0348
0.0727
0.0116
0.1
0.115
0.132
0.1885
0.2544
0.1093
0
0.2198
0.3202
0.2399
0.0636
0.0169

0.0575
0.1852
0.1737
0.0379
0.1983
0.1763
0.0414
0.116
0.082
0.042
0.208
0.556
0.208
0.11
0.256
0.14
0.2559
0.1304
0.1987
0.1997
0.1923
0.2185
0.1292
0.068
0.209
0.0845
0.0749
0.1499
0.0236
0.202
0.179
0.27
0.3292
0.38
0.2087
0.396
0.4153
0.6027
0.4533
0.2
0.0599

0.0264
0.0204
0.0184
0.0042
0.0209
0.0187
0.0045
0.0102
0.0085
0.0045
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0.0214
0.0065

1
1
1
1
1
1
1
1
1
1
1.078
1.069
1
1
1.032
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1.068
1
1
1
1
1

130
130
65
130
130
65
90
70
130
32
65
32
65
65
65
65
32
32
32
16
16
16
16
32
32
32
32
32
32
16
16
16
16
16
16
65
16
16
16
32
32

R

X

[4] J. Hetzer, D. C. Yu, and K. Bhattarai, “An economic dispatch model
incorporating wind power,” vol. 23, no. 2, pp. 603–611, 2008.
[5] J. G. Vlachogiannis and K. Y. Lee, “Economic load dispatch—a comparative study on heuristic optimization techniques with an improved
coordinated aggregation-based pso,” vol. 24, no. 2, pp. 991–1001, 2009.
[6] M. R. AlRashidi, M. F. AlHajri, and M. E. El-Hawary, “Enhanced particle swarm optimization approach for solving the non-convex optimal
power flow,” World Academy of Science, Engineering and Technology,
vol. 62, pp. 651–655, 2010.
[7] J. Kennedy and R. Eberhart, “Particle swarm optimization,” in Proc.
Conf. IEEE Int Neural Networks, vol. 4, pp. 1942–1948, 1995.
[8] M. A. Abido, “Optimal power flow using particle swarm optimization,”
International Journal of Electrical Power Energy Systems, vol. 24, no. 7,
pp. 563–571, 2002.
[9] P. N. Biskas, N. P. Ziogos, A. Tellidou, C. E. Zoumas, A. G. Bakirtzis,
V. Petridis, and A. Tsakoumis, “Comparison of two metaheuristics with
mathematical programming methods for the solution of opf,” in Proc.
13th Int Intelligent Systems Application to Power Systems Conf, 2005.
[10] Z.-L. Gaing, “Particle swarm optimization to solving the economic
dispatch considering the generator constraints,” vol. 18, no. 3, pp. 1187–
1195, 2003.
[11] J. Peschon, D. S. Piercy, W. F. Tinney, and O. J. Tveit, “Sensitivity in
power systems,” no. 8, pp. 1687–1696, 1968.
[12] P. R. Gribik, D. Shirmohammadi, S. Hao, and C. L. Thomas, “Optimal
power flow sensitivity analysis,” vol. 5, no. 3, pp. 969–976, 1990.
[13] N. D. Ghawghawe and K. L. Thakre, “Application of power flow
sensitivity analysis and ptdf for determination of atc,” in Proc. Int. Conf.
Power Electronics, Drives and Energy Systems PEDES ’06, pp. 1–7,
2006.
[14] Y. Tamura, H. Mori, and S. Iwamoto, “Relationship between voltage
instability and multiple load flow solutions in electric power systems,”
Power Apparatus and Systems, IEEE Transactions on, no. 5, pp. 1115–
1125, 1983.
[15] C. L. DeMarco and T. J. Overbye, “An energy based security measure
for assessing vulnerability to voltage collapse,” Power Systems, IEEE
Transactions on, vol. 5, no. 2, pp. 419–427, 1990.

[16] P. Löf, T. Smed, G. Andersson, and D. Hill, “Fast calculation of a voltage
stability index,” Power Systems, IEEE Transactions on, vol. 7, no. 1,
pp. 54–64, 1992.
[17] B. Gao, G. Morison, and P. Kundur, “Voltage stability evaluation using
modal analysis,” Power Systems, IEEE Transactions on, vol. 7, no. 4,
pp. 1529–1542, 1992.
[18] X. Wei, Y. Liu, and D. Zhao, “Review of dynamic dispatch research
considering intermittent power generation,” in Proc. Asia-Pacific Power
and Energy Engineering Conf. (APPEEC), pp. 1–5, 2011.
[19] A. I. Selvakumar and K. Thanushkodi, “A new particle swarm optimization solution to nonconvex economic dispatch problems,” vol. 22, no. 1,
pp. 42–51, 2007.
[20] M. A. Abuella, “Study of particle swarm for optimal power flow in ieee
benchmark systems including wind power generators,” 2012.
[21] W.-N. Lee, Y.-W. Jeong, J.-B. Park, J.-R. Shin, and K. Y. Lee, “Development of an educational simulator for particle swarm optimization and
economic dispatch applications,” in Proc. Int. Conf. Intelligent Systems
Applications to Power Systems ISAP 2007Lee2007, pp. 1–6, 2007.
[22] B. Birge, “Psot - a particle swarm optimization toolbox for use with
matlab,” in Swarm Intelligence Symposium, 2003. SIS ’03. Proceedings
of the 2003 IEEE, pp. 182 – 186, april 2003.
[23] O. Alsac and B. Stott, “Optimal load flow with steady-state security,”
Power Apparatus and Systems, IEEE Transactions on, vol. PAS-93,
pp. 745 –751, may 1974.
BIOGRAPHIES
Mohamed A. Abuella (IEEE student member). He is a PhD student at
University of North Carolina at Charlotte. He received his Bachelor’s degree
from College of industrial Technology, Misurata, Libya in 2008, and MSc
degree from Southern Illinois University Carbondale in 2012. His research
interests are in planning and operations of power systems.
C.J. Hatziadoniu (IEEE M’87, SM’06). He is a professor of electrical and
computer engineering at Southern Illinois University Carbondale. His research
interests include power systems control and protection and application of
power electronics to power systems.

