Interval arithmetic in current injection power flow analysis

Electrical Power and Energy Systems 43 (2012) 1106–1113

Contents lists available at SciVerse ScienceDirect

Electrical Power and Energy Systems journal homepage: www.elsevier.com/locate/ijepes

Interval arithmetic in current injection power flow analysis q L.E.S. Pereira 2, V.M. da Costa ⇑,1, A.L.S. Rosa 2 Department of Electrical Engineering, Federal University of Juiz de Fora, Campus Universitário, Bairro Martelos, 36036-330 Juiz de Fora, MG, Brazil

a r t i c l e

i n f o

Article history: Received 19 April 2011 Received in revised form 17 May 2012 Accepted 18 May 2012 Available online 15 July 2012 Keywords: Current injection power flow Data uncertainty Electrical power systems Interval arithmetic Krawczyk algorithm

a b s t r a c t This paper incorporates interval arithmetic into current injection method to solve power flow problem under both load and line data uncertainty. In this new methodology, the resulting interval nonlinear system of equations is solved using Krawczyk method. The proposed methodology is implemented in the Matlab environment using the Intlab toolbox. Results are compared with those obtainable by Monte Carlo simulations. IEEE test systems and a South – southeastern Brazilian network are used. Ó 2012 Elsevier Ltd. All rights reserved.

1. Introduction Power flow [1,2] is the most frequently performed study in electric power systems, and deals with the calculation of voltages and line flows, in a large sparse electrical network, for a given load and generation schedule. The conventional power flow solution comprises power equations expressed in terms of polar or rectangular voltage coordinates. Over the last years, the current injection power flow has been applied to different electric power system areas and remarkable results have been published in literature, as follows: Current-versions power flow in voltage rectangular coordinates are proposed in Refs. [3,4]. Based on [4], a sparse augmented formulation for solving a set of control devices in the current injection power flow problem is described in Refs. [5,6]. Many different flexible AC transmission system devices and controls are incorporated into the power flow problem by using this augmented methodology. Based on [4], a second order power flow methodology by using current injection equations expressed in voltage rectangular coordinates is presented in Ref. [7]. This proposed method leads to a substantially faster second order power flow solution, when com-

q This paper presents a new mathematical model regarding power flow solutions under load and line data uncertainties by using current injection equations expressed in rectangular voltage coordinates and the interval arithmetic. ⇑ Corresponding author. Tel.: +55 32 21023461; fax: +55 32 21023442. E-mail addresses: [email protected] (L.E.S. Pereira), vander.costa@ ufjf.edu.br (V.M. da Costa), [email protected] (A.L.S. Rosa). 1 Author has published other papers in the International Journal of Electrical and Energy Systems using the same kind of equations and coordinates. 2 Tel.: +55 32 32325352; fax: +55 32 21023442.

0142-0615/$ - see front matter Ó 2012 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.ijepes.2012.05.034

pared to the conventional method expressed in terms of power mismatches and written in voltage rectangular coordinates. The current injection rectangular power flow [4] is extended in [8] for solving an unbalanced distribution three-phase power flow problem. This methodology presents an expressive mathematical robustness and converges with a reduced number of iterations. A step size optimization factor to be used for solving unbalanced three phase distribution current injection is presented in Ref. [9]. On the other hand, regarding uncertainty in electric power systems, some comments become necessary. All loads are provided by measurement devices which are frequently inaccurate. This uncertainty in input data can be enlarged due to both round off error accumulation and truncating processes that occur in numerical computation. Additionally, electric power systems evolve through time, and it is reasonable to evaluate the range of all plausible conditions that might be encountered as a result of expected uncertainties [10]. Some technical publications can be pointed out as follows: The input data are subject to uncertainty. Interval analysis can be used for taking into account the effects of errors on numerical analysis. It is based on interval operations including interval arithmetic. The technique calculates the interval between the upper and lower bounds regarding variables under uncertainty. To yield solutions in a mathematical sense and to express the variable under uncertainty as an interval variable are the main advantages. On the other hand, the main limitation is sometimes to overestimate the interval between the upper and lower bounds. Thus, loads and other parameters can be characterized not by a single number but rather by a range of real values or a real interval [10–13]. In Ref. [10], the Newton operator is used to solve interval nonlinear

L.E.S. Pereira et al. / Electrical Power and Energy Systems 43 (2012) 1106–1113

