# Identification of Critical Nodes and Accident Chains of Power Systems under the Effect of Geomagnetic Storms

Authors: Qing Liu; Nan Ding; Jin-Long Li; Xiao-Hu Yi
Journal: Mathematical Problems in Engineering
 (2022)

---

## Abstract

Geomagnetically induced currents (GICs) generated by geomagnetic disturbances (GMD) during strong magnetic storms may trigger chain failures in power systems. In this paper, a critical node identification method based on the improved LeaderRank algorithm considering the effects of geomagnetic storms is proposed. A chain fault identification process based on the fragility theory of complex systems is proposed to determine the initial fault line by critical nodes and simulate the fault development process to obtain the accident chain set. The effectiveness of the method is verified by using the IEEE RTS79 node system as an example. The research results show that the method can reflect the critical nodes and accident chain set of power systems under different geomagnetic storm effects, and the research results can provide a reference for geomagnetic storm disaster quantification, assessment, and prevention.


---

## Body

## 1. Introduction

Solar activity causes rapid changes in the ionosphere in extraterrestrial space resulting in strong geomagnetic disturbances (GMD). GMD generates low-frequency geomagnetically induced currents (GIC) in power systems. The GIC may lead to half-cycle saturation of transformers, producing secondary effects such as temperature increase, harmonics, and reactive power loss [1, 2]. Not only high-latitude regions such as North America and Northern Europe have been affected by GMD [3, 4], but also mid- and low-latitude regions such as Brazil, Japan, South Africa, and China have been attacked with serious consequences [5–9]. Power systems are becoming increasingly complex and interconnected. In recent years, a number of major outages have occurred worldwide [10, 11], which have caused huge economic losses and attracted widespread attention. Through the study of previous cases, it can be found that major outage accidents are often a series of chain failures caused by the failure of some components. These components play an extremely important role in the stable operation of the power grid, called critical nodes, and their dismantling can lead to large-scale shifts in power system flows, causing other nodes or branches in the system to operate at ultralow voltage or overload, which can lead to power system chain failures. Therefore, it is important to accurately identify the critical nodes of the power system and the chain of accidents they cause to ensure the safe and reliable operation of the power system. The current research on grid disaster risk assessment under geomagnetic storm conditions focuses on the GIC level and its derived effects [12–15] and has not yet conducted in-depth studies on critical nodes and accident chains.There are two main types of research on chain failures triggered by critical nodes in power systems in academia. The first type is based on complex system theory, which considers the dynamic process of chain failures to assess the criticality of nodes. There are mainly cascade models [16, 17], such as branching process models [18, 19], optimal tidal chain outage models [20, 21], and manchester models [22]. These types of models consider the effect of operating state on cascading fault propagation and ignore the effect of topology. The second type is based on complex network theory, which considers the grid topology and evaluates the node criticality by degree and mesonumber. The main models are the small-world network model [23], scale-free network model [24], and traditional grid chain fault accident chain model. These models are more efficient but ignore the dynamic characteristics of the grid. Based on the two types of methods, Wu et al. presented a comprehensive index for the overall evaluation of the performance of the power system in terms of both grid structure and state, combining line meshes and tidal entropy [25]. Li et al. proposed a method to evaluate the importance of grid nodes based on the PageRank algorithm. Based on the establishment of a directed graph, the importance of node loads, load capacity, and network topology are considered comprehensively [26]. Zhou et al. proposed a power LeaderRank method for identifying important nodes in complex grids by combining the uncertainty of renewable energy, system topology, transmission flow, and load loss due to cascading faults [27].The above studies do not take geomagnetic storms into account. Based on this, we propose an improved LeaderRank algorithm that combines the grid topology with the system flow to accurately identify the critical nodes of the grid under the action of geomagnetic storms and to identify the set of possible accident chains caused by node failures. In this paper, we mainly accomplish the following work:(1)
An improved LeaderRank algorithm has been proposed which considers the interactions between nodes and represents the background nodes as the effects of geomagnetic storms. Based on this method, the critical nodes of the IEEE RTS79 system under the influence of geomagnetic storms have been analyzed.(2)
Based on the complex system fragility theory, the search method of accident chain under the condition of geomagnetic storm has been proposed. The IEEE RTS79 system is simulated to verify the validity and feasibility of the proposed accident chain identification method.The rest of the paper is organized as follows. Section2 of this paper presents a model for calculating GIC and transformer reactive power losses in power systems and discusses the effects of variations in the angle and strength of the geomagnetic field on GIC and transformer reactive power losses. Section 3 proposes an improved LeaderRank algorithm and identifies the critical nodes under different geomagnetic storm conditions after verifying its effectiveness. Section 4 simulates the fault development process of the IEEE RTS79 test system based on the vulnerability theory of complex systems and obtains the system accident chain under different geomagnetic storm conditions. A summary and conclusions are presented in Section 5.


