Adjacent Graph Based Vulnerability Assessment for Electrical Networks Considering Fault Adjacent Relationships Among Branches

Security issues related to vulnerability assessment in electrical networks are necessary for operators to identify the critical branches. At present, using complex network theory to assess the structural vulnerability of the electrical network is a popular method. However, the complex network theory cannot be comprehensively applicable to the operational vulnerability assessment of the electrical network because the network operation is closely dependent on the physical rules not only on the topological structure. To overcome the problem, an adjacent graph (AG) considering the topological, physical, and operational features of the electrical network is constructed to replace the original network. Through the AG, a branch importance index that considers both the importance of a branch and the fault adjacent relationships among branches is constructed to evaluate the electrical network vulnerability. The IEEE 118-bus system and the French grid are employed to validate the effectiveness of the proposed method.


L
Set of branches (i.e. branches, transformers) in a transmission network, L = · · · , L j , · · · , card|L| = N L . B Set of nodes (i.e. buses) in a transmission network, card|B| = N B . N W Number of generations in a transmission network. L i Set of branches in fault chain i, L i = · · · , L i j , · · · , L i ⊆ L, card L i = n i . X i Set of adjacent fault relationships among branches in fault chain i, X i = · · · , X i j , · · · , X i j = L i j → L i j+1 , card X i = n i − 1. α i j Loading assessment index of branch j during contingency x in fault chain i generation process, L j ∈ L.
Set of edges in g g g i , E i = · · · , E i j , · · · , card E i = n i − 1. ψ(·) Sequential mapping function. r Power exponent of cumulative distributions. R 2 Fitting effect of power law. Generally, R 2 ≥ 80% has a satisfactory fitting effect.

LOAD REDISTRIBUTION MODEL BASED COMPLEX NETWORK k a
Degree of vertex a in the AG. k in a In-degree of vertex a in the AG. ρ, τ Scale factor for initial load in the AG. ϕ a Initial load of vertex a in the AG. σ Vertex tolerance parameter, which measures the ability to tackle the extra load in the AG.

I. INTRODUCTION
Electrical network vulnerability (ENV) assessment, which is also called critical branch identification, is necessary for operators to improve the security and reliability of electrical networks.
At present, ENV assessment is studied mainly from attack and defense perspectives. On the attack side [1], the aim is to cause large-scale power outage as much as possible, which results in any property loss directly or social disorder indirectly by attacking the critical branches (or nodes) at minimal cost. On the defense side [2], focusing on protecting the critical branches is conducive to reducing the risk of network faults, especially, cascading failures [3], [4].
Although their research motivations are different, their essence is to effectively detect the weak branches (or nodes) of transmission networks. In the traditional methods, the operational indices based on steady and transient operation of power systems have been studied [5], [6]. The dynamic index based transient stability is proposed to identify the critical branches in reference [5]. The assessment index considering the transient impact is employed to analyze the system vulnerability [6]. In addition, machine learning [7], [8] and random matrix [9] as tools of analysis are also applied to reveal the features of system vulnerability.
As the size of electrical networks increases, the electrical networks are becoming one of complex artificial networks and have some universal topological features with other complex networks [11], such as communication network. Therefore, the complex network theory applied to assess ENV has been recognized as an important way. Many scholars have employed pure topological metrics based on CNT, such as, average path length [12], betweenness [13], centroid [1] and degree [14].
However, the pure topological metrics neglect some physical characteristics of the transmission networks; therefore the extended topological metrics which integrate electrical quantity into the topological metrics have attracted wide interests. For instance, in reference [15], the centrality index is redefined by considering the maximum flow from the generator nodes to the load nodes to evaluate the networks. Reference [16] studies the structural vulnerability by introducing the system bus admittance matrix into the topological model. Taking the actual path of power flow into consideration, reference [17] redefines the betweenness. Reference [18] uses the maximal demand of load and the capacity of generators to define the electric betweenness.
Although the extended topological metrics can reflect some special physical characteristics of electrical networks, there are still some problems. From the perspective of electrical network, the majority of metrics only reflect the physical characteristics under normal operating conditions but the features, such as fault adjacent relationships among branches, under fault operating conditions, are still not considered. Moreover, the metrics still focus on the topological vulnerability of the network although the metrics are improved by integrating some electrical quantity. From the perspective of CNT, the vulnerable metrics belong to basic static statistical methods and CNT-based other models are rarely employed to study ENV. The main reason is that the models [19], [20] are studied based on the topological relationships among vertices or edges but electrical networks must comply with physical and operational rules; therefore the topological relationships cannot comprehensively reveal the correlation among branches or nodes in an electrical network.
To overcome the limitations, reference [21] proposes an identification method of critical nodes via single Laplacian matrix inversion, which goes beyond pure graph theoretical approaches. In this paper, from the different perspective, we construct an adjacent graph (AG) based on the ideas of statistical graphs [22] by comprehensively considering topological features, physical features and operational features of an electrical network. We employ AG instead of the original electrical network to indirectly assess ENV. Through AG, we construct a branch importance index (BIM) based on load redistribution model (LRM) of CNT, because the topological structure of AG can reveal not only the importance of a branch, but also the adjacent fault relationships among branches intuitively and clearly. To validate the proposed method, simulation analysis has been accomplished in IEEE 118-bus system and French grid.
The remainder of this paper is organized as follows. Related works and existing problems of LRM are discussed in Section 2. Section 3 gives the AG generation method based on fault chain theory. In Section 4, LRM is analyzed and the construction method of BIM based on AG is introduced. In Section 5, the proposed method is validated by employing the IEEE 118-bus system and French grid described in Tab. 1. Finally, the conclusions are presented in Section 6.