equations. The Gauss–Seidel method is used to solve the interval linear equations. The input data for power flow problem can be described by random variables taking into account the probabilistic nature of loads, generation and networks. It is assumed that the parameter varies according to some probability distribution. The effects of uncertainty on the steady-state behavior of electric power systems can be assessed by a stochastic or probabilistic power flow [14–18]. Another approach reported in literature for dealing with uncertainties considers that the vague or imprecise information is represented by a fuzzy number. Loads and generations are represented using possibility distributions. As a result, power flow state and output variables have also possibility distributions [19–22]. There is a fundamental difference between probabilistic approaches and possibilistic ones including interval arithmetic. Both approaches are not competing but cooperative, and the type of information that each one provides is very different, because one is quantitative, though random, and the other is qualitative [10]. As can be noted, the current injection power flow formulation has been widely applied in several studies related to the power systems area and final results have been highly encouraging. In addition, when the input conditions are uncertain, numerous scenarios need to be analyzed and reliable solution algorithms that incorporate the effect of data uncertainties into the power flow analysis must be developed. Therefore, the main objective of this paper is to develop a new model of power flow under uncertainty by incorporating interval arithmetic into the current injection formulation. Besides, the computational performance of this new model in comparison with the interval power flow in polar coordinates is presented. The notations adopted in the paper are the conventional ones whenever possible. Matrices are shown in bold. The over scripts d and i refer to deterministic and interval quantities, respectively.

According to Ref. [10], interval addition and multiplication are associative and commutative. However, the distributive law does not always hold for interval arithmetic. The failure of the distribution law often causes overestimation. In some especial cases, the distributive law remains valid. In the Matlab toolbox Intlab [26,27] real intervals may be stored by either infimum and supremum or by midpoint and radius. Intlab enables basic interval operations to be performed on real and complex interval scalars, vectors and matrices. Interval functions such as trigonometric and exponential are also available. 3. Krawczyk method One of the most used approaches for solving a set of nonlinear equations is the Krawczyk method, which stems from Newton method. Let f be a nonlinear function such that

f ðxÞ ¼ 0

f ðyÞ ¼ f ðxÞ þ JðbÞðy  xÞ

JðbÞðy  xÞ ¼ f ðxÞ

JðXÞðX  xÞ ¼ f ðxÞ

½I  JðXÞðx  yÞ ¼ f ðxÞ þ x  y

diamðXÞ 2

ð3Þ

radðXÞ ¼

ð12Þ

I represents the identity matrix. Solving (12) for y

Kðx; XÞ ¼ x  f ðxÞ þ ½I  JðXÞðX  xÞ

ð2Þ

ð11Þ

Adding term (x  y) to both sides of Eq. (11)

Interval mathematics [23,24] considers a set of methods for handling intervals that approximate uncertain data. These methods are based on the definition of both interval arithmetic and optimal scalar product [25]. An interval number [x1, x2] is the set of real numbers x such that x1 6 x 6 x2. x1 is the infimum and x2 is the supremum. The interval X is defined by X ¼ ½x1 ; x2  ¼ f~ x2 Rjx1 6 ~ x 6 x2 g. The midpoint, the diameter and the radius of an interval X are given by

diamðXÞ ¼ x2  x1

ð10Þ

Defining the interval [x, y] e X then

Since [x, y] e Xy can be replaced with interval X

ð1Þ

ð9Þ

J represents the Jacobian matrix. Assuming f(y) = 0

y ¼ x  f ðxÞ þ ½I  JðXÞðy  xÞ

1 ðx1 þ x2 Þ 2

ð8Þ

Let y be a incremental value from x and b be a value between x and y. By applying the mean value theorem to Eq. (8)

2. Interval arithmetic

midðXÞ ¼

1107

ð13Þ

ð14Þ

K(x, X) is called Krawczyk operator, and corresponds to interval solution in Eq. (13). The insertion of both a preconditioning matrix C and iteration h into Eq. (14) yields [11]

KðxðhÞ ; X ðhÞ Þ ¼ xðhÞ  Cf ðxðhÞ Þ þ ½I  CJðX ðhÞ ÞðX ðhÞ  xðhÞ Þ

ð15Þ

X ðhþ1Þ ¼ X ðhÞ \ KðxðhÞ ; X ðhÞ Þ

ð16Þ

Eq. (16) means that the interval Krawczyk method provides the solution through the intersection of two interval sets. Besides, the main advantage of this operator is that no interval linear equations have to be solved at any iteration. C is a preconditioning matrix given by the midpoint inverse of J(X). 4. Interval current injection power flow – proposed method

Interval arithmetic operations are defined such that the interval encloses all possible real results [23] in order to guarantee the reliability of interval methods. The elementary operations, such as addition, subtraction, multiplication and division, are defined as follows

X þ Y ¼ ½x1 þ y1 ; x2 þ y2 