## 2. GIC and Transformer Reactive Power Loss Calculation

The calculation of GIC is divided into two steps [28]. First, the plane wave method is used to solve the induced geoelectric field based on the monitoring data of the ground geomagnetic station. Second, the induced geoelectric field is equated to the power supply in a power system and the system GIC is solved by the nodal admittance method.According to Faraday’s law of electromagnetic induction, the line-induced voltage is the integral of the electric field along the line. When the induced geoelectric field is considered to be a constant, the value depends on the geographic location at both ends of the line. The induced voltage values for the line in the north-south and east-west directions can be calculated separately and then superimposed, represented as(1)U=ENLN+EELE,where EE and EN are the east and north geoelectric field values (V/km), respectively, and LN and LE are related to the actual geographic location information of the two substations. With reference to the Earth’s WGS-84 model [29], LN and LE can be presented as(2)LN=111.103−0.56cos2φ⋅Δlat,LE=111.5056−0.1872cos2φ⋅Δlong,where φ is the average of the latitude of the transformers connected at both ends of the transmission line. Δlat and Δlong are the latitude and longitude differences between substations.According to the nodal admittance method, for any nodei, we have(3)Ji=uiyi+∑k=1nuiyki−∑k=1nukyki,where Ji is the sum of the current sources of the branch connected to i, yi is the ground admittance of node i, and yki is the branch admittance between node k and node i. By calculating the node voltages ui and uk, the GIC of the substation into the ground can be obtained.The IEEE RTS79 test system is used to calculate GIC. The system consists of 10 generators, 24 buses, and 38 transmission lines. The topology is shown in Figure1.Figure 1
IEEE RTS79 topology diagram.Since the GIC amplitude is the same when the electric field angles differ by 180° [30], 0° is set as the east direction. Four typical angles of 0°, 45°, 90°, and 135° are chosen as the electric field angles. The substation GIC for different electric field angles at an electric field strength of 1 V/km is given in Figure 2. Table 1 gives the results of GIC calculations when the electric field angle is 0°.Figure 2
GIC of substation under the action of the geoelectric field at different angles. (a) E field 0°, (b) E field 45°, (c) E field 90°, and (d) E field 135°.
(a)(b)(c)(d)Table 1
GIC of substations under the action of geoelectric fields of different strengths.
NameEastward E field GIC (A)E = 1 V/kmE = 3 V/kmE = 5 V/kmE = 9 V/kmSUB1−7.529−22.588−37.647−67.764SUB2−7.616−22.849−38.082−68.547SUB37.36622.09736.82966.292SUB4−1.435−4.305−7.176−12.916SUB5−2.006−6.017−10.027−18.049SUB61.9675.8999.83217.698SUB7−5.374−16.123−26.871−48.368SUB8−5.570−16.709−27.849−50.127SUB9−0.060−0.181−0.302−0.544SUB103.61010.82918.04932.487SUB112.8558.56614.27725.698SUB12−7.893−23.679−39.464−71.036SUB13−0.367−1.101−1.835−3.303SUB140.0000.0000.0000.000SUB156.39519.18431.97457.552SUB160.4991.4972.4944.489SUB172.3897.16611.94421.499SUB185.43716.31127.18548.934SUB192.1516.45410.75619.361SUB205.18315.54825.91346.643Combining the data in Figure2 and Table 1, it can be seen that the substation GIC varies irregularly from the angle at a determined electric field strength. The substation GIC varies linearly from the strength at a determined geoelectric field angle. In other words, the relative magnitude of GIC of each substation is fixed.The GIC is linearly related to the reactive power loss of the transformer caused by it, as shown in the following equation [31, 32]:(4)QGIC=KVpuIGIC,where QGIC is the equivalent reactive power loss generated by the GIC flowing through the transformer winding, K is the reactive power loss coefficient of different types of transformers obtained through experimental measurements, Vpu is the actual operating end voltage of the transformer, and IGICis the single-phase GIC flowing through the neutral point of the transformer.Figure3 shows the effect of variation in the angle and strength of the geoelectric field on the reactive power loss in the substation.Figure 3
Calculation results of reactive power loss in the substation. (a) The effect of angle on reactive power loss and (b) the effect of strength on reactive power loss.
(a)(b)From Figure3, it is easy to find that the reactive power loss of the substation has a strong sensitivity to the geoelectric field angle. At a certain geoelectric field angle, the reactive power loss of the substation increases with the enhancement of the geoelectric field strength, which is consistent with (4).