II. RELATED WORKS AND EXISTING PROBLEMS
At present, there are many complex networks (e.g. electrical networks, computer networks, and traffic networks) in the real world. In these networks, the nodes or branches carry certain load. The load could be concrete or abstract substance, information, or energy in different networks. Generally, in a complex network, the nodes' or branches' load capacity is finite. When a node or branch fails (or is removed) in the network, the load of the entire network will be redistributed. Consequently, the load redistribution can cause other nodes or branches overloaded, which can trigger a cascading failure. In references [23]- [26], many scholars have taken the nodes as research objects to construct LRM and then the  network vulnerability is assessed based on LRM. In LRM, the research mainly includes two aspects: load model and load redistribution strategy. In the load model, initial load and load capacity of nodes must be considered. In a complex network, degree (or betweenness) of a node can reflect its importance, therefore the degree (or betweenness) is universally employed to construct the load model of nodes. When a node fails in the load redistribution strategy,, the increment of load of a nonfaulty node, which is allocated by the faulty node, depends on the shortest path between the faulty and the non-faulty node [25], [27].
However, compared with general networks, after a branch or node fails in an electrical network, the load redistribution of the entire network depend on physical and operational rules not only topological structure; therefore CNT is not very suitable to construct LRM of electrical networks. We use the IEEE 9-bus system as an example as shown in Fig.1. When branch 8-9 is removed, the power flow of the adjacent branches 8-7 and 2-8 do not increase compared to the other branches. In contrast, the power flow of some non-adjacent branches such as branch 8-7 decreases. In summary, we can conclude that only considering the topological structure is improper during the construction of LRM for electrical networks; therefore CNT-based LRM cannot be employed to assess ENV.
To solve the above-mentioned problem, inspired by statistical graphs [22], [28]- [30], we propose an adjacent graph (AG) to construct LRM of the electrical network and then we can assess ENV based on LRM. The framework of ENV can be shown in Fig. 2. First, we can abstract  the importance of a branch and adjacent fault relationships among branches based on the topological features of AG. Furthermore, we employ both the importance of the branch and adjacent fault relationships to construct LRM. At final, we construct a branch importance index (BIM) based on LRM (LRM) of CNT to assess ENV. In addition, we intend to reveal whether the network vulnerability is related to the fault adjacent relationships among branches; therefore we introduce the adjustable parameters into BIM.