ð4Þ

X  Y ¼ ½x1  y2 ; x2  y1 

ð5Þ

X  Y ¼ ½minðx1 y1 ; x1 y2 ; x2 y1 ; x2 y2 Þ; maxðx1 y1 ; x1 y2 ; x2 y1 ; x2 y2 Þ

ð6Þ

  X 1 1 if ¼ ½x1 ; x2   ; Y y2 y1

ð7Þ

0 R ½y1 ; y2 

This paper proposes a new methodology to handle load and line data uncertainty in electric power systems by modeling the power flow problem through current injection equations written in rectangular voltage coordinates. Additionally, this paper proposes to calculate, in an interval manner, active and reactive generations, active and reactive power flows and losses. 4.1. Initialization of iterative process The interval current injection power flow (ICIPF) method is run after convergence of deterministic power flow. Its initialization is carried out based on deterministic voltage profile and on definition of load variations as follows:

1108

L.E.S. Pereira et al. / Electrical Power and Energy Systems 43 (2012) 1106–1113

Pidk ¼ ½Pddk  ð1  aPk Þ;

P ddk  ð1 þ aPk Þ

Q idk ¼ ½Q ddk  ð1  aQ k Þ;

Q ddk  ð1 þ aQ k Þ

ð19Þ

out the interval power flow, and it is evaluated from the current injection Jacobian matrix at deterministic power flow solution. f(x) refers to interval current mismatches vector, that is, f ðxÞ ¼ ½DIir DIim t , and it is calculated only once. The state variable vector x corresponds to the real and imaginary voltage components which stem from the deterministic power flow, that is, x ¼ ½V dr V dm t . The vector X corresponds to the interval power flow solution. The over script t denotes transposed vector. After calculating the Krawczyk operator, a new interval voltage solution is obtained by using (16). To check convergence of the proposed method, the difference between radii at iteration (h + 1) and radii at iteration (h) is calculated. If the difference is greater than a specified tolerance, denoted by s, then Krawczyk method must be employed to calculate new interval voltages. Otherwise, the iterative process is stopped.

ð20Þ

4.4. Calculation of interval output variables

ð17Þ ð18Þ

where P dk þ jQ dk is the complex load power at bus k. aPk and aQ k are factors which denote active and reactive load variations at bus k. Variations of resistances , reactances and shunt susceptances of each branch are defined by equations similar to (17) and (18), with the exception that Pd and Qd are replaced by corresponding line data. Interval voltages are initialized by using the deterministic voltage profile as midpoint and the largest load or line data variation factor as radius of interval. Thus

V irk ¼ ½V drk  ð1  amax Þ; V imk ¼ ½V dmk  ð1  amax Þ;

V drk  ð1 þ amax Þ V dmk  ð1 þ amax Þ

where V rk þ jV mk are the real and imaginary voltage components at bus k. The strategy adopted is to consider amax as the largest factor between all aPk and aQ k in order to ensure a good initial condition for convergence of iterative process. 4.2. Calculation of interval current mismatches The real and imaginary components of current injection equations written in rectangular coordinates, regarding bus k, are given respectively by [7,9]

ð21Þ

X V m P k  V rk Q k ðGki V mi þ Bki V ri Þ  k 2 V rk þ V 2mk i2/k

ð22Þ

Pik  V drk þ Q ik V dmk

DIimk ¼ Idmk 

ðV dk Þ2 Pik  V dmk  Q ik V drk ðV dk Þ2

ð27Þ

The increment regarding g can be expressed as nonlinear functions of real and imaginary voltage components at buses k and m. As a consequence, it is possible to linearize Eq. (27) by using the Taylor series around the corresponding state variables calculated through deterministic power flow program. Therefore

@g @g @g @g DV r k þ DV m k þ DV r m þ DV mm @V rk @V mk @V rm @V mm

ð28Þ

where DV rk þ jDV mk is the complex voltage variation at bus k. The interval increment of g can be expressed as follows:

where Irk þ jImk is the complex injected current at bus k; Gki þ jBki is the k  i element of bus admittance matrix; /k denotes the set of buses directly connected to bus k, including itself; Pk + jQk is the net complex injected power at bus k. The real and imaginary components of interval current mismatches derived from (21) and (22) are given by

DIirk ¼ Idrk 

g ¼ gðV rk ; V mk ; V rm ; V mm Þ

Dg ¼

X V r P k þ V mk Q k Irk ¼ ðGki V ri  Bki V mi Þ  k 2 V rk þ V 2mk i2/k I mk ¼