## 3. Improved LeaderRank Algorithm considering the Influence of Geomagnetic Storm

### 3.1. Improved LeaderRank Algorithm Model

Based on the PageRank algorithm, the LeaderRank algorithm adds a background node. This node establishes a bidirectional connection with other nodes in the network, solving the problem of nodes jumping to each other in the PageRank algorithm. While ensuring the strong connectivity of the network, the convergence speed of the algorithm is accelerated. The iterative formula of the algorithm can be presented as(5)LRit+1=∑j=1N+1aijkjLRjt,where the nodes i and j are connected, aij=1; otherwise, it is 0. kj is the out-degree of node j and N is the number of nodes. The initial value of the background node is 0, and the initial values of the remaining nodes are all 1. After several iterations, the stable value of each node is finally obtained. The stable values of network nodes and background nodes are LRi and LRg, respectively. The final value of each node is obtained by(6)LRi=LRit+LRgtN.In the above iterative process, the LR value of each node is uniformly assigned to other outgoing chain nodes. However, applying it to a real power system can cause a large error. In addition, the grid topology and system flow distribution affect the node importance. The nodes that are more closely connected to other nodes and transmit more power are more important. Therefore, in this paper, an improved LeaderRank algorithm is proposed according to the characteristics of power networks. The effects of node transmission power, network topology, and geomagnetic storms are considered comprehensively to accurately identify the critical nodes of the system.The linking behavior between nodes can be represented by the probabilistic transition matrixM. The improved LeaderRank iterative computation is essentially the computation of the transition matrix M. Equations (5) and (6) can be rewritten as(7)LRt+1=M×LRt.The link matrixM is represented as follows:(8)Mij=1Xij1Di×Pij∑j=1NPij+1Dj×Pij∑i=1NPij,0<i,j<N,i≠j,0,i=j,Mij=1Xave1Dg×Qj−GIC∑j=1NQGIC,0<j<N,i=N+1,1Xave1Dg×Qi−GIC∑i=1NQGIC,0<i<N,j=N+1,where Xij is the line reactance between the nodes i and j, Xave is the average line reactance, Di and Dj are the output degrees of nodes i and j, respectively, Dg is the background node output degrees, Pij is the line power between the nodes i and j, and Qi−GIC and Qj−GIC are the reactive power loss of the transformer at nodes i and j, respectively.By normalizing the elements of theM matrix, the new matrix obtained is an irreducible nonperiodic matrix with random moments. The M matrix with this property ensures that the improved LeaderRank algorithm must converge.


### 3.2. Example Analysis

