Review on Energy Resilience

PhysicaA388 (2009)4259–4266

Contents lists available at ScienceDirect

PhysicaA

journal homepage: www.elsevier.com/locate/physa

An improvedmodel for structural vulnerability analysis of powernetworks GuoChena,c,∗, ZhaoYangDongb, David J.Hillc, GuoHuaZhangc aSchool of InformationTechnologyandElectrical Engineering, TheUniversity ofQueensland,QLD4072, Australia bDepartment of Electrical Engineering, TheHongKongPolytechnicUniversity,HungHom,Kowloon,HongKong cResearch School of InformationSciences andEngineering, TheAustralianNationalUniversity, ACT0200, Australia

a r t i c l e i n f o

Article history: Received18November2008 Received in revised form31March2009 Availableonline25 June2009

Keywords: Complexnetworks Powernetworks Admittancematrix Networkefficiency

a b s t r a c t

Electric power networks have been studied as a typical example of real-world complex networks. Traditional models for structural vulnerability analysis appear to be all based on physical topological structure. In this paper, we depict a typical power network as a weighted graph based on electrical topology by introducing its bus admittance matrix, whichembodies the importantcharacteristicsofpowernetworks inamuchmorerealistic structure.Furthermore, thenumericalsimulationforboththetraditionaldynamicalmodel and the proposed electrical topological model are investigated based on the IEEE 300 bus system respectively. The comparison demonstrates that the improvedmodel ismore precise andhighlyefficient for theanalysis of structural vulnerabilityofpowernetworks.

© 2009ElsevierB.V.All rights reserved.

1. Introduction

Power systems play an indispensable role in modern society. However, the frequency of large scale blackouts all over theworldhasnotdecreased, in spiteof technological progress andhuge investments in systemreliability andsecurity. For the past decade, many countries have suffered from serious blackouts. For instance, in August 1996, more than 4 million people in severalwesternstatesof theUSAwereoutofpower service [1]. InAugust2003, ahistoricblackoutwas triggered inthepowergridof theUnitedStatesandCanada,whichdisconnected61,800MWofpowertoanareaspanningmostof the north-eastern states ofUSA and twoprovinces of Canada, totally containingmore than50million people [2,3]. This event astonishedthewholeworldandevenremindedmanypeopleof9/11.Besides, in thesummerandautumnof theyear2003, several largeblackouts in theworldhappened, suchasLondonblackout inUK,Sweden–DemarkblackoutandItalyblackout etc [3]. The series of blackouts exposed existing problems of the current dynamic security assessment (DSA) and monitoring

systems inpreventingoutages,which stimulates researchers to seek solutions fromalternativemeans. Recently, advances instatisticalphysics,modellingandcomputationalmethodshaveattractedthe interestof thescientificcommunitytostudy powergridsas complexnetworks. Initially, complexnetworks theoryand its structural vulnerabilityanalysismethodology were proposed by physicists and they mainly focused on complex abstract networks, such as Erdös–Rényi (ER) random networks, Barabási–Albert (BA) scale-free networks and so on.Much literature has beenpublished [4–12]. Albert et al. [4] demonstrated that scale-free networks have the robust-yet fragile property, which means that they are robust against random failures of nodes but fragile to intentional attacks. Latora et al. [5] proposed the concept of network efficiency to characterize small-world networks, which is a quantity of how efficiently it exchanges information. Crucitti et al. [6]

∗ Correspondingauthor. E-mail address:[email protected] (G.Chen).

0378-4371/$– see frontmatter © 2009ElsevierB.V.All rights reserved. doi:10.1016/j.physa.2009.06.041

4260 G.Chenet al. / PhysicaA388 (2009)4259–4266