Let g denote any output power flow variable such as power flow and losses, and let k  m be the branch under analysis. Thus

ð23Þ

ð24Þ

Dg i ¼

Pik ¼ Pigk  Pidk

ð25Þ

Q ik ¼ Q igk  Q idk

ð26Þ

Vk is the voltage magnitude at bus k; DIrk þ jDImk is the complex current mismatch at bus k; Pgk þ jQ gk is the generated complex power at bus k. Interval powers are constant during the interval power flow problem. Therefore, current mismatches are calculated only once. 4.3. Iterative process Eq. (15) refers to the Krawczyk method. For the sake of clarity, some comments about this equation need to be presented. The interval Jacobian matrix is calculated using interval voltages and its structure is shown in Ref. [7]. The matrix C is constant through-

ð29Þ

Active and reactive power generations are calculated in the same way. Let k be the bus generation under analysis. The function g depends on the real and imaginary voltage components in all buses directly connected to k, including itself. Therefore, the number of partial derivatives presented in Eq. (29) depends on the number of adjacent buses. The interval voltages calculated at the end of the iterative process may be substituted in Eq. (29). However, new interval operations need to be carried out and this procedure may lead to a large and inaccurate diameter regarding Dig . In order to improve the accuracy of intervals, Eq. (29) should be written in terms of interval current mismatches which are calculated from the input data of power flow program. Therefore

2

where

@g @g @g @g DV irk þ DV imk þ DV irm þ DV imm @V rk @V mk @V rm @V mm

3

.. .

2 3 .. 7 6 . 7 6 i 6 6 DV r 7 6 7 7 k 7 6 X 7 6 7 6 7 6 DV irm 7 6 7" 6 # Y 7 6 7 DI i 6 . 7 6 6 m 6 . 7 ¼ 6 .. 7 6 . 7 6 . 7 7 DI i 7 6 r 7 6 DV i 7 6 6 Z 7 6 mk 7 6 7 7 6 7 6 DV i 7 6 4W 5 6 mm 7 .. 5 4 .. . .

ð30Þ

where X, Y, Z and W are the rows of inverse current injection Jacobian matrix evaluated after convergence of deterministic power flow program. By substituting Eq. (30) in (29)

Dg i ¼



@g @g @g @g Xþ Zþ Yþ W @V rk @V mk @V rm @V mm

"

DIim DIir

# ð31Þ

1109

L.E.S. Pereira et al. / Electrical Power and Energy Systems 43 (2012) 1106–1113

The term between brackets is solved through simple algebraic operations. The corresponding interval of g is computed as follows:

g i ¼ g d þ Dg i :

ð32Þ

Table 1 Voltage magnitudes – IEEE 14 bus. Bus

Method

Lower voltage (pu)

Upper voltage (pu)

Deterministic voltage (pu)

4

ICIPF MCS

1.01724 1.01718

1.01827 1.01816

1.01775

9

ICIPF MCS

1.05536 1.05514

1.05666 1.05673

1.05601

10

ICIPF MCS

1.05044 1.05021

1.05172 1.05176

1.05108

13

ICIPF MCS

1.04988 1.04996

1.05107 1.05079

1.05048

14

ICIPF MCS

1.03505 1.03457

1.03635 1.03644

1.03570

4.5. Solution methodology The proposed algorithm can be summarized in the following steps: Step 1: Run a deterministic power flow program. Step 2: Calculate both load and line data variations using (17) and (18). Step 3: Initialize the interval voltage profile using (19) and (20). Step 4: Calculate the current mismatches using (23) and (24). Step 5: Apply the Krawczyk operator according to Eq. (15). Step 6: Update the interval voltages vector using Eqs. (33) and (34).

V rðhþ1Þ ¼ V rðhÞ \ Kðx; XÞ ðhþ1Þ ðhÞ ¼ Vm \ Kðx; XÞ Vm

ð33Þ ð34Þ

Step 7: Check convergence through s. The value of s adopted in this paper is 104. If convergence is not reached, then go back to Step 5. Otherwise go to Step 8. Step 8: Calculate power flows, losses and power generations. 5. Results

Table 2 Phase angles – IEEE 14 bus. Bus

Method

Lower angle (°)

Upper angle (°)

Deterministic angle (°)

4

ICIPF MCS

10.32885 10.51332

10.21161 10.11254

10.27021

9

ICIPF MCS

14.94749 15.21633

14.82033 14.66093

14.88390

10

ICIPF MCS

15.10555 15.37550

14.97857 14.81613

15.04204

13

ICIPF MCS

15.16411 15.44432

15.03773 14.87507

15.10091

14

ICIPF MCS