In this section, the effectiveness of the critical node identification method based on the improved LeaderRank algorithm is verified, and the critical nodes of the IEEE RTS79 system under the action of different geomagnetic storms are identified.The complex network method is used to identify critical nodes. We refer to the node electrical centrality (NEC) from [33] and add the geomagnetic storm parameter considering the effect of geomagnetic storms. The NEC considering geomagnetic storms can be presented as(9)NECi=μ1Bei+μ2Ei+1−μ1−μ2Gi,Gi=Qi−GIC∑j=1NQj−GIC,where Bei, Ei, and Gi are the electrical betweenness, eigenvector centrality, and geomagnetic storm parameters, respectively. μ is the distribution coefficient indicating the weight of the indicators. i and j are the nodes of the network, and N is the number of nodes in the whole network. The definitions of Bei and Eiare detailed in [33].Since the geomagnetic field strength is linearly related to the reactive power loss of the transformer and the geomagnetic field strength does not change the node size order, only the effect of the change of the electric field angle on the critical nodes is considered. The above two methods are used to identify the nodes of the IEEE RTS79 system under the action of eastward and northward geoelectric fields, and the top 10 nodes are selected as critical nodes. The results are shown in Table2.Table 2
Critical node order.
RankEastward E fieldNorthward E fieldNECLRNECLR115181623218151822321242318411111516524232215612121724723111138222420912114111014142017Table2 shows that the critical nodes obtained by both the improved LeaderRank algorithm and the complex network method are the same and differ only in the order of the nodes, which indicates that the critical nodes in the power network under geomagnetic storm conditions can be effectively identified by using the improved LeaderRank algorithm. The results of the node LR calculations for a typical geoelectric field angle are given in Figure 4, and the substations where the critical nodes are located are marked in Figure 5.Figure 4
Node LR values.Figure 5
Critical node locations.Figure4 shows that geomagnetic storms change the importance of the original power system nodes. The change of the geoelectric field angle caused the change in the importance of the power system nodes. Comparing Figures 2 and 5, it can be found that the GIC of the substation is larger regardless of the change in the geoelectric field angle, which indicates that the larger the GIC of the substation during the geomagnetic storm, the greater the possibility of its node becoming a critical node of the system.


## 3.1. Improved LeaderRank Algorithm Model

Based on the PageRank algorithm, the LeaderRank algorithm adds a background node. This node establishes a bidirectional connection with other nodes in the network, solving the problem of nodes jumping to each other in the PageRank algorithm. While ensuring the strong connectivity of the network, the convergence speed of the algorithm is accelerated. The iterative formula of the algorithm can be presented as(5)LRit+1=∑j=1N+1aijkjLRjt,where the nodes i and j are connected, aij=1; otherwise, it is 0. kj is the out-degree of node j and N is the number of nodes. The initial value of the background node is 0, and the initial values of the remaining nodes are all 1. After several iterations, the stable value of each node is finally obtained. The stable values of network nodes and background nodes are LRi and LRg, respectively. The final value of each node is obtained by(6)LRi=LRit+LRgtN.In the above iterative process, the LR value of each node is uniformly assigned to other outgoing chain nodes. However, applying it to a real power system can cause a large error. In addition, the grid topology and system flow distribution affect the node importance. The nodes that are more closely connected to other nodes and transmit more power are more important. Therefore, in this paper, an improved LeaderRank algorithm is proposed according to the characteristics of power networks. The effects of node transmission power, network topology, and geomagnetic storms are considered comprehensively to accurately identify the critical nodes of the system.The linking behavior between nodes can be represented by the probabilistic transition matrixM. The improved LeaderRank iterative computation is essentially the computation of the transition matrix M. Equations (5) and (6) can be rewritten as(7)LRt+1=M×LRt.The link matrixM is represented as follows:(8)Mij=1Xij1Di×Pij∑j=1NPij+1Dj×Pij∑i=1NPij,0<i,j<N,i≠j,0,i=j,Mij=1Xave1Dg×Qj−GIC∑j=1NQGIC,0<j<N,i=N+1,1Xave1Dg×Qi−GIC∑i=1NQGIC,0<i<N,j=N+1,where Xij is the line reactance between the nodes i and j, Xave is the average line reactance, Di and Dj are the output degrees of nodes i and j, respectively, Dg is the background node output degrees, Pij is the line power between the nodes i and j, and Qi−GIC and Qj−GIC are the reactive power loss of the transformer at nodes i and j, respectively.By normalizing the elements of theM matrix, the new matrix obtained is an irreducible nonperiodic matrix with random moments. The M matrix with this property ensures that the improved LeaderRank algorithm must converge.