discussed network efficiency on scale-free networks and modelled the cascading failures of complex networks [7]. Zhao et al. investigated cascade defence and control in scale-free networks [11]. Sun et al. analysed statistical properties of the evolving networks and the responses of these networks under random errors and intentional attacks [12]. Furthermore, thosephysicists tried to find the linksbetween thestructural vulnerabilityof thoseabstractnetworksandpowernetworks because mathematically, power grids can be described as a graph of nodes connected by edges [13–18]. Following their previouswork,Motteretal. [13]discussedcascade-basedattacksonrealcomplexnetworksandpointedoutthattheInternet and power grids were vulnerable to important node attacks but evolved to be quite resistant to random failure of nodes. Casals et al. [14] analysed topological vulnerability of Europeanpower grid and found that power grids displaypatterns of reactiontonode losssimilar to thoseobserved inscale-freenetworks.Similar resultshavealsobeenfoundinCrucitti etal.’s work[15,16] inwhichtheymadestructuralvulnerabilityanalysis for the ItalianelectricpowergridandtheNorthAmerican powergrid respectively. The abovework is a good start to analyse vulnerability of power grids and complexnetwork theory inaugurates a new

directionforpowernetworksresearch.Thosemodels [4–12]were initiallyproposedforcomplexabstractnetworksanalysis and then were used in power networks [13–18]. However, those physicists’ work neglected some concrete engineering features. Therefore, there are good prospects for researchers to further investigate the complex problems by considering powersystemcharacteristicsandcomplexnetworktheorytogether.Particularly,electricpowernetworksarequitedifferent from those abstract networks. They are governed by Kirchoff’s Laws which feature in the system by the bus admittance matrix. The special characteristics result in a unique pattern of interaction between nodes in power grids. Therefore, for betterexplainingcomplexblackoutsofpowersystems, in thispaper,weproposean improvedmodelwhich isbasedonthe systembusadmittancematrix, representingthespecialelectrical topological structure.Moreover, thecomparisonbetween the traditional dynamicalmodel and theproposedmodelhasbeen investigated throughnumerical simulationswhichalso verify the superiorityof the improvedmodelboth inprecisionandefficiency. The rest of this paper is organizedas follows: Section2 introduces traditional physical topologicalmodels. In Section3,

firstly the characters of electrical structure of power grids are presented and then they are embedded into the improved model. Thenumerical simulationhasbeendisplayed inSection4andconclusionsare followed inSection5.

2. Traditional physical topological models

Traditional models were developed to the structural vulnerability analysis via error and attack resilience of complex networks andhavebeen employed in abstract networks and real-world networks (the Internet andpowernetworks etc.). Usually,agenericnetworkcanbedescribedasagraphGwithN nodesandK edges [19].G isdenotedbyanN×N adjacency matrix {wij}. In traditionalmodels, if there isnotanedgebetweennodes iand j, thenwij is set to0,otherwisewij equals to1 whichrepresentsthe lengthofeachedge. In fact, {wij}alsodescribesphysical topologicalconnectionsof thenetworkandall the relevantpropertiesof thenetworkdependonanalysisof theadjacencymatrix. Basedon thephysical topological graph structure, initially, traditionalmodels adopted a static analysismethod [4,8–10,14,17] for vulnerability analysis: a certain percentageofnodesaredeletedandthenthenewnetwork isevaluatedbytheperformance.Refs. [4,8–10,14,17] illustrated that the removal of a sizable group of nodes can have serious consequences. Actually, in most real networks, such as the Internet [18]andpowernetworks [20,21], themalfunctionofasingleorofaverysmall sizegroupofnodescancauseacas- cadingfailurewhichwouldleadthewholesystemtocollapseduetothedynamicsofredistributionofflowsonthenetworks. Thus, to take intoaccount thisphenomenon, adynamicalphysical topologicalmodelhasbeendeveloped [7,11,15,18]. Thedynamicalmodel iterativelyappliesa rule for the timeevolutionofG thatmimics thecascading failure followingby

removinganode randomlyor intentionally. Inorder toevaluatehowwell a systemworksbeforeandafter thebreakdown, the average efficiencyof anetwork is proposed in themodel. Refs. [7,11,15,18] assumed that the flowbetween twonodes i and j takes the shortest path dij connecting them. eij is used to denote the efficiency between nodes i and j. In short, the efficiencybetweennodes iand j isdefinedasthe inverseof theshortestdistance, i.e.eij = 1/dij,whichmeansthat the larger the dij is, the less efficiency that the information canbepropagatedbetweennodes i and j. If there is nopath between the nodes i and j, namely dij = +∞, eij = 0. Once the efficiency is known, the average efficiency of a network can then be calculatedas