III. ADJACENT GRAPH
AG is a directed graph derived from fault chains, which can reveal the adjacent fault relationships among branches. In this section, we focus on the construction method of AG and its topological properties.

A. FAULT CHAINS 1) BASIC DEFINITION OF FAULT CHAINS
After L 1 is removed in the current contingency as shown in Fig. 3, we can choose one of three candidate branches L 2 , L 3 and L 4 as the fault branch to trip in the next contingency. Compared with L 2 and L 4 , if L 3 can lead to the maximum damage of network function/ structure after removal, we can define the removed sequence as a fault adjacent relationship between L 1 and L 3 . To investigate the adjacent relationships among branches, for an electrical network with N L branches, we need to calculate N k L in the n-k contingency. It is difficult to complete the work for a large-scale electrical network.
To overcome the problem, we can trace the adjacent relationships among branches by constructing an index to reduce the computational burden and then we take every branch as an initial point to construct the fault chains as shown in Fig.4.
For a fault chain, we propose a branch loading assessment index (BLAI) to trace the adjacent relationships among branches [31], [32]. BLAI can reflect the loading burden of a branch and the possibility of failure of the branch under current contingency. BLAI is represented as We adopt the load shedding percentage [22] to measure the scale of power blackout and mark the end of a fault chain. The load shedding rate is defined as where P i dx is the active power of load bus during contingency x in fault chain i; δ i x is the load shedding percentage during contingency x in fault chain i. To determine the end of the fault chain, we define a threshold [30], [33]. When i ≥ , we terminate the fault chains generation process.

2) GENERATION ALGORITHM OF FAULT CHAINS
To simplify the calculation, we employ the DC power method [31] to calculate the transmitted power over branches and the DC optimal power flow (DC OPF) algorithm [32], [34] to calculate the i . The flow diagram of generation algorithm of fault chain C i = L i , X i , n i is shown in algorithm1, where L i represents the set of branches in fault chain i; X i represents the set of adjacent fault relationships among branches in fault chain i; n i represents the number of branches in fault chain i.

B. AG GENERATION METHOD
For a fault chain C i [1], we employ mapping function : [35], i.e. g g g i = (C i ), where V i is the set of vertices and E i is the set of edges. The g g g i = V i , E i is defined as follows: Set of vertexes: According to : The set of vertices is defined as  Set of edges: According to : The set of edges is defined as Mapping relationships among vertexes and edges: Employ the function ψ : For a given electrical network with N L branches, we can take each branch as a starting point to develop N L fault chains, i.e., C 1 , C 2 , . . . , C N L . Furthermore, we map the N L fault chains to N L graphs, i.e., g g g 1 , g g g 2 , . . . g g g N L , according to Equations (4)- (6). Finally, we merge N L graphs to construct the AG of the electrical network. The AG is represented as

C. GRAPH PROPERTIES OF AGs
To analyze the graph properties of AGs, we construct AGs of IEEE 118-bus system (Fig.5) and French grid. It is noted that AG of the French grid is not given due to the space limitation. Here, we establish the threshold = 20%. 20% power loss is a big enough blackout event for an electrical network REF_Ref494789773 \r \h [6]. As can be seen in Fig.5, compared with electrical networks that are the spatial association networks, AGs are the temporal correlation networks that can reveal the fault adjacent relationships among branches. We further analyze the cumulative distribution of vertex degree in-degree and out-degree in AGs. The cumulative distributions of the vertex degree, in-degree and out-degree P (K > k) = K >k P (K > k) are from power laws family P (K > k) ∝ k −(r−1) ; therefore we can conclude that AGs are scale-free graphs. In addition, it is noted that the values of all exponents r exceed 2, which satisfies the features of power law distributions [37]. Similar conclusions can be also drawn in other systems, for example the IEEE 39-and 300-bus system.
(1) for the IEEE 118-bus system The scale-free properties indicate most of vertices have small degree (in-degree or out-degree), but there are a few vertices with high degree (in-degree or out-degree); therefore two text benchmarks have high robustness under random attacks, but there is high vulnerability under intentional attacks. In other words, for a vertex (i.e. branch) with high degree, its importance for improving the robustness of a network is larger than that with small degree. The vertexes can be further subdivided. Some with high in-degree can be affected easier than those with small indegree once other vertices fail. For some vertexes with high out-degree, once they fail, they are easier to result in fault propagation.