## 3.2. Example Analysis

In this section, the effectiveness of the critical node identification method based on the improved LeaderRank algorithm is verified, and the critical nodes of the IEEE RTS79 system under the action of different geomagnetic storms are identified.The complex network method is used to identify critical nodes. We refer to the node electrical centrality (NEC) from [33] and add the geomagnetic storm parameter considering the effect of geomagnetic storms. The NEC considering geomagnetic storms can be presented as(9)NECi=μ1Bei+μ2Ei+1−μ1−μ2Gi,Gi=Qi−GIC∑j=1NQj−GIC,where Bei, Ei, and Gi are the electrical betweenness, eigenvector centrality, and geomagnetic storm parameters, respectively. μ is the distribution coefficient indicating the weight of the indicators. i and j are the nodes of the network, and N is the number of nodes in the whole network. The definitions of Bei and Eiare detailed in [33].Since the geomagnetic field strength is linearly related to the reactive power loss of the transformer and the geomagnetic field strength does not change the node size order, only the effect of the change of the electric field angle on the critical nodes is considered. The above two methods are used to identify the nodes of the IEEE RTS79 system under the action of eastward and northward geoelectric fields, and the top 10 nodes are selected as critical nodes. The results are shown in Table2.Table 2
Critical node order.
RankEastward E fieldNorthward E fieldNECLRNECLR115181623218151822321242318411111516524232215612121724723111138222420912114111014142017Table2 shows that the critical nodes obtained by both the improved LeaderRank algorithm and the complex network method are the same and differ only in the order of the nodes, which indicates that the critical nodes in the power network under geomagnetic storm conditions can be effectively identified by using the improved LeaderRank algorithm. The results of the node LR calculations for a typical geoelectric field angle are given in Figure 4, and the substations where the critical nodes are located are marked in Figure 5.Figure 4
Node LR values.Figure 5
Critical node locations.Figure4 shows that geomagnetic storms change the importance of the original power system nodes. The change of the geoelectric field angle caused the change in the importance of the power system nodes. Comparing Figures 2 and 5, it can be found that the GIC of the substation is larger regardless of the change in the geoelectric field angle, which indicates that the larger the GIC of the substation during the geomagnetic storm, the greater the possibility of its node becoming a critical node of the system.


## 4. Accident Chain Identification Based on the Fragility Theory of Complex Systems

### 4.1. An Accident Chain Model Based on the Fragility Theory of Complex Systems

The scale and complexity of power systems are increasing and the whole power system can be considered as a complex system. The cascading failure of a power system can be considered as a chain reaction of power outages caused by a node or line failure, which can be regarded as a fragile transmission process [34, 35]. Suppose the system blackout accident chain is(10)L=T1,T2,……,Tn,where Ti=i=1,2,…,n is the element that constitutes the accident chain, T1 is the initial failure of the accident chain, and T2,…,Tn is the intermediate link. The system blackout accident is caused by Tn−1 induced Tn.If any event does not occur in the accident chain, the accident chain will not be generated. Therefore, the blackout accident can be avoided by predicting the accident chain in advance.The power system accident chain search is divided into two parts, namely, initial failure selection and intermediate link setting.


#### 4.1.1. Initial Failure Selection

The LR value of the line is determined by the LR value of the head-end node. The LR of the line is defined as(11)LRa,b=LRa+LRb,where a and b represent the power system nodes.


#### 4.1.2. Intermediate Link Setting