16.03751 16.33681

15.91053 15.73436

15.97401

5.1. Initial considerations In order to perform this study some simulations were accomplished by using IEEE 14 and 300 bus systems and a practical Brazilian network composed of 1768 buses, 2527 branches and 119 generation buses. For the three test systems, the total amount of buses, branches and generation buses is 2082, 2958 and 171, respectively. The tolerance adopted for convergence of the iterative process, related to both deterministic and interval power flow methods, is 104 pu. The Monte Carlo simulation (MCS) method validates the proposed methodology. One million of Monte Carlo simulations were performed for IEEE test systems, and three hundred thousand for 1768 bus system. For each simulation, different values of power injections and line data, within the intervals previously defined, are selected and conventional power flow solutions are performed. Interval solutions are obtained by monitoring the largest and the smallest values of both state and output variables calculated during all simulations. For the lack of space, results are displayed for only buses and branches that present the five largest relative errors (e) in comparison with MCS method. These errors are expressed in terms of absolute values. The radius of an interval associated with any input and output variable is defined around its respective deterministic value. In this paper, the radius of 5% is considered for series resistances and reactances and shunt susceptances of all IEEE 14 bus branches. Radii of 3% and 2% are considered for all active and reactive loads of 300 bus and 1768 bus systems. In addition, radii of 3% and 2% are considered for only 20% of the branches of IEEE 300 bus and 1768 bus, respectively. These branches were randomly chosen. The radius of 1% is assumed for all active power generations.

Table 3 Voltage magnitudes – IEEE 300 bus. Bus

Method

Lower voltage (pu)

Upper voltage (pu)

Deterministic voltage (pu)

17

ICIPF MCS

1.03297 1.06133

1.09686 1.06835

1.06492

120

ICIPF MCS

0.92968 0.95378

0.98719 0.96300

0.95844

139

ICIPF MCS

0.98136 1.01065

1.04206 1.01279

1.01171

192

ICIPF MCS

0.90933 0.93249

0.96558 0.94131

0.93746

234

ICIPF MCS

1.00754 1.03765

1.06987 1.03987

1.03871

Table 4 Phase angles – IEEE 300 bus. Bus

Method

Lower angle (°)

Upper angle (°)

Deterministic angle (°)

17

ICIPF MCS

13.81965 13.52391

12.30700 11.87730

13.04295

120

ICIPF MCS

9.23880 10.20123

8.20891 7.46745

8.70912

139

ICIPF MCS

3.72098 3.96595

3.30109 3.14304

3.50478

5.2. Calculation of state variables

192

ICIPF MCS

11.59565 11.95755

10.31392 10.57828

10.93694

Tables 1–4 present the state variables yielded by ICIPF and MCS methods for IEEE 14 and 300 bus systems. For example, the ICIPF

234

ICIPF MCS

21.81043 21.27168

19.54066 19.52708

20.64997

1110

L.E.S. Pereira et al. / Electrical Power and Energy Systems 43 (2012) 1106–1113

Fig. 1. Voltage magnitude error – IEEE 14.

method calculates that the voltage magnitude at bus 10, IEEE 14 bus, is included in the interval [1.05044, 1.05172] pu. Therefore, the voltage magnitude at bus 10 is not smaller than 1.05044 pu and not larger than 1.05172 pu. On the other hand, the phase angle at bus 10 is not smaller than 15.10555° and not larger than 14.97857°. The same reasoning can be extended to all output variables presented throughout the paper. The largest possible voltage variation around its deterministic value is 2.99% at bus 139 (IEEE 300), and the smallest possible voltage variation is 0.05% at bus 8 (IEEE 14). Regarding phase angles, the largest and the smallest variations are 6.08% at bus 120 (IEEE 300) and 0.19% at bus 11 (IEEE 14), respectively. By comparing with MCS, Figs. 1–4 depict the largest relative errors of voltage magnitudes and phase angles for some buses including those ones presented in Tables 1–4. NSV denotes the amount of each state variable under consideration. The proposed formulation calculates 4164 voltage magnitudes and the same amount of phase angles taking into account the lower and upper bounds. It can be observed that 98.82% of the voltage magnitudes yielded by ICIPF method present an error less than 1% and only 1.18% present an error greater than 1%. Regarding phase angles, 96.47% present an error less than 1% and only 0.47% an error greater than 5%.

5.3. Calculation of output variables

Fig. 2. Phase angle error – IEEE 14.

Tables 5 and 6 present the active and reactive losses yielded by ICIPF and MCS methods for IEEE 14 bus. The largest possible active loss variation around its deterministic value is 15.86% at branch 13–14 (IEEE 14), and the smallest possible loss variation is 0.53% at branch 6–11 (IEEE 14). Regarding reactive losses, the largest