IV. CONSTRUCTING LRM AND BIM OF ELECTRICAL NETWORK USING AG
Thanks to the scale-free properties of AG, we can employ the topological structure of AG to construct LRM of electrical network indirectly. Before constructing LRM, we introduce some concepts based on CNT. We define the load called the vulnerability flow that exists in AG. Each vertex of CG carries a certain proportion of the vulnerability flow. The proportion of vulnerability flow over vertex can reflect importance of the corresponding branch in the original electrical network. The edges of CG reflect the transmission path for the vulnerability flow among vertices, which can reveal the fault adjacent relationships among branches in the original electrical network.

A. LRM BASED ON AG
In LRM, we need to construct the load model and load redistribution strategy. In AG, the vertex with higher degree can reflect the larger importance of the corresponding branch in the original electrical network; therefore we can employ the degree k a to define the initial load (i.e. vulnerability flow) of vertices. The initial load ϕ a of vertex a is represented as where ρ, τ are scale factor for initial load. Furthermore, let the load capacity [38] of vertex a be proportional to ϕ a as follows: where σ is the vertex tolerance parameter; σ ϕ a represents capacity margin. When a vertex b fails, its vulnerability flow will be transferred to the other non-faulty vertexes [39] as follows: where ϕ b a is represented as In Equation 10, φ b a l ab , θ, k in a , β is the proportional coefficient which is defined as follows: where Q b adjusts the range of load distribution of the faulty vertex b and proportion control factors θ ∈ [0, +∞) , β ∈ [0, +∞) ; l ab is the shortest path between vertices a and b, and k in a is the in-degree of vertex a. In Equation 11, φ b a l ab , θ, k in a , β depends on l ab and k in a . l ab adjusts the distance of load redistribution, which reflects the adjacent fault relationship between a and b. In addition, when a non-faulty vertex has higher in-degree, the non-faulty vertex can be easily affected by faulty vertex, which indicates that the nonfaulty vertex can receive more load from the faulty vertex; therefore we introduce k in a into φ b a l ab , θ, k in a , β . In addition, for the non-faulty vertex a, if a ∈ Q b , Equation (12) is satisfied.
As is shown in Equations (11) and (12), load redistribution strategy is decided by l ab and k in a (a ∈ Q b ).