In the power system, other lines will be affected by the initial failure, which depends on the fragility relationship between the other lines in the system and the brittle source. The fragility correlation between lines and line overload coefficientλ is defined as(12)Iij=EijPi−EijPiμj,μj=Pj.max−Pj.0Pj.max,λi=EijPiμj,where Eij is the power change of line j after the fault of line i, Pi is the active power borne by the initial faulted line itself, Pj⋅0 is the active power of line j in actual operation, and Pj⋅max is the thermal stability limit of the active power of line j. μj is the thermal stability margin of other circuits.The propagation process of grid cascade faults during geomagnetic storms is as follows: firstly, calculate the line LR value, and set the largest line as the initial fault; secondly, take the initial fault line out of operation to determine whether the power flow has converged. If it converges, calculate the overload coefficient of the remaining line. If there is an overload coefficientλ<0, set the line as the next-level fault; otherwise, calculate the remaining line fragility correlation I, and select the line with the largest fragility correlation I as the next-level fault. Until the system power flow does not converge, stop searching and generate an accident chain. The process is shown in Figure 6.Figure 6
Accident chain search process.


### 4.2. IEEE RTS79 System Accident Chain Identification

Any geoelectric field can be composed of superimposed eastward and northward geoelectric fields. The total reactive power loss of the substation is calculated to be the smallest for a geoelectric field angle of 0° and the largest for a geoelectric field angle of 90°. For 0° and 90° to study the accident chain under different geoelectric field strengths. The top ten lines of line LR values are selected as the initial faults and the accident chains that may lead to system outages are found according to the flow in Figure4. The results are shown in Table 3.Table 3
Accident chain collection.
RankEastward E fieldE = 0E = 1 V/kmE = 3 V/kmE = 5 V/kmE = 9 V/km1L28-L24-L7-L11L27-L6L27-L6L27L272L30-L28-L24-L7-L11L32-L11L32-L11L32-L11L32-L113L23-L29-L7L25-L28-L26L25-L28-L26L25-L28-L26L254L24-L28-L7L21-L23-L22-L11L21-L23-L22-L11L21-L23L215L31-L38L24-L28-L7L24-L28-L7-L11L24-L28L246L27-L6L30-L28-L24-L7-L11L30-L28-L24-L7-L11L30-L28-L24L30-L287L29-L23-L7L19-L23L19-L23L19-L23L198L32-L11L1-L11L1-L11L1-L11L1-L119L19-L23L18-L23-L7-L11L18-L23-L11L18-L23L18-L2310L7-L23-L27L22-L23-L21-L11L22-L23-L21-L11L22-L23L22RankNorthward E fieldE = 0E = 1 V/kmE = 3 V/kmE = 5 V/kmE = 9 V/km1L28-L24-L7-L11L22-L23-L21-L11L22-L23-L21-L11L22-L23-L21L222L30-L28-L24-L7-L11L36-L37-L29L36-L37-L29L36-L37L36-L373L23-L29-L7L24-L28-L7L24-L28L24-L28L244L24-L28-L7L27-L6L27-L6L27L275L31-L38L31-L38L31-L38L31-L38L31-L386L27-L6L21-L23-L22-L11L21-L23-L22L21-L23-L22L217L29-L23-L7L30-L28-L24-L7-L11L30-L38-L24-L7L30-L38-L24-L7L30-L388L32-L11L28-L24-L7-L11L28-L24-L7-L11L28-L24-L7L289L19-L23L23-L29-L7L23-L29-L7L23-L29-L7L2310L7-L23-L27L18-L23-L7L18-L23L18-L23L18-L23Comparing the accident chains in Table3, it can be found that(1)
Different initial faults triggered different accident chains under the eastward and northward geoelectric fields. This is because the change of geoelectric field angle causes the change of GIC and transformer reactive power loss, which leads to the change of power system tide and triggers different accident chains.(2)
At the same geomagnetic field angle, the accident chain shortens with increasing geomagnetic field strength. This is because during weak geomagnetic storms, the power system is kept in balance by flow distribution, and as the strength of the magnetic storm increases, the power system flow imbalance causes the accident chain to shorten. This means that under strong magnetic storm conditions, the system is more likely to collapse or experience a major outage.(3)
L7/L11/L23/L24/L27/L28 are the vulnerable lines of the IEEE RTS79 system, where L7/L23/L24/L28 may not only be the initial failure of the accident chain but also propagate chain failures as intermediate links, which should be focused on prevention.