Table 5 Active losses – IEEE 14 bus. Branch

Method

Lower active loss (MW)

Upper active loss (MW)

Deterministic active loss (MW)

4–5

ICIPF MCS

0.48756 0.49086

0.54246 0.53876

0.51230

1–5

ICIPF MCS

2.64202 2.64390

2.88102 2.88714

2.76284

9–10

ICIPF MCS

0.01097 0.01138

0.01478 0.01437

0.01283

9–14

ICIPF MCS

0.10455 0.10717

0.12782 0.12533

0.11563

13–14

ICIPF MCS

0.04590 0.04791

0.06232 0.06040

0.05379

Fig. 3. Voltage magnitude error – IEEE 300.

Table 6 Reactive losses – IEEE 14 bus.

Fig. 4. Phase angle error – IEEE 300.

Branch

Method

Lower reactive loss (MVAr)

Upper reactive loss (MVAr)

Deterministic reactive loss (MVAr)

4–5

ICIPF MCS

1.53793 1.54833

1.71107 1.69942

1.61594

7–9

ICIPF MCS

0.75321 0.76604

0.85208 0.83898

0.80001

9–10

ICIPF MCS

0.02915 0.03024

0.03925 0.03817

0.03408

9–14

ICIPF MCS

0.22240 0.22797

0.27189 0.26659

0.24597

13–14

ICIPF MCS

0.09345 0.09755

0.12688 0.12299

0.10953

1111

L.E.S. Pereira et al. / Electrical Power and Energy Systems 43 (2012) 1106–1113 Table 8 Reactive power flow – IEEE 300 bus. Deterministic reactive power flow (MVAr)

71.67347 71.81591

195.43647 194.07225

132.23509

ICIPF

120.57547

5.67565

62.45981

MCS

120.33121

5.15433

ICIPF

5.57547

18.87626

Method

15–14

ICIPF MCS

120– 116 172– 139 225– 192

Fig. 5. Active loss error – IEEE 14.

Upper reactive power flow (MVAr)

Branch

234– 228

Lower reactive power flow (MVAr)

MCS

5.81241

17.48281

ICIPF

380.75735

589.13778

MCS

381.33703

587.10320

ICIPF

261.55547

150.87546

MCS

260.06179

153.95327

7.21049

88.08206

57.57514

Tables 7 and 8 present the active and reactive power flows yielded by ICIPF and MCS methods for IEEE 300 bus. The largest possible active power flow variation around its deterministic value is 295.44% at branch 225–192 (IEEE 300), and the smallest possible active power flow variation is 10.35% at branch 90–92 (IEEE 300). Regarding reactive losses, the largest and the smallest variations are 568.85% at branch 225–192 (IEEE 300) and 16.28% at branch 86–102 (IEEE 300), respectively. By comparing with MCS, Figs. 7 and 8 depict the largest relative errors of power flows for some branches including those ones presented in Tables 7 and 8. The proposed formulation calculates 11,832 active power flows and the same amount of reactive power flows taking into account the lower and upper bounds. It can be

Fig. 6. Reactive loss error – IEEE 14.

and the smallest variations are 15.20% at branch 9–10 (IEEE 14) and 1.62% at branch 10–11 (IEEE 14), respectively. By comparing with MCS, Figs. 5 and 6 depict the largest relative error of losses for some branches including those ones presented in Tables 5 and 6. NOV denotes the amount of output variable under consideration. The proposed formulation calculates 5916 active losses and the same amount of reactive losses taking into account the lower and upper bounds. It can be observed that 98.93% of the active losses yielded by ICIPF method present an error less than 1% and only 0.22% an error greater than 5%. Regarding reactive losses, 98.02% present an error less than 1% and 0.46% an error greater than 5%.

Table 7 Active power flow – IEEE 300 bus.

Fig. 7. Active power flow error – IEEE 300. Deterministic active power flow (MW)

Branch

Method

Lower active power flow (MW)

Upper active power flow (MW)

15–14

ICIPF MCS

239.76552 238.18281

8.35758 8.04065

123.39609

120– 116

ICIPF

22.23347

144.65247

65.69656

MCS

22.52047

143.59711

172– 139

ICIPF

21.65324

15.76347

MCS

21.66069

15.51364

225– 192

ICIPF

36.86554

19.32554

MCS

38.78419

19.98552

ICIPF

55.43678

52.13547

MCS

56.25572

51.82317

234– 228

19.46229

9.88830