E(G) = 1

N(N −1)

∑ i6=j∈G

eij. (1)

It is used tomeasure the performance ofG at a given iteration step. To characterize the load distribution in networks, the conceptof loadorbetweenness isused [7,11,15,18]. The loadatanode i isdefinedas the totalnumberof theshortestpaths passing throughthisnode.Thecapacityof anode is themaximumload that thenodecanhandle. Fora real-worldnetwork, the capacity is severely limited by cost. It can be assumed that the capacity Ci of node i is proportional to its initial load carriedby i

Ci = αLi(0) i= 1,2, . . . ,N, (2)

whereα >= 1isatoleranceparameterofthenetwork[13]andLi(0) is theinitial loadofthenode iat iterationstep t = 0[7, 11,15,18].Withsuchadefinitionof capacity, thenetwork is inastationarystate inwhich itoperateswithan initial average efficiencyE. The removalofanodesimulates themalfunctionof it and thentriggers thedynamicsof redistributionof flows

G.Chenet al. / PhysicaA388 (2009)4259–4266 4261

on thenetwork. In fact the removal of anodechanges the shortestpathsbetweennodesandconsequently thedistribution of loads,whichwouldcreateoverloadsonsomenodes.Ateach iterationstep t, the following iterative rule isadopted [7,11, 15,18]

wij(t +1) =

 wij(0)

Li(t) Ci

if Li(t) > Ci wij(0) if Li(t) ≤ Ci

(3)

where j extends to all the closest neighbours of i. In this way if at each iterative step, a node i is overloaded, the length of all the edges passing through it is increased, which can alter the shortest paths between vertexes, resulting in a new redistribution of the loads and then some nodes may be overloaded. The process will continue and produce a dynamic evolutionofnetworks, namelya cascading failure,whichcancause theaverageefficiencydegradation.

3. The proposed electrical topological model

Modelingphysical topologyofpowergridsasagraphandusingcomplexnetworkstheoryisagoodstart forvulnerability analysis of powernetworks.However, traditionalmodels ignore some important traits of powergrids. Inorder to improve accuracyandefficiency, theelectrical structureandpowernetwork featuresneed tobeconsidered.

3.1. Bus admittancematrix

Mathematically, we can form the power networks equations in the bus (or nodal) frame of reference in which the performance is described by n linear independent equations for n+ 1 number of nodes (the reference node, which is at groundpotential, is alwaysneglected) [22]. In theadmittance forms, theperformanceequationcanbedescribedas

IB = YBVB (4) where IB is thevectorof injectionbuscurrents.Theusualconventionfortheflowofcurrent isthat it ispositivewhenflowing toward thebus andnegativewhen flowingaway thebus [22].VB is the vector of busornodal voltagesmeasured fromthe slackbusandYB ={Ykl}n×n is thebusadmittancematrix. ExpandingEq. (4)∣∣∣∣∣∣∣

I1 I2 ·

In

∣∣∣∣∣∣∣ = ∣∣∣∣∣∣∣ Y11 Y12 · Y1,n Y21 Y22 · Y2,n · · · ·

Yn,1 Yn,2 · Yn,n

∣∣∣∣∣∣∣ ∣∣∣∣∣∣∣ V1 V2 ·

Vn

∣∣∣∣∣∣∣ (5) where

Ykl =

  −(Gkl + jBkl) k 6= l∑ k6=l

(Gkl + jBkl) k= l. (6)

The definition of the Y bus matrix used here [23] captures both the real and reactive (imaginary) portions of the line admittances (include resistances, capacitors, inductors and so on). For pairs of nodes k and l that do not share a direct physical connection,Ykl = 0, otherwiseYkl is theadmittanceof the linebetweennodeskand l.

3.2. The electrical topologicalmodel

Accordingtotheadmittancematrix, anewweightedundirectedgraphGcanbereconstructedwithN nodesandK edges to describe power networks by a updated N ×N adjacency matrix