## 4.1. An Accident Chain Model Based on the Fragility Theory of Complex Systems

The scale and complexity of power systems are increasing and the whole power system can be considered as a complex system. The cascading failure of a power system can be considered as a chain reaction of power outages caused by a node or line failure, which can be regarded as a fragile transmission process [34, 35]. Suppose the system blackout accident chain is(10)L=T1,T2,……,Tn,where Ti=i=1,2,…,n is the element that constitutes the accident chain, T1 is the initial failure of the accident chain, and T2,…,Tn is the intermediate link. The system blackout accident is caused by Tn−1 induced Tn.If any event does not occur in the accident chain, the accident chain will not be generated. Therefore, the blackout accident can be avoided by predicting the accident chain in advance.The power system accident chain search is divided into two parts, namely, initial failure selection and intermediate link setting.


### 4.1.1. Initial Failure Selection

The LR value of the line is determined by the LR value of the head-end node. The LR of the line is defined as(11)LRa,b=LRa+LRb,where a and b represent the power system nodes.


### 4.1.2. Intermediate Link Setting

In the power system, other lines will be affected by the initial failure, which depends on the fragility relationship between the other lines in the system and the brittle source. The fragility correlation between lines and line overload coefficientλ is defined as(12)Iij=EijPi−EijPiμj,μj=Pj.max−Pj.0Pj.max,λi=EijPiμj,where Eij is the power change of line j after the fault of line i, Pi is the active power borne by the initial faulted line itself, Pj⋅0 is the active power of line j in actual operation, and Pj⋅max is the thermal stability limit of the active power of line j. μj is the thermal stability margin of other circuits.The propagation process of grid cascade faults during geomagnetic storms is as follows: firstly, calculate the line LR value, and set the largest line as the initial fault; secondly, take the initial fault line out of operation to determine whether the power flow has converged. If it converges, calculate the overload coefficient of the remaining line. If there is an overload coefficientλ<0, set the line as the next-level fault; otherwise, calculate the remaining line fragility correlation I, and select the line with the largest fragility correlation I as the next-level fault. Until the system power flow does not converge, stop searching and generate an accident chain. The process is shown in Figure 6.Figure 6
Accident chain search process.


## 4.1.1. Initial Failure Selection

The LR value of the line is determined by the LR value of the head-end node. The LR of the line is defined as(11)LRa,b=LRa+LRb,where a and b represent the power system nodes.


## 4.1.2. Intermediate Link Setting

In the power system, other lines will be affected by the initial failure, which depends on the fragility relationship between the other lines in the system and the brittle source. The fragility correlation between lines and line overload coefficientλ is defined as(12)Iij=EijPi−EijPiμj,μj=Pj.max−Pj.0Pj.max,λi=EijPiμj,where Eij is the power change of line j after the fault of line i, Pi is the active power borne by the initial faulted line itself, Pj⋅0 is the active power of line j in actual operation, and Pj⋅max is the thermal stability limit of the active power of line j. μj is the thermal stability margin of other circuits.The propagation process of grid cascade faults during geomagnetic storms is as follows: firstly, calculate the line LR value, and set the largest line as the initial fault; secondly, take the initial fault line out of operation to determine whether the power flow has converged. If it converges, calculate the overload coefficient of the remaining line. If there is an overload coefficientλ<0, set the line as the next-level fault; otherwise, calculate the remaining line fragility correlation I, and select the line with the largest fragility correlation I as the next-level fault. Until the system power flow does not converge, stop searching and generate an accident chain. The process is shown in Figure 6.Figure 6
Accident chain search process.


## 4.2. IEEE RTS79 System Accident Chain Identification