Construction of BIM:
When vertex b fails, we can introduce the load increment to analyze the fault probability of nonfaulty vertices under the next contingency as follow: When γ > 1, i.e., load increment ϕ b a is greater than capacity margin σ ϕ a , vertex a is overloaded, which could cause vertex a fails under the next contingency. When γ ≤ 1, the opposite conclusion occurs. To simplify the calculation, we use Equation (14) instead of Equation (13).
Plug Equations (4), (7) and (8) into Equation (14). γ can be further represented as γ can qualify the degree of overload of vertex a after vertex b fails [25], which can reflect the impact of faulty vertex b on vertex a. Furthermore, we employ the average value of γ of all non-faulty vertices to define BIM of vertex b as follows: where |Q b | is the number of vertexes in the Q b . If χ b ≥ 1, faulty vertex b causes other vertices to be overloaded with high probability; If χ b < 1 |Q b |, the opposite occurs; Here, we call χ b = 1 the dangerous threshold. We call χ b = 1 |Q b | the secure threshold. In addition, it is noted that χ b is also dependent on tolerance parameter σ . We define two thresholds of σ : dangerous critical value σ d and secure critical value σ s . When σ ≤ σ d , χ b of all vertices is larger than the dangerous threshold, i.e., any faulty vertex can cause other vertices to be overloaded. When σ ≥ σ s , χ b of all vertices is smaller than the secure threshold, i.e., any faulty vertex cannot cause other vertices to be overloaded.
1. When θ = ∞, according to Equation (12), |Q b | = k b , i.e., Q b = L a |∃e q ∈ E, e q = L a L b ∨ e q = L b L a . Therefore the load redistribution is the nearest neighbor redistribution rule(NNR), i.e., the faulty vertex only transfers its load to its adjacent non-faulty vertices. In the corresponding original electrical network, θ = ∞ demonstrates that we only consider the fault adjacent relationships among branches under N -1 secure criterion.
2. When θ = 0, according to Equation (12), Therefore the load redistribution is the uniform sharing redistribution rule (USR), i.e., the faulty vertex transfers its load to all non-faulty vertices. In the corresponding original electrical network, θ = 0 demonstrates that we consider the fault adjacent relationships among branches under N -k secure criterion where k = max n i |i = 1, 2, . . . ., N L . In summary, in GA, θ decides the range of load redistribution of faulty vertex. In the corresponding electrical network, θ decides the fault adjacent relationships among branches under different secure criterions.

V. SIMULATION AND ANALYSIS
Using the IEEE 118-bus system and French grid as an example, we calculate BIMs of branches to rank the top critical branches shown in Tab. 2. Tab. 2 shows the BIMs of branches under NNR(θ = ∞) and USR (i.e. θ = 0), respectively, where β = τ = 1, σ = 0.4. In the IEEE 118-bus system, when load redistribution is NNR, branches 7, 38,9,125,94,127,126,37,137 and 92 are the 10 most critical branches. By contrast, when load redistribution is USR, branches 7,9,31,137,104,61,30,96,105 and 90 are the 10 most critical branches. It demonstrates that BIMs and rankings of branches will also change as load redistribution rule changes.

A. ANALYSIS OF ENV BASED ON THE TOLERANCE PARAMETER σ
We analyze two types of thresholds (dangerous critical value σ d and secure critical value σ s ) of tolerance parameter σ of a network. Figs.6-7 show that the regularity of changes of two types of thresholds as τ and β increase under different redistribution rules in two text benchmarks. For σ s , whether load redistribution is the NNR or USR in two text benchmarks, σ s increases as β increases. By contrast, σ s decreases as τ increases.
For σ d , on the IEEE 118-bus system, under NNR, 1) when β is set at a fixed value, τ has little effect on the σ d . 2) When τ is set at a fixed value, σ d decreases with β increasing. Under USR, 1) when β is set at a fixed value and β ∈ (0, 1], σ d decreases as τ increases. When β is set at a fixed value and β ∈ (1, 2], τ has little impact on σ d . 2) When τ is set at a fixed value and τ ∈ [0, 0.8], σ d first increases and then decreases with increasing β. When τ is set at a fixed value and τ ∈ (0.8, 2], σ d decreases with increasing β.
For the French grid, under NNR, τ has no impact on σ d when β is set at a fixed value; σ d decreased with increasing β when τ is set at a fixed value. Under USR, σ d decreases with increasing τ or β.
In addition, comparing the IEEE 118-bus system with the French grid, it is not difficult to find that French grid is less likely to trigger a cascading failure because σ d of the French grid is smaller than that of the IEEE 118-bus system but σ s of the French grid is larger than that of the IEEE 118-bus system when β and τ of two systems are set at the same value, respectively. This demonstrates that the French grid have lower vulnerability than the IEEE 118-bus system.