{ wij } . Different from traditional models, the value of

wij is equal to |Ykl|, that is to say, if there is an edge between nodes i and j, the wij is an absolute value according to the admittance between i and j, otherwise it is 0.

{ wij } represents the efficiency of each edge, meaning that the larger the

admittance between two nodes, the more power can flow between them. Actually, the initial adjacency matrix {wij} in traditionalmodels also depicts the efficiency of each edge, because the lengths of edges are all equal to 1 [16]. Therefore, Ref. [16]proposedanequivalent formforcalculating theefficiencyeij basedontheadjacencymatrix {wij},whichneedsnot to calculate the shortest paths first. In theproposedelectricalmodel, thenewformof eij is employedandwemakea slight alteration forbetter explaining thedynamical behaviourofpowernetworksafter the removal of anode. Wedefine the electrical efficiencyof a pathbetween twovertices i and j as theharmonic compositionof efficiencies of

edgeson thepath [18]. Theharmonic compositionofN numbers x1,x2,x3, . . . ,xn is denotedas [16] n∑ i=1

(1/xn) −1 . (7)

Usually, different paths between nodes i and j have different electrical efficiencies. Thus, the efficiency eij is denoted as the maximum electrical efficiency and the corresponding path is regarded as an efficient path. To characterize the load distribution in the improvedmodel, the conceptof loadorbetweenness is still adopted. The load in theproposedmodel at

4262 G.Chenet al. / PhysicaA388 (2009)4259–4266

Table 1 Efficienciesof edges (bus i, bus j) in traditionalmodels.

Edges Efficiency Edges Efficiency Edges Efficiency Edges Efficiency

(1, 2) 1 (1, 5) 1 (2, 3) 1 (2, 4) 1 (2, 5) 1 (3, 4) 1 (4, 5) 1 (4, 7) 1 (4, 9) 1 (5, 6) 1 (6,11) 1 (6,12) 1 (6,13) 1 (7, 8) 1 (7, 9) 1 (9,10) 1 (9,14) 1 (10,11) 1 (12,13) 1 (13,14) 1

Table 2 Efficienciesof edges (bus i, bus j) in theproposedelectrical topologicalmodel.

Edges Efficiency Edges Efficiency Edges Efficiency Edges Efficiency

(1, 2) 16.0609 (1, 5) 4.3575 (2, 3) 4.9147 (2, 4) 5.3865 (2, 5) 5.4654 (3, 4) 5.4440 (4, 5) 22.6370 (4, 7) 4.8895 (4, 9) 1.8555 (5, 6) 4.2574 (6,11) 4.5369 (6,12) 3.5235 (6,13) 6.8445 (7, 8) 5.6770 (7, 9) 9.0901 (9,10) 11.0755 (9,14) 3.3471 (10,11) 4.7879 (12,13) 3.3566 (13,14) 2.5791

a node i is redefined as the total number of efficient paths passing through this node,which is an equivalent formof load defined inSection2. In the dynamical physical topological model, the average efficiency E(G) of a network is employed as a measure for

performance. However, for a power network, if the representation of adjacency matrix is different, the initial average efficiency is different. In order to easily compare the performance in differentmodels after a breakdownof networks, the averageefficiencycanbenormalizedby

NE(G) = E(G) E(G0)

(8)

where E(G) is the average efficiency of current networks at a given iteration step after occurrence of anymalfunction and E(G0) is the initial average efficiency of networks. Comparing with normalized average efficiency, another corresponding measure isdamagedefficiencyDwhich is thenormalizedaverageefficiency loss [16]

D= E(G0)−E(Gf) E(G0)

(9)

whereE(Gf) is the final efficiencyafter theendof the iterationprocess. Furthermore, once considering the adjacencymatrix as containing the efficiencies of edges, the iterative rule shouldbe

adjustedas

wij(t +1) =

 wij(0)

Ci Li(t)

if Li(t) > Ci

wij(0) if Li(t) ≤ Ci. (10)