Any geoelectric field can be composed of superimposed eastward and northward geoelectric fields. The total reactive power loss of the substation is calculated to be the smallest for a geoelectric field angle of 0° and the largest for a geoelectric field angle of 90°. For 0° and 90° to study the accident chain under different geoelectric field strengths. The top ten lines of line LR values are selected as the initial faults and the accident chains that may lead to system outages are found according to the flow in Figure4. The results are shown in Table 3.Table 3
Accident chain collection.
RankEastward E fieldE = 0E = 1 V/kmE = 3 V/kmE = 5 V/kmE = 9 V/km1L28-L24-L7-L11L27-L6L27-L6L27L272L30-L28-L24-L7-L11L32-L11L32-L11L32-L11L32-L113L23-L29-L7L25-L28-L26L25-L28-L26L25-L28-L26L254L24-L28-L7L21-L23-L22-L11L21-L23-L22-L11L21-L23L215L31-L38L24-L28-L7L24-L28-L7-L11L24-L28L246L27-L6L30-L28-L24-L7-L11L30-L28-L24-L7-L11L30-L28-L24L30-L287L29-L23-L7L19-L23L19-L23L19-L23L198L32-L11L1-L11L1-L11L1-L11L1-L119L19-L23L18-L23-L7-L11L18-L23-L11L18-L23L18-L2310L7-L23-L27L22-L23-L21-L11L22-L23-L21-L11L22-L23L22RankNorthward E fieldE = 0E = 1 V/kmE = 3 V/kmE = 5 V/kmE = 9 V/km1L28-L24-L7-L11L22-L23-L21-L11L22-L23-L21-L11L22-L23-L21L222L30-L28-L24-L7-L11L36-L37-L29L36-L37-L29L36-L37L36-L373L23-L29-L7L24-L28-L7L24-L28L24-L28L244L24-L28-L7L27-L6L27-L6L27L275L31-L38L31-L38L31-L38L31-L38L31-L386L27-L6L21-L23-L22-L11L21-L23-L22L21-L23-L22L217L29-L23-L7L30-L28-L24-L7-L11L30-L38-L24-L7L30-L38-L24-L7L30-L388L32-L11L28-L24-L7-L11L28-L24-L7-L11L28-L24-L7L289L19-L23L23-L29-L7L23-L29-L7L23-L29-L7L2310L7-L23-L27L18-L23-L7L18-L23L18-L23L18-L23Comparing the accident chains in Table3, it can be found that(1)
Different initial faults triggered different accident chains under the eastward and northward geoelectric fields. This is because the change of geoelectric field angle causes the change of GIC and transformer reactive power loss, which leads to the change of power system tide and triggers different accident chains.(2)
At the same geomagnetic field angle, the accident chain shortens with increasing geomagnetic field strength. This is because during weak geomagnetic storms, the power system is kept in balance by flow distribution, and as the strength of the magnetic storm increases, the power system flow imbalance causes the accident chain to shorten. This means that under strong magnetic storm conditions, the system is more likely to collapse or experience a major outage.(3)
L7/L11/L23/L24/L27/L28 are the vulnerable lines of the IEEE RTS79 system, where L7/L23/L24/L28 may not only be the initial failure of the accident chain but also propagate chain failures as intermediate links, which should be focused on prevention.


## 5. Conclusion

In this paper, considering the impact of geomagnetic storms on power systems, a critical node identification method is proposed by integrating topology and system flow. A chain fault accident chain simulation model is established based on the complex system fragility theory. The model selects the initial fault lines based on the critical nodes, predicts the intermediate links by calculating the line overload coefficients and correlation degrees, and finally obtains the accident chain set. The feasibility and effectiveness of the method are verified by using the IEEE RTS79 system as an example. The obtained research results are as follows:(1)
During geomagnetic storms, critical nodes are mostly concentrated in substations with large GICs and their strong sensitivity to the angle of the geoelectric field.(2)
As the strength of magnetic storms increases, the accident chain of the power system shortens, meaning that major outages are more likely to occur.(3)
The intersections of the accident chains are the weak link of the power system. For IEEE RTS79 systems, appropriate measures should be taken and focused on prevention for L7/L11/L23/L24/L27/L28.



---

*Source: 2804781-2022-02-09.xml*