2.55217 Fig. 8. Reactive power flow error– IEEE 300.

1112

L.E.S. Pereira et al. / Electrical Power and Energy Systems 43 (2012) 1106–1113

Table 9 Reactive power generation – 1768 bus. Bus

Method

Table 11 Computation times relationships – 1768 bus. Upper reactive power generation (MVAr)

Lower reactive power generation (MVAr)

Deterministic reactive power generation (MVAr)

10

ICIPF MCS

803.11232 802.96909

566.24378 567.50238

143.9106

12

ICIPF MCS

1596.57575 1597.05173

667.65247 667.21511

433.0160

14

ICIPF MCS

84.13478 84.29319

130.03576 129.55984

25.3061

16

ICIPF MCS

1069.34268 1069.05045

226.13537 225.03091

422.6598

18

ICIPF MCS

507.75434 503.58183

1007.87645 1008.67623

222.4069

Tasks

IPPF

ICIPF

C f(x) J(X) K(x, X) Total time

1.14947 1.04351 2.83300 1 2.80120

1 1 1 1 1

that 96.82% of the generations yielded by ICIPF method present an error less than 1% and only 1.53% an error greater than 5%. 5.4. Computational Performance

Fig. 9. Reactive power generation error – 1768 bus.

observed that 98.62% of the active power flows yielded by ICIPF method present an error less than 1% and only 0.19% an error greater than 5%. Regarding reactive power flow, 98.20% present an error less than 1% and 0.49% present an error greater than 5%. The 1768 bus system is used to calculate the active and reactive power generations at slack bus. The interval active generation yielded by ICIPF method is [451.252, 470.487] MW and by MCS is [451.304, 469.704] MW. On the other hand, the reactive generations yielded by both methods are [1615.244, 1607.195] MVAr and [1616.652, 1608.332] MVAr, respectively. Therefore, a good agreement between results can be observed. The largest relative error of 0.166% refers to the upper bound of active power generation. Table 9 presents the reactive power generations yielded by ICIPF and MCS methods. The largest possible reactive generation variation around its deterministic value is 458.06% at bus 10 (1768 bus), and the smallest possible reactive generation variation is 113.09% at bus 55 (1768 bus). By comparing with MCS, Fig. 9 depicts the largest relative errors for some generation buses including those ones presented in Table 9. The proposed formulation calculates 342 power generation taking into account the lower and upper bounds. It can be observed

Table 10 Computation times – 1768 bus.

The ICIPF method achieves convergence with 3, 2 and 2 iterations for IEEE 14, IEEE 300 and 1768 bus, respectively. Table 10 displays computation times per iteration, in seconds, required by ICIPF to calculate the power flow solution regarding the 1768 bus system. In addition, there is a column for interval polar power flow (IPPF) with the purpose of comparing the performance of both methodologies. Computation times were obtained when using a personal computer AMD Athlon II X4 630 processor and 4 GB RAM. Table 11 displays the relationships between computation times required by both voltage coordinates for each one of the tasks related to interval power flow calculation. The main tasks considered are the assembling of preconditioning matrix C and the calculation of f(x), J(X) and K(x, X). ICIPF times are taken as reference. According to Table 10, the largest time is associated with the assembling of both interval Jacobian matrices. The times for the remaining tasks can be neglected. The assembling of interval current injection Jacobian matrix is faster than the polar version because the majority of its elements are constant, whereas the polar Jacobian matrix has all elements calculated through trigonometric functions. The ICIPF method presents a computational gain around 180% in comparison with the IPPF approach. Therefore, the use of current injection formulation requires a smaller computation time and, therefore, leads to a substantially faster power flow analysis subjected to uncertainty. 6. Conclusion This paper presents a new methodology for handling uncertainty in electrical power systems by using the interval arithmetic incorporated into the current injection power flow formulation. If input data vary within relatively small ranges, interval arithmetic yields good results that include all possible solutions. This methodology is very simple and reliable, and converges with a few iterations. Convergence is not an issue even in large systems. In general, ICIPF method not only presents good results in comparison with Monte Carlo simulations, but it also propitiates a remarkable computational gain in comparison with polar power flow coordinates. Since electric system data are uncertain, this paper demonstrates that the proposed method can be regarded as a powerful tool in power flow analysis under uncertainty and it can be indeed recommended for general use. References

Tasks

IPPF

ICIPF

C f(x) J(X) K(x, X) Total time

10.46428 1.01703 1790.26936 0.00001 1802.83546

9.10356 0.97462 631.93411 0.00001 643.59364