Therule isalsoapplicable for {wij} in the traditionaldynamicalmodel, because{wij}depictsboth lengthsandefficienciesof edges. In thisway if at each iterative step, anode i is overloaded, theefficienciesof all linespassing through it isdecreased, which can adjust efficient paths, leading to a new redistribution of the loads and then some nodes may be overloaded. Therefore, the proposed electrical topologicalmodel can also cause a dynamic evolution of power networks after removal of anode.

3.3. A case study

In this section, a simple power network model (the IEEE 14 bus system shown in Fig. 1) is used to illustrate the characteristics of physical topological structure and electrical topological structure. The system has 14 buses and 20 branches, which forms a small network with 14 nodes and 20 edges. The relative smaller side of this 14 bus system also enables a clear illustrationof theoriginal powernetwork (Fig. 1) and the correspondingphysical topological graph (Fig. 2). Traditional models [13–18] only consider the physical connections. Thus, the adjacency matrix {wij} is calculated by the physical topological graph (Fig. 2) according to the rule that if there is an edge between nodes i and j then wij is set to 1, otherwise wij equals to 0. Table 1 displays the efficiencies of edges in traditional models which show that all edges have equal function forpropagate flow. In fact, the real situation in the IEEE14bus system isnot like this case.Wecannotneglect the influenceof resistances,

capacitors, inductors etc., which have been considered in the proposed electrical topological model by calculating bus admittance matrix. Table 2 illustrates the efficiencies of edges in the proposed model which shows that the function of eachedge isdifferent.

G.Chenet al. / PhysicaA388 (2009)4259–4266 4263

G

C

GENERATORS

SYNCHRONOUS COMPENSATORS

12

G

G

C

5

6

11

13

14

10

2

C

3

1

THREE WINDING

TRANSFORMER EQUIVALENT

9

7 8

4

C

4

7 8

9 C

Fig. 1. The IEEE14-bus system.

1 2

3

4

5

6

7

8 9

10

11

12

13

14

Fig. 2. Physical topological graphof the IEEE14bus system.

4. Numerical analysis

In this section, the standard IEEE300bus test system is employed. The systemhas300buses (nodes) and411branches (transmission lines),whichcan represent a commoncomplexnetworkwithcomplex features [23]. Theproposedelectrical topologicalmodelandthetraditionaldynamicalphysical topologicalmodel (wecall itphysical topologicalmodel for short) are implementedbyMATLABR2007. In thesimulation,weadoptsinglenoderemoval strategy from[7,11,15,18], inwhicha singlenode is removedat iterationstep t = 0andthenetwork isevolvedbyiterationaccordingtoEq. (10).Thestrategycan simulate theeffectof externalperturbationofpowernetworks, because in fact, theprobabilityofmultiplenodal failuresat thesametimeisquitesmall. Followingthescheme,wechoosenodeseither randomly(randomfailure)orselectivelybythe highest load (load-based attack) andonce removed, thenormalized efficiency of thenetwork and the loads of nodeswere continually recalculatedby iteration. Fig. 3 displays the final normalizedefficiency (100 iterations or iteration step t = 100) after removal of a randomnode

andahighest loadnode for thephysical topologicalmodel and theproposedelectrical topologicalmodel. Fig. 4 shows the correspondingdamagedefficiency.Therandomfailurecurveswereobtainedbyaveraging100individualremovals. Itcanbe seenfromFig.3 that, roughly, thetendencies inthetwomodelsaresimilar.Bothof themcanreveal that thepowernetwork is robust to randomperturbations, yet disturbances affecting key nodes greatly reduce its ability to function. However, in detail thedifference for the twomodels is quite clear. For example,when the toleranceparameterα is 1, under load-based attack, the normalized efficiency of the physical topological model decreases to below 0.1, which means that almost the wholepower gridwill experience ablackout. Itmaybenot practical from theknowledgeof apower engineer. It is known thatthemostseriousblackout inhumanhistory,namelythe2003US–CanadaBlackout, involvedalargeareaofAmericaand a little part of Canada,whichdidnot break thewholepowernetwork. As a result the correspondingnormalized efficiency

4264 G.Chenet al. / PhysicaA388 (2009)4259–4266