B. RELATIONSHIP BETWEEN BIM AND DEGREE
The relationship between BIM and degree with different values of τ : Fig. 8 shows the relationship between BIM and degree under different redistribution rules in two test benchmarks when τ = 1.4 and 0.1 respectively. Under NNR, there is a negative correlation between BIM and degree generally VOLUME 7, 2019  when τ = 0.1. Conversely, when τ = 1.4, there is a positive correlation between BIM and degree generally. Therefore we can conclude that τ decides the positive or negative correlation between BIM and degree. Under the USR, regardless of the values of τ , there is a positive correlation between BIM and degree.
The relationship between BIM and degree with different values of β: Figs.9-10 show the relationships between BIM and degree under different redistribution rules. β affects  the steep extent of curves but has no effect on positive or negative correlation between BIM and degree.
By changing τ and β, we conclude that τ has a larger impact on the relationship between BIM and degree than β. Especially under NNR, values of τ decide the correlation between BIM and degree.
In summary, through the above analysis, BIMs of branches are affected mainly by θ and τ which decide the vulnerability rankings of the branches.

C. ANALYSIS OF CRITICAL BRANCHES
Critical branch attacking: We investigate the percentage of remaining load of systems when attacking the critical branches ranked by BIMs. We respectively remove (1) the top 10 and 20 critical branches on IEEE 118-bus system, (2) the top 60 and 120 critical branches on French grid, ranked by BIMs under the different values of θ and τ as shown in Fig.11. In the figure, with changing θ and τ , the percentage of load remaining load changes obviously. In addition, as the number of attacking branches increases, we find that the percentage of load remaining load decreases generally as θ and τ simultaneously increase on the IEEE 118-bus system while the results are opposite on the French grid.
Furthermore, we investigate the percentage of remaining load with increasing critical branches as shown in Fig.12. The ranking method of critical branches have five ways: 1) random, 2) degree of AG, 3) LIMs under the different values of θ and τ , 4) topological metrics (including betweenness, electrical betweenness, network efficiency [47] and network ability [30]) based on topological structure of electrical network and 5) state metric based operational state of electrical network [48]. In Fig.12, the percentage of remaining load drops faster when removing the critical branches ranked by BIMs than those ranked by random and degree. In addition, taking the IEEE 118-bus system as an example ( Fig. 12(a)), when the critical branches ranked by BIMs under different values of θ and τ excluding θ = ∞ and τ = 0.1, the results are better or equal to those ranked by degree of AG; therefore we can infer that only considering the importance of a branch itself cannot comprehensively reveal the vulnerability of the branch. In the future, we can employ some optimization algorithm (e.g., genetic algorithm, particle swarm optimization) to adjust these parameters so that the accuracy of identification of critical branches can be improved.
In addition, we compared the proposed method with the existing methods (i.e., topological metrics and state metrics) shown in Fig. 12. In the figure, we can find that our proposed method is better than the topological metrics and state metric after the removal of suggested branches. It is noted that even if we do not consider the fault adjacent relationships among branches, i.e., we only employ the degree of AG to identify the critical branches, the attacked results is still better than the topological metrics and state metric. In summary, it demonstrates that employing BIMs to identify critical branches of electrical networks and to analyze electrical network vulnerability is valid.
In summary, we can conclude that 1) it is a reasonable idea to employ AG instead of the original electrical network to assess the ENV, which can effectively overcome the limitations of CNT applied to the ENV; 2) the vulnerability of branches is related to not only the own importance of branches but also the fault adjacent relationships among branches.

VI. CONCLUSIONS
In this paper, we employ AGs instead of original electrical networks from the perspective of CNT to evaluate the ENV. The topological, physical and operational features of an electrical network is considered during the construction of AG. Then, AG is adopted to construct LRM which considering both the importance of branches and fault adjacent relationships among branches. Furthermore, we employ LRM to construct BIMs. The simulation analysis demonstrates that employing BIMs to evaluate the ENV is valid and reasonable. Meanwhile, considering the fault adjacent relationships in BIMs can obviously improve the accuracy of identification of vulnerable branches. In addition, our works provide a new view of electrical network vulnerability in this field.