[1] Kulworawanichpong T. Simplified Newton–Raphson power flow solution method. Int J Electr Power Energy Syst 2010;32(6):551–8. [2] Mallick S, Rajan DV, Thakur SS, Acharjee P, Ghoshal SP. Development of a new algorithm for power flow analysis. Int J Electr Power Energy Syst 2011;33(8):1479–88. [3] Exposito AG, Ramos ER. Augmented rectangular load flow model. IEEE Trans Power Syst 2002;17(2):271–6.

L.E.S. Pereira et al. / Electrical Power and Energy Systems 43 (2012) 1106–1113 [4] Da Costa VM, Martins N, Pereira JLR. Developments in the Newton–Raphson power flow formulation based on current injections. IEEE Trans Power Syst 1999;14(4):1320–6. [5] Da Costa VM, Pereira JLR, Martins N. An augmented Newton–Raphson power flow formulation based on current injections. Int J Electr Power Energy Syst 2001;23(4):305–12. [6] Variz AM, Da Costa VM, Pereira JLR, Martins N. Improved representation of control adjustments into the Newton–Raphson power flow. Int J Electr Power Energy Syst 2003;25(7):501–13. [7] Ferreira CA, Da Costa VM. A second order power flow based on current injection equations. Int J Electr Power Energy Syst 2005;27(2):254–63. [8] Garcia PAN, Pereira JLR, Da Costa VM, Martins N. Three phase power flow calculations using the current injection method. IEEE Trans Power Syst 2000;15(2):508–14. [9] Da Costa VM, Oliveira ML, Guedes MR. Developments in the analysis of unbalanced three-phase power flow solutions. Int J Electr Power Energy Syst 2007;29(2):501–13. [10] Wang Z, Alvarado FL. Interval arithmetic in power flow analysis. IEEE Trans Power Syst 1992;7(3):1341–9. [11] Barboza L, Dimuro G, Reiser R. Towards interval analysis of the load uncertainty in power electric systems. In: 8th international conference on probabilistic methods applied to power systems. Ames, Iowa; 2004. [12] Vaccaro A, Canizares CA, Villacci D. An affine arithmetic-based methodology for reliable power flow analysis in the presence of data uncertainty. IEEE Trans Power Syst 2010;25(2):624–32. [13] Zhang P, Li W, Wang S. Reliability-oriented distribution network reconfiguration considering uncertainties of data by interval analysis. Int J Electr Power Energy Syst 2012;34(1):138–44. [14] Dimitrovski A, Tomsovic K. Boundary load flow solutions. IEEE Trans Power Syst 2004;19(1):348–55.

1113

[15] Stefopoulos GK, Meliopoulos AP, Cokkinides GJ. Probabilistic power flow with non-conforming electric loads. Int J Electr Power Energy Syst 2005;27(4):627–34. [16] Su CL. Probabilistic load-flow computation using point estimate method. IEEE Trans Power Syst 2005;20(4):1843–51. [17] Hu Z, Wang X. A probabilistic load flow method considering branch outages. IEEE Trans Power Syst 2006;21(2):507–14. [18] Usaola J. Probabilistic load flow with wind production uncertainty using cumulants and Cornish–Fisher expansion. Int J Electr Power Energy Syst 2009;31(3):474–81. [19] Ramaswamy M, Nayar KR. On-line estimation of bus voltages based on fuzzy logic. Int J Electr Power Energy Syst 2004;26(1):681–4. [20] Bijwe PR, Hanmandlu M, Pande VN. Fuzzy power flow solutions with reactive limits and multiple uncertainties. Electr Power Syst Res 2005;76(2):145–52. [21] Cortés-Carmona M, Palma-Behnke R, Jiménez-Estévez G. Fuzzy arithmetic for the DC load flow. IEEE Trans Power Syst 2010;25(1):206–14. [22] Kalesar BM, Seifi AR. Fuzzy load flow in balanced and unbalanced radial distribution systems incorporating composite load model. Int J Electr Power Energy Syst 2010;32(1):17–23. [23] Moore RE. Methods and applications of interval analysis. Philadelphia: SIAM; 1979. [24] Keafort RB, Kreinovich V. Applications of interval computations. Boston: Kluwer; 1996. [25] Kulish UW. Advanced arithmetic for the digital computer design of arithmetic units. Electron Notes Theor Comput Sci 1999;24. [26] Rump SM. Intlab-interval laboratory in developments in reliable computing. In: Csendes T. editor. Dordrecht. Kluwer; 1999. p. 77–104. [27] Hargreaves GI. Interval analysis in Matlab, Manchester Centre for Comp Math, Manchester, Numerical Analysis Report N 416; 2002.