1.3 1.6 1.9 2.2 2. 0.

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

random failure of electrical model

load-based attack of electrical model random failure of physical model

load-based attack of physical model

1

1

N o rm

a liz

e d E

ff ic

ie n cy

E

1 5

Tolerance parameter α

Fig. 3. FinalnormalizedefficiencyE of the twomodelsunderdifferenta.

1.3 1.6 1.9 2.2

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8

9 random failure of electrical model

load-based attack of electrical model random failure of physical model

load-based attack of physcial model

0

0.

D a m

a g e d E

ff ic

ie n cy

D

2.51

Tolerance parameter α

Fig. 4. FinaldamagedefficiencyDof the twomodelsunderdifferenta.

10 20 30 40 50 60 70 80 90 5

0.55

0.6

0.65

0.7

0.75

0.8

0.85

0.9

0.95 a=1

a=1.3

a=1.6

a=1.9

a=2.2

a=2.5

0.

1

N (E

)

0 100

Iteration t

Fig. 5. Randomfailure in thephysical topologicalmodel.

fallsto0.35intheelectricaltopologicalmodelismorepreciseandreasonable.Furthermore,whenα increasesto2.5meaning that there is about 1.5 times for securitymargin, in thephysical topologicalmodel, thenormalizedefficiency is only about 0.5 under load-based attack. That is to say still about half of the power network will be out of function, which is also not practical under suchhigh securitymargin. About 0.8 in the improvedmodel ismuchmore appropriate. Similar results can be observed fromother parts of Figs. 3 and 4. Thus, the proposedmodel ismore precise and reasonable than the physical topologicalmodel.

G.Chenet al. / PhysicaA388 (2009)4259–4266 4265

10 20 30 40 50 60 70 80 90 0.

0.7

0.75

0.8

0.85

0.9

0.95

a= 1

a= 1.3

a= 1.6

a= 1.9

a= 2.2

a= 2.5

65

1

N (E

)

0 100

Iteration t

Fig. 6. Randomfailure in theelectrical topologicalmodel.

10 20 30 40 50 60 70 80 90 0.

0.2

0.3

0.4

0.5

0.6

0.7

0.8

0.9

a= 1

a=1. 3

a=1. 6

a=1. 9

a=2. 2

a=2. 5

1

1

N (E

)

0 100

Iteration t

Fig. 7. Load-basedattack in thephysical topologicalmodel.

10 20 30 40 50 60 70 80 90

0.4

0.5

0.6

0.7

0.8

0.9

1

N (E

)

a= 1

a= 1. 3

a= 1. 6

a= 1. 9

a= 2. 2

a= 2. 5

0 100

Iteration t

Fig. 8. Load-basedattack in theelectrical topologicalmodel.

Figs. 5–8 illustrate thedynamicalbehaviourofnormalizedefficiencyafteroccurrenceofbreakdownofanodewithboth iterationstep t andthetoleranceparameterα increasing.Thedynamicalbehaviour istriggeredbyremoval(randomorload- based) at t = 0. The results showdynamical normalized efficiency for six values of α, namely α = 1, 1.3, 1.6, 1.9, 2.2 and 2.5 respectivelyunder iteration step t from0to100. FromFigs. 5–8, it is very clear that the improvedmodel ismuch faster

4266 G.Chenet al. / PhysicaA388 (2009)4259–4266

to converge to stable points. For example, it can be observed fromFig. 5, under random failure,when α = 1 the physical topologicalmodel needs about 90 iteration steps (t = 90) to obtain a stable point,while the electrical topologicalmodel only needs about t = 20 (Fig. 6). For a = 1.3, the physical one is about t = 80 (Fig. 5), the electrical one is about t = 35 (Fig. 6).As for load-basedattack, the situation ismuchmoreobvious. For all toleranceparametervaluesofα, the improved modelcanconvergeat lessthan20iterations(Fig.8).However, inthephysical topologicalmodel,convergencesneedtotake about 60–70 iterations (Fig. 7). Therefore, it can be concluded that the electrical topologicalmodel ismuchmore efficient than thephysical one.

5. Conclusions and further work

Complexnetworkstheorywas initiallyproposedbyphysicistsandthenwasemployedindifferent fields, includingsome initial explorationwithpowernetworks.However, thosephysicists’workhadnot considered theunique characteristics of power networks. Consequently, such simplified physical-based approach often leads to unrealistic resultswhich basically limit theirpotential andsignificance forengineeringapplications. Theyneedtobe improvedbyembeddingpowernetwork features inorder toobtainmorerealisticandreliablevulnerabilityassessment result. In thispaper, anelectrical topological modelisproposedforstructuralvulnerabilityanalysisofpowergrids, inwhichanewweightedundirectedgraphstructureis presentedbasedonKirchoff’sLawsandbusadmittancematrix.Themodel isacloserapproximationtorealpowernetworks. Furthermore, the numerical tests on the IEEE 300 bus system also display that the proposed model is more precise and efficient. Theproposedelectrical topologicalmodel canbeuseful tomakeavulnerability assessmentand todesign specific action to reduceweaknessesofpowergrids. Sincepowernetworksblackouts arediverse andcomplicated, it is impossible to consider all factors in this paper alone.

Instead, we demonstrate that complex networks theory provides a new direction for complex power networks research. However, at present this work is still in its early age and a more comprehensive complex network theory based power networkvulnerability analysis is our futurework.

References

[1] D.N.Kosterev, C.W.Tayor, IEEETransactionsonPowerSystems14 (1999)967–979. [2] U.S–Canada power system outage task force, Final Report on the August 14, 2003 Blackout in the United States and Canada: Causes and Recommendations,April 2004.

[3] D.Q.Gan, J.Y.Hu,AutomationofElectric PowerSystem28 (2004)1–5. [4] R.Albert,H. Jeong,A.-L. Barab’asi,Nature406 (2000)378. [5] V. Latora,M.Marichiori, PhysicalReviewLetters87 (2001)198. [6] P. Crucitti, V. Latora,M.Marchiori, PhysicaA320 (2003)622. [7] P. Crucitti, V. Latora,M.Marchiori, PhysicaA340 (2004)388–394. [8] R. Cohen,K. Erez,D. ben-Avraham,S.Havlin, PhysicalReviewLetters85 (2000)4626. [9] R. Cohen,K. Erez,D. ben-Avraham,S.Havlin, PhysicalReviewLetters86 (2001)3682. [10] A.E.Motter, T.Nishikawa,Y. Lai, PhysicalReviewE66 (2002)065103. [11] H. Zhao, Z.Y.Gao, TheEuropeanPhysical JournalB67 (2007)95–101. [12] S. Sun, Z. Liu, Z. Chen, PhysicaA373 (2007)851–860. [13] A.E.Motter, Y.C. Lai, PhysicalReviewE66 (2002)065102. [14] M.R. Casals, S. Valverde, R. Sole, International Journal ofBifurcationandChaos17 (7) (2007)2465–2475. [15] P. Crucitti, V. Latora,M.Marchiori, PhysicaA338 (2004)92–97. [16] R.Kinney, P. Crucitti, R.Albert, V. Latora, TheEuropeanPhysical JournalB46 (2005)101–107. [17] R.Albert, I. Albert,G.L.Nakarado, PhysicalReviewE69 (2004)025103 (R). [18] P. Crucitti, V. Latora,M.Marchiori, PhysicalReviewE69 (2004)045104. [19] D.J.Hill, G.R. Chen, in: IEEE International SymposiumonCircuits andSystems, 2006. [20] R. Baldick, B. Chowdhury, I.Dobson, et al., in: IEEEPowerEngineeringSocietyGeneralMeeting, PittsburghPA,USA, July2008. [21] IEEEPESComputerandanalyticalmethodssubcommittee, in: IEEEPowerandEnergySocietyPowerSystemsConferenceandExposition2009,Seattle,

WA,USA. [22] J.Das, PowerSystemAnalysis,MarcelDekker, Inc., 2002. [23] P.Hines, S. Bluemsack, in: Proceedingsof the41stHawaii InternationalConferenceonSystemSciences, 2008.