Theory and Method of Power System Integrated Security Region Irrelevant to Operation States: An Introduction

Yixin Yu*, Yanli Liu, Chao Qin, Tiankai Yang

Key Laboratory of Smart Grid of Ministry of Education, Tianjin University, Tianjin 300072, China

Keywords:Security region Composition Power flow security Static voltage stability small-disturbance stability Transient stability Dynamics nature Topological and geometric characteristics Hyperplane

A B S T R A C T How to comprehensively consider the power flow constraints and various stability constraints in a series of power system optimization problems without affecting the calculation speed is always a problem.The computational burden of probabilistic security assessment is even more unimaginable. In order to solve such problems, a security region (SR) methodology is proposed, which is a brand-new methodology developed on the basis of the classical point-wise method. Tianjin University has been studying the SR methodology since the 1980s, and has achieved a series of original breakthroughs that are described in this paper. The integrated SR introduced in this paper is mainly defined in the power injection space,and includes SRs to ensure steady-state security, transient stability, static voltage stability, and smalldisturbance stability.These SRs are uniquely determined for a given network topology(as well as location and clearing process for transient faults) and given system component parameters, and are irrelevant to operation states. This paper presents 11 facts and related remarks to introduce the basic concepts, composition, dynamics nature, and topological and geometric characteristics of SRs. It also provides a practical mathematical description of SR boundaries and fast calculation methods to determine them in a concise and systematic way. Thus, this article provides support for the systematic understanding, future research,and applications of SRs.The most critical finding on the topological and geometric characteristics of SRs is that, within the scope of engineering concern, the practical boundaries of SRs in the power injection space can be approximated by one or a few hyperplanes. Based on this finding, the calculation time for power system probabilistic security assessment (i.e., risk analysis) and power system optimization with security constraints can be decreased by orders of magnitude.

1. Introduction

For a long time, the classical method used to analyze the security and stability of power transmission and distribution systems[1,2] is the point-wise method, i.e., the conclusion that a system is secure/unsecure or stable/unstable is drawn by simulation calculation in one or more fault modes of a specific scenario (operation mode; i.e., operation point). This method still plays an important role in power system analysis. However, it is difficult to use this method to quickly and directly obtain an overall evaluation of power system operation states,such as how far the current operation point is from the stability boundary, how large the stability margin is,and so on.How to comprehensively consider power flow constraints and various stability constraints in a series of power system optimization problems without affecting the decisionmaking speed is always a difficult issue.The computational burden of power system probabilistic security assessment is even more unimaginable.

Wu et al. [3] introduced the concept of probabilistic security assessment and Refs. [4–6] put forward the concept of steadystate security region (SSSR) and dynamic security region (DSR) in the injection space. Security region (SR) methodology is a new methodology developed on the basis of the point-wise method. It addresses problems from the perspective of region and describes a region in which a power system can operate securely as a whole(the green area in Fig.1[7]is the SR).The relative relation between system operation points and SR boundary can provide information on the security margin and optimal control direction, which will enable power system operators to perform online real-time security monitoring, assessment, and control in a more scientific and efficient manner. The SRs introduced in this paper are mainly defined in the power injection space (or decision space). They are uniquely determined for a given network topology and given system component parameters (as well as location and clearing process of faults for DSR), and do not vary with operation states.Therefore, it is only necessary to calculate the SR boundaries (i.e.,hyperplane (HP) coefficients, etc.) once and store them in a database for use in future analysis and calculation,without adding calculation burden to online application. In view of this advantage,this paper only introduces the SRs that are irrelevant to the operation states and does not cover the method of determining SR boundaries based on operation states obtained through real-time measurement.

Fig. 1. Schematic diagram of a SR [7].

The mathematical description[3]of power system security and stability problem not only includes electromechanical dynamics,but also has ultrahigh dimensional and nonlinear properties,which cause difficulties in imagining and grasping the topological and geometric characteristics of the SRs.

Regarding the study of the topological nature of SRs, many properties have been clearly understood through previous extensive and in-depth research on the power flow stability region(with the boundary of saddle-node bifurcation (SNB)) [8] and the smalldisturbance stability region (of the system expressed by differential algebraic equations)[9]defined in the state space.For example,each of the stability regions is not necessarily connected; within a restricted domain of interest (compact set), the number of nonconnected parts in the stability region is finite;and there are some sufficient conditions to ensure that the power flow solution of the power system is unique under normal conditions, and is stable within the power flow stability region. Here, the so-called state space is composed of node voltage magnitudes and relative angles between node voltages. In order to facilitate engineering application, it is preferable to give SRs the following characteristics:They are defined in the power injection space (decision space); for a given network topology (as well as location and clearing process for transient faults)and given system component parameters,each SR is connected and uniquely determined, and the SRs are irrelevant to operation states; there is no void inside them; and the SRs’ boundary surfaces are piecewise smooth.

Early representative research on the geometric characteristics of SRs has been presented in Refs. [6,8–10]. These works propose the use of an internal truncated super-cuboid of the static SR in order to approximate the static SR. The internal truncated supercuboid approximation presented in Refs. [6,8–10] is very convenient to use, but is very conservative. Often, they cannot include some of the secure operation points of interest. Therefore, it is necessary to develop an SR expression with higher precision and more convenient for power system security analysis, assessment,and control.

At the same time, there is a lack of systematic research on the application methodologies of SRs, especially in the areas of probabilistic security assessment (risk analysis) and the optimization of power systems with security constraints.

Tianjin University has been studying SR methodology since the 1980s [7,11–14], and has made systematic research progress,including: insight into the dynamics nature and topological and geometric characteristics of SRs; an approximate HP description approach of SR boundaries and corresponding fast calculation methods; and encouraging applications of SRs in practical largescale power systems. The research has been supported by the National Natural Science Foundation of China, the Research Fund for the Doctoral Program of Higher Education of China, and the US Electric Power Research Institute project. Some achievements have been extensively cited by Ref.[15].The purpose of the present paper is to provide a concise and systematic introduction of these achievements and offer their original sources.

This paper is organized as follows: Section 2 takes the steadystate operation of the grid as an example, and then explains what the so-called SR method is,in comparison with the classical pointwise method.It further defines the mathematical description space of SRs. This paper mainly uses the injection space [3] (although it sometimes also uses the decision space) as the mathematical description space of SRs, where each point (each vector) is composed of all independent variables in the system under certain assumptions. Defining the mathematical description space of the SR is one of the most important conditions to ensure that the SR under study is uniquely determined and connected. At the same time, the so-called critical cut-set space that dispatchers usually care about (although in such a case, the SR is not uniquely determined) is also discussed.

Section 3 introduces the composition of an integrated SR(IGSR)and its dynamics nature.Two basic assumptions(Hypothesis 1 and Hypothesis 2)are given and will be used later,and two facts(Fact 1 and Fact 2)are further given,based on the two assumptions;these are the basis for ensuring that the various SRs described in this paper are uniquely determined and connected. The IGSR of the power system in this paper is further defined as the intersection set of the following four SRs: the SSSR to ensure power flow security; the DSR to ensure transient stability; the SR to ensure static voltage stability (SVSR); and the SR to ensure small-disturbance stability (SDSR). These SRs are defined one by one, and some dynamics theories related to SR research are briefly introduced.Fact 3 is derived from Refs. [16–18], while Facts 4–7 are new discoveries to meet the needs of SR research. Facts 3–7 construct a theoretical foundation for subsequent research on the topological and geometric characteristics of SRs.

Section 4 is the core of this paper.Four important facts(Facts 8–11) are proposed in this section, which are used to describe the topological properties(e.g.,the SR is uniquely determined and connected; the SR does not change with operating states; there is no void inside the SR; the SR boundary consists of a finite number of smooth sub-planes) and the geometric characteristics (e.g., the composition of the SR or its boundary, and the property that the boundary of the SR can be approximated by one or a few HPs in the power injection space,within the acceptable range of the engineering application) of the SRs—namely, the SSSR, DSR, SVSR, and SDSR—in practical power systems. These facts are concluded from the research results obtained by our research group over a period of many years. The details of and references for the four facts are provided in the remarks,which also involve fast calculation methods for practical SR boundaries.

Based on the HP boundary of the SRs, new methods to address the power system optimization problem and probabilistic security analysis are systematically established in Section 5.By formulating all the security constraints as the linearly combined inequalities of decision variables (i.e., node power injections) in the objective functions,these methods have two major benefits:①A large class of power system optimization problems with security constraints,which used to be difficult to be solved, are greatly simplified; and②the n-fold integration problem of the probability density function of n-dimension variables in the probabilistic security assessment problem can be mathematically converted into the threshold comparison problem of a one-dimensional (1D) probability distribution function. Thus, the online computational speed of these two types of problems can be greatly improved by orders of magnitude.In order to provide readers with the concepts of the computational burden of the SR itself and the computational burden that can be saved in practical application,an application example is given in this section.This section also demonstrates that it is easy to realize the visualization of an SR and quickly determine the security margins based on HP boundaries of an SR, providing a powerful tool for the situation awareness of power systems. The paper is concluded in Section 6.

2. What is the SR methodology?

The SR methodology will be explained by taking the steadystate operation of power systems as an example. For this purpose,the point-wise method [7] needs to be presented first.

2.1. The point-wise method

Let a power system be composed of n + 1 nodes and nbbranches,with 0 to ngrepresenting generator nodes,0 representing the reference node, and ng+ 1 to n representing load nodes. Let G ?｛0,1,2,...,ng｝ represent the set of generator nodes,L ?｛ng＋1,...,n｝ represent the set of load nodes,N ?｛0,1,2,...,n｝ represent the set of all nodes,B ?｛1,2,...,nb｝ represent the set of all branches, and black italics represent vectors.Hence,the power flow equations for a power system are as follows:

where Piand Qidenote the active and reactive power injections of node i, respectively; Viand Vjare the voltage magnitudes of node i and node j;θiand θjare the voltage phase angles of node i and node j;θijis the branch angle,which can be calculated by θij=θi–θj;and Gijand Bijrepresent the real part and the imaginary part of the component at the ith row and jth column of the node admittance matrix. Each node has four variables: Pi, Qi, Vi, and θi. Therefore,the power flow equations involve 4n + 4 variables; but only 2n + 2 equations can be formulated according to Eqs. (1a) and(1b). To make the power flow equations solvable, in general, two variables are specified for each node, and the other two variables are to be calculated.For example,for a load node i,Piand Qiare usually given,while Viand θiare to be calculated.Considering that only the branch angle θijis involved in the power flow equations, there are only 2n + 1 independent variables in Eqs. (1a) and (1b): That is, the voltage phase angle of one node in the system can be specified arbitrarily,and if and only if one of the voltage phase angles of nodes is specified (e.g., the voltage phase angle of the reference node is specified as θ0= 0), the others of the voltage phase angles of nodes can be determined.

In power transmission systems, the branch admittance Gij≈0 and branch angle θijare very small, making sinθij≈ θijand cosθij≈1.Eqs.(1a)and(1b)can be simplified as decoupled power flow equations, shown as follows:

In addition, there are some operational constraints, such as the voltage magnitude of each node, the current of each branch, the active power and reactive power of generators and loads, and the branch angle of each branch,which should be within certain limits,as shown in Eqs. (3a)–(3c):

The constrained power flow problem is to find solutions satisfying the equality constraint F（x）＝0 and the inequality constraint G（x）≤0.Figs.2(a)and(b)[7]are schematic diagrams showing feasible solutions and no feasible solutions, respectively, in a simple two-dimensional (2D) case of x=(x1,x2)T. In either diagram, the curve is the solution of F(x)＝F(x1,x2)＝0,and the rectangle shows the bounding range of xm1≤x1≤xM1and x2m≤x2≤x2Mthat corresponds to G(x) ≤0.

To deal with this kind of problem, the traditional variable separation method may be used,which includes the following four steps:

Fig. 2. Schematic diagrams of constrained power flow problems: (a) with a solution satisfying the security constraints; (b) without a solution satisfying the security constraints [7].

Step 3: Solve F（x）＝0 and obtain xα.

Step 4:Check whether xαmeets xmα≤xα≤xMα.If yes,xαis a feasible solution to the constrained power flow problem, and xβis identified as being steady-state secure.

Specifically, the constrained power flow calculation consists of the following four steps:

where R is the set of all real numbers,xβis the given variable vector of the power flow equations, xα1is the solution of the power flow equations, and xα2is the vector of the variables that can be calculated directly by xβand xα1based on the power flow equations.xα3is also the vector of the variables that can be calculated by xβand xα1. However, those variables are not included in the power flow equations, and xα3may be calculated by Eq. (5):

If Eqs.(7a)–(7c)can be met,the given xβ(point a in Fig.3[7])is identified as being steady-state secure; in contrast, if any constraint of Eqs. (7a)–(7c) is not met, the given xβ(point b in Fig. 3[7]) is identified as being steady-state unsecure. Since the given xβis only a point in the space where it is located, this kind of analysis method is called the point-wise method. This method is widely used in the steady-state security analysis of power systems.

2.2. The concept of SR

An SR[7]describes a region in the space of xβ,as shown by the shaded part in Fig. 3 [7]. An xβcorresponding to any point inside this region (e.g., point a) is secure when verified by the pointwise method.An xβcorresponding to any point outside this region(e.g.,point b)is unsecure when verified by the point-wise method.

2.3. The SR definition space

In the abovementioned power flow calculation, the space defined in Eq.(4a),where xβis located,is called the decision space.As the decision space is identical to the given variable space of the power flow calculation,it is usually used in power system situation awareness and in dispatchers’ decision-making. The power injection space,which is usually used in SR research,is the space where the following vector is located:

where P and Q are the vector of the active and reactive power injections, respectively.

The current injection space, which is also used to define the SR in some research, is the space where the following vector is located:

where IPis the vector of the node’s active current injections; IQis the vector of the node’s reactive current injections; IPiand IQiare the active and reactive injection currents of node i, respectively,corresponding to the scenario of Vi≈1,?i ∈N. In power system operation optimization,the SR defined in the power injection space shown in Eq.(8a)is often more convenient than that defined in the decision space, shown in Eq. (4a). The current injection space defined in Eq. (8b) is only used for qualitative research.

Fig. 3. Schematic diagram of an SSSR [7].

In the following paragraphs,the spaces defined by Eqs.(4a)and(8a) are further explained in combination with practical engineering application requirements:

(1) When the complex voltage （V0,θ0） of the reference node is specified to be constant—as is often done in power flow calculation—the power injection space and the decision space can be expressed with the following two vectors, respectively:

where PGis the vector of the generator’s active power injections;QGis the vector of the generator’s reactive power injections; PLis the vector of the load’s active power injections; QLis the vector of the load’s reactive power injections; VGis the vector of the generator buses’ voltage magnitudes.

In the following paragraphs of this paper,the spaces defined by Eqs.(8c)and(8d)are mainly used to define SRs.It should be noted that （V0,θ0） has been specified for them.

(2) If more variables are specified to be constants, the dimension of the SR definition space may be further reduced, which can facilitate the calculation, analysis, and visualization of SRs.However,it should also be noted that these SRs are under specified conditions. For example, for a high voltage alternating current(HVAC) power system, it can be assumed that the reactive power is locally balanced;thus,a change in the active power will have little effect on the voltage magnitude,and only the SR defined in the active power injection space needs to be studied. Then, as the transmission loss is very small and ∑ni＝0Pi≈0, only n active power injections in the grid are independent of each other. Hence, the injection space shown in Eq. (8c) can be transformed as follows:

Eq. (8e) is often used in SR research related to the transient power angle stability problem of large-scale transmission systems.

It should be emphasized that the solution of power flow equations will be uniquely determined for a given injected power vector only when the complex voltage of the reference node is specified(e.g., V0= 1, θ0= 0).

(3)Considering the physical characteristics of HVAC power systems and the characteristics of power flow equations—that is,that the power flow equations can be decoupled based on the relationship between the active power injections and the voltage phase angles of the nodes, and based on the relationship between the reactive power injections and the voltage magnitudes of the nodes—Eq. (8c) is further divided into the active power injection space and the reactive power injection space. Therefore, we can study the SSSR in the reactive power injection space when the active power injections or voltage phase angles of the nodes are specified, and we can study the SSSR in the active power injection space when the reactive power injections or voltage magnitudes of the nodes are specified.In this way,not only can the space dimensionality be reduced, but also some very clear physical concepts can be obtained.

It should be noted that each x mentioned above refers to a vector composed of all independent variables in a system under certain assumptions, which is the basis for regarding the space where they are located as the space for defining an SR. In this paper, the spaces where the vectors shown in Eqs. (8a), (8c), and(8e) and the vectors with further dimensionality reduction are located will be collectively called the power injection space,while the space corresponding to Eq. (8d) will be called the decision space.

Furthermore, as power system dispatchers are also often concerned about the power transmission limits of certain critical sections(cut-sets)of power systems,in addition to the decision space and the power injection space, this paper also adopts the critical cut-set space as needed. The so-called ‘‘cut-set” refers to the minimum branch set for cutting a connected graph into two subgraphs.In subsequent paragraphs,there will be two typical critical cut-sets,as follows:①the critical cut-set in transient power angle stability research, which refers to the critical cut-set of a pre-fault system consisting of the critical cut-set of the post-fault network(i.e.,the system power angle splitting cut-set)and the branch that is cut due to protective action(on the premise that these two parts can form a cut-set);and ②the critical cut-set for describing an SR to ensure static voltage stability in the cut-set power space(CVSR),which refers to a complete cut-set in a system that splits the system into two disconnected parts,i.e.the part with voltage stability weak nodes and the part without voltage stability weak nodes, as shown in Fig. 4 [7]. In transient stability research, Region 1 indicates a non-critical node group and Region 2 indicates a critical node group.In static voltage stability research,Region 1 indicates a non-weak node group of system and Region 2 indicates a weak node group of system.It should be noted that there might be multiple critical cut-sets in a large-scale power system.

Since the definition variables of an SR in a power injection space and those of an SR in a decision space are controllable,such regions will bring great convenience to power system optimization,probabilistic security analysis, and risk assessment (as described below). However, from the perspective of power system monitoring and dispatching, the dispatchers of interconnected systems are particularly likely to monitor the transmission power on several sections (cut-sets) of a system, because the dimension of those cut-sets is very low and it is clear at a glance.

3. Composition and dynamics nature of SRs

The following hypotheses will be adopted in the subsequent paragraphs of this paper:

Hypothesis 1.The voltage magnitude and phase angle of the reference node are specified; that is, V0≈1,θ0＝0.

Hypothesis 2.During research on SRs in the space where the abovementioned xβis located, only consider the region enclosed by the SR boundaries encountered for the first time during a slow outward continuous extension from the initial normal operation point xoβin a quasi-static form.

Fact 1.Under Hypothesis 1, there is a one-to-one correspondence between the system steady operation states（xα1,xα2）and xβ.

Fact 2.Under Hypothesis 2, in the space of xβ, the boundary of the SR is uniquely determined and connected for a given network topology and given system component parameters.

As shown in Fig. 5 [18], the IGSR(symbolized by Ω)of a power system is the intersection set of the SSSR (symbolized by ΩSS) to ensure power flow security, the SVSR (symbolized by ΩSV)to ensure static voltage stability, the SDSR (symbolized by ΩSD)to ensure small-disturbance stability, and the DSR (symbolized by Ωd) to ensure transient stability; that is

Ω ?ΩSS∩ΩSV∩ΩSD∩Ωd（9）

Fact 3.If,and only if,the operation point xβis within Ω,the system is secure [17,18].

3.1. SSSR (ΩSS) to ensure power flow security

For a given network topology and given system component parameters,ΩSSis the set of all xβthat can satisfy the power flow equation F（x）＝0and the constraint G （x ）≤0in the power injection space or the decision space. It is expressed as follows:

Fig. 4. Schematic diagram of a critical cut-set [7].

Fig. 5. Schematic diagram of the IGSR (dotted grid region) [18].

3.2. SDSR (ΩSD) to ensure small-disturbance stability

For a given network topology and given system component parameters, the SDSR is the set of all operation points defined in the power injection space, each of which can ensure the smalldisturbance stability of the power system.

Power system dynamic models can be generally represented by the following differential-algebra equations (DAEs) [7,19–21]:

where xs∈Rmstateis the state variable vector; xa∈Rlalgeis the algebraic variable vector;and xβ:＝（PTG,QTG,PTL,QTL）T∈R2nis the control variable vector (as given in Eq. (8c)). For a given xβ, the set of system equilibrium points (EPs)—that is, EPs（xβ）—can be defined as follows:

Note 1:Based on Eq.(10d)and nonlinear system theory,we can have the following important concepts:

(1) For the system defined by Eq. (10), when matrix D(xβ) is non-singular, if and only if all eigenvalues of matrix A～（xβ） have negative real parts, the system is small-disturbance stable at（xs0,xa0）; when matrix D（xβ） is non-singular and the eigenvalues of A～（xβ） vary continuously with xβ, if any real eigenvalue λ of A～（xβ） is changed from negative to positive, the point x′βcorresponding to λ ＝0 is called a SNB of the system and,after this point,the system will lose stability in a monotonic mode when there are some small disturbances.

(2) If the real parts of a pair of conjugate eigenvalues λ ＝a±jb （b ≠0） of D（xβ） change from negative to positive, the point x′βcorresponding to λ ＝±jb is called a Hopf bifurcation(HB)of the system defined by Eq.(10)and,after this point,the system will lose stability in a continuous oscillation mode with increasing magnitude when there are some small disturbances.

(3)The matrix D（xβ）is not always non-singular,and may be singular in some cases. If D（xβ） is singular, it will be impossible to eliminate Δxafrom Eq. (10c); thus, it will be impossible to obtain Eq. (10d). In this case, the point x′βthat makes D（xβ） singular is called as singularity-induced bifurcation (SIB) of the system defined by Eq. (10). At this point, D（xβ） has a real eigenvalue μ,which changes its sign at the SIB. When μ crosses the zero point,A～（xβ）will have an eigenvalue λ that changes its sign suddenly from infinity at one end to infinity at the other end (i.e., +∞→-∞or-∞→+∞). If matrix A～（xβ） has an eigenvalue λ that changes suddenly from-∞to+∞,the system will lose stability in a monotonic mode when there are some small disturbances.

Hence, in the power injection space (i.e., parameter space) R2n,the SR to ensure the small-disturbance stability of the system defined by Eq. (10) is defined as follows:

The research on the boundary of ΩSD(denoted by ?ΩSD) also involves the chaos phenomenon [19–22] in addition to the SNB,HB, and SIB. About the boundary of SDSR and chaos in power systems, Fact 4 was found.

Fact 4.Chaos occurs outside the bifurcation boundaries of the SDSR.

Remark:According to Ref. [21], the chaos phenomenon in power systems is formed by further development of singleperiod bifurcation, and is likely to be an intermediate stage of large-disturbance instability. Several ways of instability and collapse induced by chaos are shown in Fig.6.If(based on Hypothesis 2) the boundary of the SDSR is the boundary that is first encountered during outward extension (i.e., the gradual increase of the system load) from the normal operation point in the power injection space, the chaos phenomenon will occur only outside the HB boundary of the SDSR in power systems (i.e., chaos will not occur in advance of HB).Since HB is not allowed in the secure and stable operation of power systems, there is no need for further research on chaos.Then,if non-periodic instability caused by SNB or SIB corresponds to the boundary that is first encountered during outward extension from the normal operation point in the power injection space,and if no single-period oscillation of HB occurs,no chaos will occur.Therefore,chaos does not need to be investigated too much in SDSR research, thus greatly reducing the search range of ΩSD.

After excluding the consideration of chaos, it is also noted that the definition shown in Eq. (11) does not assume that D（xβ） is not singular,so the boundary of ΩSDcan be formulated as follows[23]:

Note 2:In recent years,static voltage stability has been of great concern in power systems.There are two definitions of voltage stability. One definition is given by the Conseil International des Grands Réseaux électriques(CIGRE)working group on voltage stability. This definition is based on the small-disturbance stability bifurcation theory, and thus belongs to the category of the abovementioned small-disturbance stability SR. The other definition,which is also commonly used,is given by the Institute of Electrical and Electronics Engineers(IEEE).Voltage stability is defined as the capability of a system to maintain its voltage. When the load admittance increases, the load power increases accordingly, and both the power and the voltage are controllable. As shown in Fig. 7 [7], for general ZIP static load models (where Z stands for constant impedance, I stands for constant current, and P stands for constant power),the SNB(i.e.,singular point of the power flow equation’s Jacobian matrix) of the system is located at the lower half branch of the P–V curve. An increase in load admittance will lead to a decrease in load power in the lower half branch of the P–V curve. According to the IEEE’s definition of voltage stability,the lower half branch belongs to the unstable region. Therefore,only the upper half branch of the P–V curve is considered, and the nose tip of the P–V curve (i.e., fold bifurcation) is taken as the critical point of the static voltage stability. This nose tip is also the system SNB point obtained in the constant power load model.Therefore,only the constant power load model needs to be considered in the calculation of the static voltage stability SR.

Fig. 6. Several ways of instability and collapse induced by chaos. Q1d: the reactive power load; KA: proportion coefficient of the excitation system; TA: time constant of the excitation system.

Fig. 7. Singular point on the P–V curve [7]. ①System P–V curve; ②characteristic curve of ZIP load model; ③characteristic curve of constant power load model.

Remark:Through bifurcation analysis and the two-step analysis method, the following conclusions are given in Ref. [24]:

(1) Considering the recoverable dynamic load while ignoring the dynamics of other components such as generators, the system voltage stability limit point obtained from small-disturbance analysis—that is, the SNB point—is the same as the fold bifurcation point obtained based on the CPF of the static constant power load model(the nose tip of the P–V curve).That is to say,the boundary of the power flow feasible region coincides with the SNB boundary of the SDSR. When the specific dynamic load and the dynamics of the generator and its regulating systems are considered,it is probable that the SNB point will occur in advance of fold bifurcation. In other words,the SNB points of DAEs might still exist in the power flow feasible region, i.e., the SR to ensure the existence of EPs.

(2)These SNB points will not necessarily cause voltage collapse in power systems. Not all SNB points are on the boundary of the SDSR. The nature of SNB points needs to be analyzed according to the specific conditions of the system. When induction motors(IMs) account for a large proportion of the load, it is essential to judge the nature of SNB points by means of a time-domain simulation when determining the boundaries of the SDSR.

Fact 6.The presence of IMs in the load might lead to no occurrence of HB in power systems before the increase of power injections causes voltage instability of the SNB type.

Remark:In theoretical research,HB generally occurs in advance of the SNB point. However, this type of voltage oscillation instability phenomenon is seldom discovered in the actual recording of voltage instability events. The voltage instability of actual power systems is often in a collapse mode with a monotonic voltage drop, which is closely related to SNB. In the opinion of some scholars, the occurrence of monotonic voltage collapse before reaching the HB of an unconstrained system might be caused by the operation limits in power systems. This is called limitinduced bifurcation (LIB); thus, no HB occurs on the SDSR boundary. However, in reality, many voltage instability events are not necessarily caused by the operation limits of the system.According to Ref.[24],the presence of IMs might lead to no occurrence of HB,and the voltage oscillation phenomenon resulting from HB before SNB-type voltage instability is caused by the increase of power injections. Therefore, when there is a high proportion of IM load in the system, SNB can be prioritized instead of HB.

Fig. 8 provides the voltage stability limit and the relationship between instability modes and the proportion of IMs in the load.On the one hand, the figure shows the particularity of the IM load model in voltage stability limit calculation, such that HB will not occur due to the increase in the proportion of IMs in the load. On the other hand, the figure indicates that there is a big difference between the result of CPF based on static models(fold bifurcation)and the more accurate result of small-disturbance voltage stability analysis based on dynamic models(SNB point)and the SNB points of the DAE system are still likely to exist in the power flow feasible region (i.e., the SR to ensure the existence of EPs).

3.3. DSR (Ωd) to ensure transient stability

The DSR [7,13] of a power system in power injection space is defined as follows:

where Fdis a given fault; idand jdare the pre-fault and post-fault network topology, respectively; xτ（xβ） represents the system state at the instant of fault clearing τ; and S（xβ） represents the stable region encircled by EPs xsin the post-fault injection space. Both S（xβ）and xsare determined by xβ.It should be noted that Ωddetermines a set of all points in the pre-fault power injection space that can ensure the transient stability of the power system.

In research on transient stability with direct methods, one injected power vector xβ(an operation point) will correspond to a transient stability region S（xβ） [19,20]. Fig. 9 [18] shows the distinction and correlation between the DSR Ωdand the transient stability region S（xβ）.

After the fault is cleared, if in the group of generators with relative accelerations that are equal to zero in the generator rotor swing equations and phase angles that are greater than π/2, there are some generators with a relative speed of changes in the phase angles that is still greater than zero,such generators will be out of synchronization with the others,and the power system will not be able to recover to synchronous operation.In general,the set of generators with phase angles greater than π/2 is called an unstable generator group, while the set of other generators is called a non-unstable generator group. The different classifications of the unstable generator group and the non-unstable generator group are called different instability modes in this paper.

It is well-known that the transient power angle instability mode of a power system is closely related to the controlling unstable EP (CUEP). The post-fault trajectory of the system evolves along the unstable manifold of the CUEP. Considering that the power system’s loss of transient power angle stability is always accompanied by a sharp and excessive increase of the branch angle in certain or some critical cut-sets of the power system, it has been realized that the critical cut-set (see Note 3) in the power system is equivalent to the unstable EP (UEP) of the system, and the cut-set stability criterion has been put forward.To enable researchers to correctly apply the stability criterion related to the critical cut-set, the following fact has been proved theoretically in Ref. [25]:

Fig. 8. (a) Voltage stability limit and (b) the relationship between instability modes and the proportion of IMs in the load.

Fig. 9. Correlation between the DSR and the transient stability region [18].

Fact 7.If the CUEP of transient power angle stability is k-type,there will be k corresponding critical cut-sets in the network.

Note 3:If there is a cut-set in a network and the absolute values of all branch angles (radian) in this cut-set are greater than π/2,this cut-set is called a critical cut-set and the branches comprising it are called saturated branches. An individual critical cut-set divides the whole network into two connected sub-graphs, which correspond to two node groups, respectively. Of these groups, the one in which the voltage phase angle of the nodes is obviously larger is defined as a critical node group,and the other one is called a non-critical node group [25]. So, a critical node group contains—and only contains—generators in the unstable generator group.

Remarks:The correspondence is as follows:

(1) Most of the UEPs in a system are 1-type. If the CUEP of a system is a 1-type hyperbola (see Note 4), there must be a unique critical cut-set in the network that will separate the system into two parts: a critical node group and a non-critical node group,which correspond to the unstable generator group and the non-unstable generator group, respectively.

(2)If the CUEP in a system is an k-type hyperbola,there must be k critical cut-sets in the network corresponding to it.These critical cut-sets separate the system into k critical node groups (i.e., k unstable generator groups) and one non-critical node group (i.e.,one non-unstable generator group). The CUEP corresponds to k instability modes.

Fig. 10 gives an example of the relationship between the transient instability mode and the critical cut-set in a system of four generators and 11 nodes. In the case of a fault on branch 6–9(which refers to the branch between bus 6 and bus 9 in Fig.10(a)),instability among three groups will occur in the system.The CUEP of this instability mode is 2-type and,in the corresponding system, there are two branches (i.e., branch 10–8 and branch 7–5) whose angle differences are sharply changed in separate direction to form two critical cut-sets.

Through Fact 7, the pure mathematical concept of the CUEP is linked with the physical properties of a network, thus providing useful information for the analysis of complex instability modes of power systems from the perspective of network topology.Based on this fact,the cut-set stability criterion can be easily extended to the multiple generator group instability. A modified cut-set stability criterion has been established to eliminate the conservatism of the original cut-set stability criterion [14].

Note 4:The basic types of EPs can be identified by the local linearized expression of the dynamic system near these points.That is, the stability type of an EP can be decided by the eigenvalues of its linearized dynamic matrix. If there is no eigenvalue whose real part is zero, the number of eigenvalues with positive real parts will be regarded as the exponent of this EP. If this exponent is equal to zero,this EP is stable in linear approximation and is generally referred to as a stable EP (SEP). If this exponent is equal to i,this EP is called an k-type UEP(UEP-k).In general,an EP can be called a hyperbolic equilibrium point only when all the eigenvalues of the local linearization of the dynamic system have no zero real parts [7].

4. Topological and geometric characteristics of SRs and methods to determine SR practical boundaries

As mentioned above, the most basic issues for power system security and stability identification include power flow stability,transient stability,static voltage stability,and low-frequency oscillation stability. Therefore, this section will introduce the topological and geometric characteristics of SRs and methods to determine their practical boundaries.

Based on Hypothesis 1 and Hypothesis 2, four important facts about the boundaries of the above SRs in real power systems have been found through a great deal of simulation research and theoretical analysis.

Fact 8.The SSSR, ΩSS, which is defined in the power injection space and the decision space to ensure the power flow security of transmission and distribution networks, is the intersection of ΩT(the thermal stability SR to ensure the thermal stability constraints (THSR)), ΩV(the steady-state voltage SR to ensure the node voltage magnitude constraints), and the region constrained by the equipment capacity limits. For a given network topology and system component parameter is, the SSSR ΩSS(is) is uniquely determined, connected, and void-free, and does not change with operation states. From an engineering perspective, both ΩTand ΩVcan be described approximately by a hyper-polyhedron, and their boundaries can be described as pairs of HPs. The region between each pair of HPs respectively corresponds to the thermal stability SR of a given line or the steady-state voltage SR of a given node. Such an SSSR is called a practical SSSR (PSSR) in this paper.

Remarks:

(1) In 1989, the thermal stability SR was studied in active current injection spaces based on the decoupled power flow model [11],in which the affine relationships between node voltage angles and branch angles (i.e., the difference of voltage phase angles at both ends of the branch) are used, and the maximum allowable value of line angles is used to approximate the maximum line current magnitude. Meanwhile, all node voltage magnitudes of the system should be known in order to generate the mathematical expression of the SR. In the three-node system shown in Fig. 11(a)[11], node 0 is the reference node and its complex voltage is V0∠0. With the given node voltage magnitudes V = (V1, V2)T, the angle differences between the two endpoint nodes (θ1–θ0, θ2–θ0,θ1–θ2) is proportional to the active power of lines 1–0, 2–0, and 1–2, respectively. On this basis, an affine transformation relation can be established between (θ1, θ2) and the node active current injections (Ip1,Ip2) (which are proportional to the node active power injections). For the affine transformation shown in Fig. 11[11], the preimage is the rectangle and the image becomes the parallelogram.Here,the image is the thermal stability SR in active current injection spaces, each pair of edges in the parallelogram corresponds to the thermal stability boundaries for the line currents in the forward and reverse directions, and the region between two parallel edges is the thermal stability SR for the line current.

Fig. 10. Diagram of (a) a four-generator 11-node system and (b) the angle difference curve after a contingency.

Fig. 11. Affine transformation diagram. (a) The three-node system; (b) the thermal stability SR in the voltage angle space; (c) the thermal stability SR in the active current injection space. TB: tree branch. The expressions C1（V）∈Rn×n and D（V）∈Rn×n are given in Ref. [11].

Then,the affine transformation shows the following characteristics [11]: ①The image and the preimage are in one-to-one correspondence, including the vertex, side, edge, and internal point;②the images of parallel edges are still parallel edges; and ③the images of parallel planes are still parallel planes. When the above conclusions are applied to a system with n1nodes (excluding the reference node), if the node voltages in the decoupled power flow model Eq. (2a) are specified, it can be seen that the thermal stability SR of power systems can be approximately expressed as the convex hyper-polyhedron enclosed by nbparallel HPs in n1-dimensional Euclidean space, each pair of HPs corresponds to the thermal stability boundaries for a line current in the forward and reverse flow directions, and the region between two HPs is the thermal stability SR of the line current. These features make the form of the SSSR boundaries concise. Even when studying the SR based on the AC power flow, these features are beneficial for understanding the entire SR. In addition, based on the ideology of affine transformation,Ref.[12]has studied the steady-state voltage SR in reactive current injection spaces.

However, because the above method simplifies the power flow model and line current constraints,and does not consider the relationship between the active and reactive line currents, there is an obvious error between the generated HP boundary and the real SR boundary.

For convenience in the description,the AC power flow equation Eq.(2d)F(x)=F(xα,xβ)=0 is rewritten as F′(xα)=xβin the following. Obviously, F′(xα) is continuous nonlinear mapping, and the preimage (i.e., the definition space of xα) is a hypercube defined by Eqs. (7a)–(7c). Under this mapping, ①the image of the vertex of the hypercube is a vertex;②the image of the edge(straight line)of the hypercube is an edge,but not a straight line;③the image of the boundary of the hypercube is a continuous (smooth) curved surface, but not a HP; and ④the images of the hypercube that is internally free of void are still free of void.Therefore,it can be concluded that, strictly speaking, the SSSR defined by the AC power flow equations is a polyhedron enclosed by several smooth hypersurfaces and is internally free of void.

(2) Ref. [26] examines the calculation method for the thermal stability SR in decision spaces. Based on the AC power flow model and the sensitivity method, the proposed method successively generates the mathematical expression of the corresponding thermal stability SR boundariesfor each line i ∈B.The boundaries are pairs of HPs, and each HP corresponds to the critical current limit in one direction. The method to determine HPs is also provided in Ref. [26]. This method involves searching for a reference critical point on the boundary,which is the first encountered point on the SR boundary when continuously extending outwards from the initial normal operation point xoβ. It can be seen that this method satisfies Hypothesis 2. Then, the thermal stability SR of the entire power system is the intersection of the thermal stability SR of all lines:

It should be noted that the HPs in pairs are no longer strict parallel planes, as described in the affine transformation in Remark (1).

Eq. (16) is expressed as many simultaneous linear inequalities,and can be easily dealt with by computers.

(3)The above references mainly focus on the SSSR of transmission networks ΩSS. For a distribution network, since it uses the same power flow model as that of transmission networks, they likely have the same abovementioned characteristics. Indeed, in Refs. [27,28], it is demonstrated that Eqs. (14)–(16) can also be used to define the steady-state voltage SR and the thermal stability SR in power injection spaces of the distribution network,and that a pair of HPs can also be used to describe the corresponding SR boundaries (such kind of characteristics are consistent with the transmission network).The difference is that Refs.[27,28]propose simple generation methods for the SSSR based only on the radial topology of distribution networks, and HP coefficients can be directly generated by using the network topology and line impedances. The simulations show that the proposed methods are fast enough to meet the needs of a distribution network with frequent changing topologies. By combining these generation methods of SRs in distribution networks with the abovementioned generation method of SRs in transmission networks,coordinated optimization of transmission and distribution networks can be executed effectively. Considering that the loads, distributed generations, and energy storages can be expressed as power injections,this research can be used in the monitoring, optimization, and security assessment of smart grids. (Specific application examples are given in Section 5.) Ref. [29] examines the thermal SR of distribution systems. Under the assumption that the node voltages are specified,SRs are used to describe the relationship of the power flows in radiant distribution networks and are applied to evaluate the distribution network security.

(4)Because the ΩSS(i)defined in this paper satisfies Hypothesis 1 and Hypothesis 2,the boundary of ΩSS(is)is uniquely determined and connected for a given network topology and given system component parameters in the space defined by xβ. Furthermore,ΩSS(is) is internally void free due to the continuity of the region restrained by Eqs. (7a)–(7c) and the power flow mapping F(x) =0.

(5) Each HP of the SSSR in the form of a hyper-polyhedron can be described by the following equation:

where αj,iand βj,iare constant coefficients,and msis the total number of SR boundary surfaces. Because the SSSR is connected and uniquely determined for a given network topology and given system component parameters,and is irrelevant with the system operation states, all HP coefficients can be calculated offline and stored for online security monitoring, assessment, and optimization.

Fact 9.For a given network topology, system component parameters, pre-fault graph id, post-fault graph jd, and fault Fd,the DSR Ωd(id, jd, Fd) defined in power injection spaces is uniquely determined and connected, and does not change with the operation states. Also, it is internally void free, and its boundaries Ωd(id,jd,Fd) may consist of a finite number of smooth subsurfaces, each of which corresponds to a specific instability mode.Each sub-surface can be approximately described by an HP within the scope of practical engineering. The DSR described in such way is called the practical DSR (PDSR) in this paper.

Fig.12 shows the diagram of the New England 10-machine 39-node system and the cross-sections of its PDSR in three 2D spaces.The PDSR corresponds to the three-phase grounding short-circuit fault on bus 26 of lines 26–29. The duration of the fault is τ ＝0.1 s, and the fault is cleared by disconnecting lines 26–29.

Remark:

(1) In 1990, by fitting a large number of transient stability critical points obtained by simulation,the authors of Ref. [13]first discovered that the PDSR boundary in the power injection space of the pre-fault system, which is shown in Eq. (8e), can be approximately described by a critical HP around the ‘‘basic operation point”in the scope of practical engineering application.The critical HP can be expressed as follows:

The absolute value of the HP coefficient represents the influence of the corresponding node power injection on the system stability.The positive or negative value of the HP coefficient represents the influence trend of the corresponding node power injection on the system stability; that is, an increase of the node power injection with a positive coefficient will deteriorate the system stability,while an increase of the node power injection with a negative coefficient will improve the system stability.

The discovery of these characteristics is of great significance because the linear combination constraint in the power injection space, which is shown in Eq. (18), is of great advantage in mathematical processing for power system assessment, operation, and control.

In Ref. [13], the potential energy boundary surface (PEBS)method [7] is used to find the critical injection vector, and then the HP coefficient αiin Eq. (18) is determined by the leastsquares fitting method. Ref. [13] recommends the use of the quasi-orthogonal method to select possible searching directions in order to achieve uniform distribution of critical points, which can reduce the number of required critical points and ensure accuracy.

(2)In Ref.[30],the 984-node system of the Central China Power Grid is used to validate the proposed method,considering the dualaxis reaction of the generators,the excitation system and the speed governing system, features of the IMs, and the static var compensator(SVC).Using transient stability simulation tools,a great number of critical points on the DSR boundaries are searched through numerical simulations. The authors further prove that the DSR boundary for transient power angle stability can be approximately described by an HP via least-squares fitting. Each critical point is the SR boundary point that is first encountered by continuous outward expansion from the initial normal operation point x0βalong different ray directions in the active power injection space, as shown in Eq. (8e).

(3) The least-squares fitting method requires approximately 2n suitably distributed critical points.It can only be calculated offline and then used online due to the great amount of calculation required, although highly accurate SR boundary coefficients can be obtained. To solve this problem, the analytical expressions of the HPs of the PDSR boundary are derived in Refs. [31,32] based on the structure-preserving model of a power system, as well as the sensitivity of the initial operation state, the operation state at fault-clearing time, and the energy function to the node power injections.Based on the fact that the Gramian matrix of the power system dynamic model at the short-circuit fault-clearing time remains almost unchanged compared with the Gramian matrix of the power system dynamic model near the initial operation point,Ref.[33]proposes a method,that combines coherency identification based on k-medoids method with the initial acceleration,to rapidly identify the transient instability mode(and thus to identify the critical cut-set) near the initial operation point.

Fig.12. The diagram of the New England 10-machine 39-node system and the cross-sections of its PDSR in three 2D spaces.(a)The diagram of the New England 10-machine 39-node system;(b)PDSR in the active power injection space of G8 and G9;(c)PDSR in the active power injection space of G9 and L29;(d)PDSR in the active power injection space of G8 and L28. The line from M to N is the transversal of the critical HP of the PDSR in the 2D space. The active power limits of generators are shown as dotted-lines.

(4) In order to provide a theoretical basis for the application of the DSR, Ref. [34] has proved the following characteristics of the DSR based on the differential topology theory for the structurepreserving model of power system transient stability analysis

[20] of a power system: ①The DSR determined by the CUEP method is internally void free, that is, showing density; ②?Ωd(id,jd,Fd),the boundary of the DSR,will not knot—that is,it is not torsional or expansive—and the local surface corresponding to the same instability mode kdis continuous; and ③the boundary of the DSR is compact—that is, it can be expressed by the union set of a finite number of sub-surfaces, and each critical sub-surface corresponds to a different instability mode.

(5) According to Fact 7, when the CUEP of the system is an i-type hyperbola, there must be i critical cut-sets in the corresponding network, that is, there must be i instability modes in the system corresponding to the CUEP. Therefore, the corresponding ?Ωd(id, jd, Fd) must have i sub-surfaces, each of which corresponds to an instability mode.

(6) Furthermore, the cut-set power space can be used to simplify the description of the DSR. However, the SR in the cut-set power space is not unique, as it is not defined in the xβ-defined space (see Fact 2). As indicated by studies on the DSR in a cut-set power space (CDSR) in Refs. [35,36], for a given network topology and given system component parameters, the total transfer capacity (TTC) in the critical cut-set of the system is relevant to the operation states (power injection vector), that is, the CDSR is not uniquely determined. For this reason, it is necessary to seek the maximum value of the TTC, by which the HP of the PDSR in the power injection space can greatly reduce the computational burden. However, the boundary of the CDSR may vary in a small range and can be used approximately in a few practical scenarios.

(7)A large number of simulation studies in Ref.[37]show that,for a large number of different (id, jd, Fd), all ?Ωd(id, jd, Fd) can be approximately expressed by the following HP in a complex power injection space defined by Eq. (8c) [7]:

If the changes in both the node active power injections and the node reactive power injections are considered, especially when involving the transient voltage stability,Eq.(19)will have a better accuracy. According to Ref. [37], the boundary of the PDSR in the pre-fault power injection space has the approximate translation property in the case of a change in the proportion of the IMs in the load.Therefore,the interpolation method can be used to determine the DSR boundary under different proportions of IMs in the load.

For the load that consists of both constant impedance loads and IMs, transient stability simulations show that the critical faultclearing time of the transient stability is closely relevant to the proportion of motors in the load,and the increase of the proportion of motors in the load will result in the decline of the critical faultclearing time. In addition, the system transient instability mode mainly manifests as transient voltage instability when the proportion of motors is high and as transient power angle instability when the proportion is low, as shown in Fig. 13.

(8) Ref. [38] shows that, for an alternation current (AC)–DC hybrid transmission system,the critical HP description of the PDSR in the decision space still has good accuracy. The boundary equation of the PDSR corresponding to a certain instability mode can be expressed as follows:

where G and L are the set of generator nodes (except equilibrium nodes) and the set of load nodes, respectively; Viis the voltage of generator node i; and αi, βi, ηj, λj, and μ are the critical HP coefficients of the PDSR.μ ＝αrb-αsb, where sb and rb are the number of AC nodes at the sending end and the receiving end of the DC line,respectively. Therefore, the DC power Pdcan be regarded as an active power variable of the critical HP equation. Under the same failure mode,instability mode,and DC control mode,all the critical HPs corresponding to different DC power have approximate parallelism, and the spatial geometric distance between them is approximately proportional to the change in DC power.

Fig. 13. Relationship between the transient instability mode and the proportion of IMs in the load.

(9) Ref.[39]examines the DSR of a power system with doublyfeed induction generators (DFIGs). The accuracy of its boundary is verified through time-domain simulation.The DSR of a power system with photovoltaic (PV) generation is calculated in Ref. [40].According to Refs. [39,40] and within the engineering perspective,the boundaries of the PDSR in the power injection space can still be approximately described by HPs for a power system with DFIG or PV generation. Ref. [39] also analyzes the influence of the integration of DFIGs on the DSR, and finds that the integration of DFIGs would cause external expansion of the DSR in a power system.

Fact 10.For a given network topology and given system component parameters, the SVSR in the power injection space corresponding to fold bifurcation is uniquely determined, connected,internally void free,and does not change with the operation states.Furthermore, its boundaries are smooth. When the load node power injections are specified, the boundary of the SVSR in the generation power injection space can be approximately expressed by an HP in the practical engineering range.When there is only one critical cut-set in the system,CVSR can be approximately described by an HP within the practical engineering range (as shown in Fig. 14). When several critical cut-sets should be considered in the system, the boundary of the CVSR corresponding to each critical cut-set can be approximately described by an HP within the practical engineering range. In addition, the union set of the boundaries of the CVSR of all cut-sets is the boundary of the CVSR under a given network topology and given system component parameters.

Remarks:

(1) In Ref. [41], since the voltage stability problem has strong local characteristics, the voltage stability boundary can be approximately described by the injected power of the nodes with weak voltage stability in the system.Its relatively complex boundaries can be visualized by an artificial neural network.In Ref.[42],a quadratic polynomial is used as an approximate analytical expression for the large-scale boundaries of the SVSR in active and reactive power injection spaces. Meanwhile, in order to reduce the dimension of the power injection space, the key nodes of the system are selected by using modal analysis.A large number of simulation results show that this method has satisfactory engineering accuracy. In Ref. [43], an approximate analytical expression of the SVSR boundary is constructed by the eigenvalue sensitivity and eigenvector sensitivity of Jacobian matrix.Although the precision of this method is higher than that of the HP-based linear approximate expressions, it is far less convenient than the HP method when applied in optimization and risk analysis.

Fig.14. Local diagram of a real power system and boundary of the CVSR on critical cut-set marked by dotted lines.(a)The diagram of a real power grid;(b)the boundary of CVSR in a 3D power injection space.

(2)Ref.[44]is based on the power flow feasible solution region.It takes whether the power flow equation has a solution or not as the decisive factor of the SSSR by using a new‘‘hybrid method”for feasible region boundary computation,proposed by the authors in Ref. [45]. With a prediction–correction ideology and optimization technique, this hybrid method takes into account the limit of the system equipment while tracking the boundary of the power flow feasible region composed of singular points of the Jacobian matrix of the power flow equation (fold bifurcation) and LIB. In Ref. [44],the power injection space is divided into the load injection space and the power generation injection space.By fitting a large number of critical points obtained by static voltage stability simulations, it can be found that the boundaries of the power flow feasible region in the above two sub-spaces have distinct geometric characteristics.In the power generation injection space with a specified direction of load growth, a part or the majority of the boundary of the power flow feasible region can be well approximated by an HP in the practical operation range, as shown in Fig. 15 [24]. Another new algorithm is put forward in Ref. [44] to track the farthest (L1norm) boundary point of the power flow feasible region in the high-dimensional power generation injection space, which avoids the frequent start of CPF calculation. The multi-solution of the results is also used to further verify that the boundary of the power flow feasible region can be approximated by an HP.

Note 5:Because the vector of power generation has a very high dimension, an SR as whole is difficult to visually represent.Therefore,the above figures show only the transects of the SR of concern in a 2D power generation space.All transects are obtained by supposing that the variables other than the given coordinate variables remain unchanged. If the edge lines of the transects all appear as straight lines in a large range of actual engineering concerns, then it can be inferred that, in a high-dimensional space, the boundary surface of SRs can be approximated by an HP in a large range of engineering concerns.

(3) Ref. [46] provides a local visualization method for the steady-state voltage stability region in the cut-set power space based on the ideology of the power flow feasible region. Refs.[47,48] clearly propose that the practical boundary of the CVSR can be approximately described by Eq.(21)with an HP expression:

In Fig. 16, a method to quickly determine the HP of the CVSR boundary is presented. First, the CPF is used to find a SNB point of the system. Then, the tangent plane of the SVSR boundary at the SNB is represented by the eigenvector at the critical point.Finally, the critical boundary of the CVSR with an HP expression is obtained by the transformation from the power injection space to the cut-set power space.

(4)When there is a large error,some power injections of nodes with weak voltage stability at the receiving end of the system can be used as additional variables of Eq. (21) to improve the fitting accuracy.

(5)Considering that the number of lines in the critical cut-set is very limited, the HP coefficients αiand βican be obtained by the least-squares fitting method after finding about 4nc(where ncis the total number of lines in the critical cut-sets)appropriately distributed critical points through simulation.Based on this method,a software to determine the CVSR has been developed and the results show that the software can meet real-time online needs in practice [48].

(6)In some systems or cases,the weak nodes of the system may be distributed in different areas, which requires multiple critical cut-sets [46]. Therefore, all the operation points that ensure static voltage stability are in the region defined below:where C1,...,CmSVare the msvcritical cut-sets that need to be considered, and msvis generally a small integer.

Fig. 15. View of 2D transects (see Note 5) of the power flow feasible region of IEEE 118-bus system in a power generation injection space [24].

Fig. 16. A method to determine the HP of the CVSR boundary. Pg: a critical point on the boundary of the SVSR.

(7)The statistical analysis in Refs.[48,49]shows that,under the actual possible load growth mode, the weak node sets in the system are concentrated in several areas and most of them are distributed at the terminals of long transmission lines without reactive power support.The graph model corresponding to the grid structure can be established according to the connection relationship and the degree of closeness between nodes.The spectral clustering algorithm can be used to determine the voltage stability area(the node set with weak voltage stability in the power system)and the critical cut-set (the potential critical cut-set)of the power system.

(8) Static voltage instability is based on the absence of power flow solutions, and the SR defined in this paper complies with Hypothesis 2; that is, the SVSR boundary is the first static voltage instability boundary encountered by continuous outward extension from the initial operation point.Then,along with the continuity of power flow functions F(x) =0, it can be concluded that the SVSR is uniquely determined, connected, and internally void free.Its surfaces are smooth in xβ-defined space. The mapping of this region to the low-dimensional cut-set power space, that is, the CVSR, can retain these properties. Therefore, the CVSR in Fact 10 is also internally void free.

Fact 11.For a given network topology and given system component parameters, the SDSR, ΩSDin the power injection space is bounded by SNB boundary (?ΩSD(SNB)) and/or HB boundary(?ΩSD(H)). It is uniquely determined, connected, and internally void free, and does not change with operation states. ?ΩSD, the boundary of ΩSD, is composed of several smooth surfaces, and the sudden change occurs at the junction between ?ΩSD(SNB)and ?ΩSD(H), or on ?ΩSD(H) when its dominant oscillation mode changes. Each smooth surface can be approximated by an HP in a large range(as shown in Fig.17[50]),which can meet the practical engineering requirements.In this paper,the SR described in such a way is called a practical SR to ensure small-disturbance stability(PSDSR).

Remarks:

(1) In the investigation of small-disturbance stability [21], the detailed evolution process from HB to chaos is demonstrated by an illustrative power system, and the law of increasing energy is revealed. It is also found that the boundary of the SDSR relevant to oscillatory instability can be described by the HB boundary?ΩSD(H), and that there is no need to account for more complex chaos.

(2) Ref. [50] shows that ?ΩSD(H) in the active power injection space is composed of several smooth curved surfaces, each of which can be approximately described by Eq.(23),where the reactive power is assumed to be locally balanced:

Fig. 17. Piecewise fitting under multiple dominant modes (New England 10-machine 39-node system) [50].

where Piis the injected active power at node i,αiis the HP coefficient of node i,and n is the total number of other nodes in the network without the reference node.

(3) ?ΩSD(H) may be composed of one or several smooth surfaces,and the reason why the sudden change occurs at the junction among surfaces is that those surfaces have different dominant oscillation modes.Based on the above knowledge,Ref.[50]applies the following piecewise fitting strategy to obtain the approximate HP for ?ΩSD(H). Eigenvalue analysis is first conducted for each critical point to determine its dominant oscillation mode. Then all the critical points are divided according to their dominant oscillation modes, and those with the same dominant mode are clustered into the same set. The critical points in each set can be fitted by an HP, and the union set of the fitting boundaries corresponding to all dominant oscillation modes constitutes ?ΩSD(H)in the power injection space (as shown in Fig. 17 [50]).

(4) The influence of DFIGs on power system electromechanical oscillation is examined in Ref. [51], and it is found that the SDSR boundary of a power system with DFIGs is composed of several smooth surfaces,including ?ΩSD(SNB)and ?ΩSD(H)under different dominant oscillation modes. Fig. 18 [51] shows the 2D crosssection of a SDSR boundary in the active power injection space of two generators. When the active power output constraints are not considered, the boundary is composed of several smooth curves.Boundaries 1 and 3 correspond to HB,and boundary 2 corresponds to SNB. There are sudden changes at the junction between boundaries 1 and 2 and at the junction between boundaries 2 and 3. When the active power limits (shown as the rectangle in Fig. 18 [51]) of the generators are considered, the actual SDSR of concern is the shaded area. The appearance of SNB on boundary 2 is due to the fact that the excitation system of generator G4 reaches its limit, which is the LIB problem in Fact 5.Therefore, it is necessary to determine the boundary of an SDSR through time-domain simulation, as is described in Fact 10.

Fig. 18. 2D section of an SDSR boundary [51].

(5)The uniqueness,connectivity,and the void-free characteristics of SDSRs can be proved in a similar way to those of Remark(8)in Fact 9, which will not be repeated here.

5. Methodology of SR with its boundaries expressed by HPs for power system analysis

5.1. Application of SR methodology in power system optimization

The power system optimization problem with security and stability constraints can be described by the following general model.Let xoptbe the optimization variable vector,which usually refers to the active/reactive output of generators or the active/reactive power of loads. fopt(xopt) is the optimization objective. Securityconstrained optimization problems of power systems can be described by the following general model, as Eqs. (25)–(30) show.The model will change with the specific scenarios of applications.For example, in security-constrained optimal power flow, the objective function usually considers the minimum power generation cost or loss, while in the security-constrained optimal control of power systems, the objective function usually considers the minimum control costs, such as the load shedding cost and so forth. hopt(xopt) is the equality constraint function, which mainly includes the power flow equations.gopt(xopt)is the inequality constraint function, including the generator output limit, the node load constraint, and so forth. Eqs. (27), (28), (29), and (30) represent the node voltage constraint, power flow security constraint,small-disturbance stability constraint, and transient stability constraint, respectively.

Therefore, based on the fact that the boundaries of SRs in the power injection space can be approximated by HPs, system security constraints such as the power flow constraint, smalldisturbance stability constraint, and transient stability constraint can be formulated as explicit linear inequality functions of decision variables (i.e., node power injection), so that it will be extremely convenient to take into account various types of security constraints in power system optimization, and the calculation efficiency will be significantly improved.

Up to now,SR methodology has been widely applied to various power system optimization problems.For example,an optimization model of power system security cost in the electricity market is proposed in Ref. [52]. The transient stability constraint, fault probability,and system instability loss are considered by means of PDSR, through which optimal preventive control can significantly improve the economic efficiency of the power system. Ref. [53]reports on experiential laws about the parallelism and superposition of the critical surface of the extended PDSR. Based on this, it provides an optimal transient stability-constrained emergencycontrol method,and quantifies the control measures as cost indicators, which solves the problem of the difficulty quantifying the effectiveness of emergency-control measures. In Ref. [54], optimal active and reactive coordination of a power transmission system is realized considering security constraints based on SRs.By means of the SSSR,SVSR,and DSR,the power flow constraint,static voltage stability constraint,and transient stability constraint are considered comprehensively. Through the affine relationship between the active power and the branch angle of the power system and the affine relationship between the reactive power and the node voltage magnitude, an algorithm for solving the optimal power flow through QP is established,which significantly improves the calculation efficiency.In Ref.[55],the power flow constraint,steady-state voltage stability constraint,and transient stability constraint are,for the first time,simultaneously taken into account in the day-ahead dispatching of power systems, providing an important approach to solve the contradiction between security and economy in power system operation. In Ref.[56], models for the pricing of active and reactive power as well as corresponding solution algorithms are proposed, which systematically consider many complicated constraints,including quantification of transient stability in node pricing for the first time. In the proposed model, the static voltage stability constraint and the transient stability constraint are considered through the CVSR and PDSR, respectively, and a decoupled optimization and iteration method for active power production cost and reactive power production cost is suggested.Based on the marginal cost theory and Karush–Kuhn–Tucker(KKT)optimality condition,active power and reactive power are priced,respectively,and the component prices related to various security constraints are derived. In the proposed method, not only is it convenient to consider the contingency set, but this is also done with more concise expressions and clear physical meaning.In Ref.[57],SR methodology is adopted to establish static security and dynamic security value-based power system extension planning. In Refs. [27,28],the SR ideology is introduced into the reactive optimization of a power distribution network for the first time, providing a critical decision-making tool for the reactive voltage control of a smart power distribution network. An example of SR application with HP boundaries in the coordinated optimal power flow of a power system is provided here,as follows.

Example:SR application with HP boundaries in the coordinated optimal power flow of a power system with a given topology.

With increasing penetration of intermittent, fluctuating, and uncertain renewable distributed energy in smart distribution grids,transmission and distribution grids can affect the tie-line power by changing their own operation states. This means that it is necessary to consider the impact of tie-line power changes on each other’s operation state in both the transmission grid and the distribution grid.Hence,it is necessary to comprehensively consider the power generation resources in the entire grid for coordinated optimal power flow.

Nowadays, the methods for solving this coordinated problem can be mainly divided into two types. One type is the centralized optimization method [58]. Although this method can yield optimization results, it is necessary to collect and process the entire grid data during the application, and the difficulty of data maintenance and computation burden are very large.The other type is the distributed optimization method, which decomposes the coordinated optimal power flow problem of the entire transmission and distribution grids into several sub-problems for optimization.The solution speed of the latter method is faster than that of the centralized method [59]. However, the constraints in each subproblem include capacity constraint inequalities for devices, node voltage constraint inequalities,line current constraint inequalities,and AC power flow equations (similar to Eqs. (25)–(28)), and the decision variables in the objective function are power generation injections, which are nonlinearly related to node voltages and line currents based on power flow equations; thus, the determination process for each sub-problem and for the KKT conditions is still complicated.

In order to solve the problems in the existing methods,Ref.[60]applied the SSSR with HP boundaries to the optimal power flow problem of the entire transmission and distribution grids for the first time, and established a coordinated optimization method for the transmission and distribution grids. The characteristics of this method are as follows:

(1)The approximate HP expressions for the SSSR boundaries are generated based on the respective topological structure and boundary information of the transmission grid and the distribution grids in their power supply area. Considering that the SSSR is uniquely determined by the given grid topology and is irrelevant to operation states,for a transmission grid whose topology is relatively stable, the HP expression coefficients can be calculated offline and recalled online. For a radial distribution network, the HP expression coefficients can be generated online, which requires a small computation burden and can be fully adapted to the requirements of frequent changes in topology.

(2) In order to minimize the overall power supply cost of the transmission and distribution grids, SSSRs in the transmission and distribution grids are applied to describe the steady-state security constraints(as shown in Eq.(32)),and to establish the distributed optimization model for coordinated optimal power flow in the transmission and distribution grids;the boundary variables are shown in Fig.19.During the execution of the proposed method,the transmission grid and distribution grids can complete the optimization only according to the boundary node voltage, the tieline power and the clearing prices(determined by KKT conditions)exchanged between the grids. Since the constraints become linear combination inequalities of decision variables, and the objective function is consistent with the variables in the constraints, the optimization process of each transmission grid and distribution grid and the determination of the KKT conditions become extremely simple.

In order to verify the effectiveness of the proposed method,Ref.

[60]first compares the proposed method with the centralized optimization method[58].In the optimization process,the centralized method establishes a nonlinear optimization model based on an AC power flow model. The overall supply cost/calculation time obtained by the two methods is shown in Table 1.The optimization results obtained by the two methods are basically the same, but there are great differences in the calculation time.

Note 6:

(1) In order to briefly show the advantages of the practical SR boundary with HP expression in power system optimization, only the SSSR is considered in this example.If other security constraints are taken into account,the advantages of the SR methodology will be even more prominent.

Fig. 19. Equivalent model of transmission and distribution grids. (a) Equivalent model of the entire grid; (b) equivalent model of transmission and distribution grids;(c) equivalent model of distribution grid k.

Table 1 Comparison of daily overall supply cost and calculation time (see Note 6).

(2)Case 1 is a modified IEEE 24-bus grid with nine distribution grids. Case 2 is a modified IEEE 118 with 20 distribution grids.Case 3 is a modified IEEE 300 with 40 distribution grids.

(3) Since the SSSR boundaries of the transmission grid are uniquely determined by the grid topology, they are independent of operation states, and the topology does not change frequently.Therefore,the coefficients of the SSSR boundaries can be calculated offline, stored in the database, and recalled according to the topology when used online. However, the coefficients of the SSSR in distribution grids need to be calculated in real time. Therefore,for the statistics of optimization time burden, the generation time of the SSSR in the distribution grid is included, and the generation time of the SSSR in the transmission grid is not.

The proposed method is then compared with the alternating direction method of multipliers method (ADMMM) [61] and the optimality condition decomposition (OCD) method [62], which are commonly used. The average optimization time consumed by the three methods in each iteration is shown in Table 2. It can be seen that since the proposed method is based on the SR method to establish the optimization models in the transmission and distribution grids,the variables in the objective function and the constraints are identical, and the constraints are linear combination inequalities of node power injections, so that the optimization speed of each sub-problem is fast.

Furthermore,the average number of iterations required to complete a coordinated optimal power flow by the three methods is shown in Table 3. Since the Lagrangian multiplier corresponding to each constraint is needed in the process of determining the KKT condition, and considering that the constraints used in the optimization model of the ADMMM method and the OCD method are nonlinear and the constraints (i.e., the boundary HPs of SSSR)used in the proposed method are linear, the number of iterations varies greatly.

In summary, the total time required to complete the coordinated optimization of the entire transmission and distribution grid is mainly determined by the average optimization time consumed by each iteration of the transmission grid multiplied by the average iteration times to complete a coordinated optimal power flow.The proposed method has obvious advantages in both of the above factors, so that the average total time burden for coordinated optimization can be decreased by several orders of magnitude (as shown in Table 4).

This example shows that the SR application with HP boundaries presented in this paper can successfully solve the problem of the security and stability constraints being difficult to deal with in a large class of power system optimization problems.

It should be noted that the average time consumed by the proposed method to generate the SSSR in the transmission grid is 1.67,6.3,and 45.1 s for Cases 1,2,and 3,respectively,showing that even if the average time for generating the SSSR of the transmission grid is included in the optimization time statistics, the superiority will not be affected when compared with the existing methods.

Table 2 Comparison of average optimization time in each iteration.

Table 3 Average iteration times to complete a coordinated optimal power flow.

Table 4 Average total time for coordinated optimization.

5.2. Application of SR methodology in power system probabilistic security assessment

Power systems often suffer from various disturbances, such as fluctuations of the node injected power,incidents,and component outages.The system is considered to be secure if it is capable of bearing any disturbance possibly occurring at the next moment.Considering the inherent uncertainty in the fluctuations of node injected power, incidents, and so forth, especially after the integration of large amounts of wind and PV power,probabilistic security assessment becomes even more important. Born at the right moment,probabilistic security assessment aims to obtain the security probability for a system to satisfy the security constraints or the probability distribution of key state variables considering the uncertainties of node injection power,incidents,and so forth,in order to indicate the system security level in the sense of probability for the near term(e.g.,for every half hour or even shorter time interval in the coming 24 h).

Probabilistic security assessment is usually conducted with the simulation method or the analytical method. The simulation method is based on the Monte–Carlo (MC) simulation approach.First, a large number of operation points are randomly generated.Next, the stability of each operation point is judged one by one,that is, one simulation computation will be carried out separately for each operation point to determine whether it is stable. Finally,a probability will be obtained based on the law of great numbers.This method can easily consider various uncertainty factors, but its calculation error is inversely proportional to the square root of the test times.Thus,more calculation time is required to reduce the error, which leads to a considerable calculation burden in probabilistic security assessment. In the analytical method, a model is set up to evaluate the probabilistic security indicators,or an analytic expression is established for probabilistic security indicators. For example, in Ref. [3], the probability distribution of the time to insecure is obtained via the two-layer Markov model,while in Ref. [63], an analytical expression for system dynamic security probability is established based on the condition probability theory. In comparison with the simulation method,the analytical method has a relatively complete theory. However,considering the uncertainty of the node injection power—especially the simultaneous consideration of the transient stability constraint for the given fault—the solution process will be extremely complex. In case of the evaluation of the system probability satisfying the transient stability constraint in a certain fault considering the uncertainty in the node injection power, i.e. the probability for the vector of injection power y to be located within Ωdcan be expressed as follows:

The calculation process is thus significantly simplified. The key to solving Eq.(34)is to calculate the probability density function of ye, that is, g(ye). Whether the node injection power is relevant or not, valid methods such as point estimation and the combination of semi invariants with series can be applied to realize a fast solution with the desired calculation accuracy.

The research results herein have provided the following advantages compared with the simulation method:

(1) Application of the PDSR improves the calculation efficiency of Eq.(33)by 1×104times(the New England system is calculated on an individual computer [64]).

(2) As Eq. (33) is the most basic computational unit in the two-layer Markov model, considering the curse of dimensionality due to varied states, application of the PDSR and the typical characteristics of the HP coefficients(offline calculation and online application) improve the calculation efficiency by 1 × 106times(the New England system is calculated on an individual computer,enumerated to N-3) [65].

With the rapid development of high-performance computing,there will be unlimited possibilities for the SR method to significantly improve the calculation efficiency of probabilistic security assessment. Meanwhile, with the integration of enormous quantities of distributed energy sources, the probabilistic analysis method is playing an increasingly important role. Since uncertainty of the node injection power needs to be considered in various security-constrained optimizations of power systems, the SR with its boundaries expressed by HPs undoubtedly has a broad application space.

5.3. SR application on the visualization of power system security monitoring

Based on the research results of the PDSR in the power injection space and the CVSR in the cut-set power space,an SR visualization system for the power system has been successfully developed[48].Since the visualization can only be implemented in 3D or 2D space,the system uses the SR to display transient stability and voltage stability boundaries in the 3D or 2D power injection subspace to improve the observability of the power system, and to help dispatchers or operation planners to identify important contingency power transfer levels defined by voltage and power angle stability constraints and system security margins. Figs. 20(a) and (b) provides visualization examples of the cross-sections of the CVSR and PDSR in the 3D power injection subspace, respectively. The transparent section in Fig. 20(a) corresponds to the CVSR HP boundary, and the transparent section in Fig. 20(b) corresponds to the PDSR HP boundary (where each transparent section corresponds to a contingency,and the intersection of the two parts they surround is the SR under two contingencies). These visual images can strengthen an operator’s insight into critical situations(e.g.,the situation of voltage instability,transient instability,or the danger of a cascading outage), which enhances the ability of situation awareness assisting operators to respond immediately to prevent the consequence of a blackout.

Fig. 20. SR visualization in the 3D power injection subspace. (a) The CVSR of a specific network; (b) intersection of the PDSR under two contingencies.

6. Conclusions

This research initiates and significantly promotes SR methodology (distinct from the traditional point-wise method). The main original research achievements are as follows:

(1)In terms of the composition and dynamics nature of SRs,we have found that:

· Chaos occurs outside the boundary of the SDSR that corresponds to HB. Several types of instability and collapse induced by chaos have also been discovered. Chaos can be ignored in the study of the IGSR since the secure operation of power systems does not allow HB to occur.

· There are some rules for the impact of models and parameters on the SR boundary. For example, the presence of IMs in the load might lead to no occurrence of HB in power systems before the increase of power injections causes voltage instability of the SNB type; not all SNB points are on the boundary of the SDSR ?ΩSD∩｛SNBs｝, and thus it is essential to determine the property of the SNB points by means of time-domain simulation in order to finally determine the boundaries of the SDSR when IMs account for a large proportion of the load.

· There is a quantitative relationship between the type of UEP of transient stability and the number of critical cut-sets of transient power angle stability. If the UEP of the system is k-type hyperbolic, there must be k corresponding critical cut-sets in the network, and these critical cut-sets divide the system into k critical node groups (i.e., k unstable generator groups) and one non-critical node group (i.e., one nonunstable generator group), corresponding to k unstable modes.

(2)In regard to the topological and geometric characteristics of SRs,the critical finding is that for a given network topology(and its evolution process) and given system component parameters,within the practical engineering scope of concern, the SSSR that ensures the power flow security in the power injection space and the decision space, the DSR that ensures the transient stability in the power injection space,the SVSR that ensures the static voltage stability (composed of SNB points) in the generation power injection space, and the SDSR that ensures the small-disturbance stability (composed of HB points) in the power injection space have the following characteristics:

· The SR defined in the power injection space and the decision space is connected, uniquely determined, and irrelevant to operation states when only the range, surrounded by the SR boundary first encountered by a slow outward continuous extension from the initial normal operation point in the form of a quasi-steady state, is considered.

· There is no void inside the SR, which means that the SR is only confined by the boundaries.

· The boundaries are piecewise smooth.The occurrence of sudden changes at the junction between the smooth surfaces is caused by the fact that each smooth curved surface has a different instability mode(critical cut-set)for the DSR boundary and a different bifurcation type or different dominant oscillation mode for the SDSR boundary.

· The SR can be approximated by one or a few HPs. The mathematical description of the corresponding security constraints is a set of linearly combined inequalities of the power injection variables.Since the SR is irrelevant to operation states, the HP coefficients of SRs can be calculated and stored offline for online use.

(3)Based on the facts that the boundaries of the practical SR in the power injection space can be approximated by HPs and that this SR is irrelevant to the operation states, great advantages of the SR method for the security monitoring, probabilistic security(risk) assessment, and optimization of power systems have been discovered, as follows:

· All security constraints can be formulated as the linear combination inequalities of decision variables (i.e., node power injection) in the optimization problem, such that it is very easy to simultaneously consider the security and stability constraints in power system optimization. This greatly simplifies the optimization algorithm and improves the online computation speed by orders of magnitude. Consequently,a successful solution to the difficulties of addressing security and stability constraints in power system optimization problems has been achieved, and is applicable in areas such as coordinated optimization in transmission and distribution networks,security cost optimization(optimal security control), emergency control, and security-constrained unit commitment.

· Based on the SR boundary with HP expression, the n-fold integration problem of the probability density function with n-dimensional variables in probabilistic security assessment has been mathematically transformed into a threshold comparison problem of a 1D probability distribution function.In this way,the computational burden of the online probabilistic security assessment of power systems can be reduced by several orders of magnitude. In the future, SR methodology will be a powerful tool for the security analysis of smart grids with the integration of massive uncertain distributed generations.

· SR methodology is also a powerful tool for situation awareness, since it can make visualization easier to implement and can enable rapid determination of the security margin.

(4) Analytical methods for the fast calculation of the HP boundaries of the practical SSSR(for both transmission and distribution systems), the DSR, the SVSR, and the CVSR have been invented.These methods greatly improve the computational speed and can achieve the required accuracy for practical engineering application,in comparison with the fitting method, which requires a large number of simulations to search out a large number of critical points.

· For transmission and distribution networks, the HP coefficients of SSSR boundaries have been given by the sensitivity matrix of AC power flow; moreover, a fast searching method for a reference critical point on the boundary has been proposed.In particular, the HP coefficients of SSSR boundaries for the distribution network with a tree topology can be directly generated based on the network topology and line impedance.The proposed algorithm is very simple and fast and can meet the requirements of the online real-time analysis of smart distribution grids, the topology of that changes frequently.Moreover, as the HP expressions are derived from the AC power flow model, the HP coefficients, that reflect the influence of active and reactive power injections of the nodes on the node voltage and line current,are more suitable for distribution networks with active and reactive power coupling characteristics.

· The HP analytical expression of the PDSR boundary has been derived, and a method for quickly searching for a reference critical point in the PDSR direct method has been provided.

· A hybrid method for calculating the SVSR boundary has been established. On this basis, a new algorithm has been put forward to track the farthest (L1norm) boundary point of the power flow feasible region in a high-dimensional power generation injection space. This algorithm avoids the frequent start of CPF calculation, and the multi-solution of its results is also used to further verify the characteristics of the power flow feasible region.

· Based on the fact that the Gramian matrix of the power system dynamic model at the short-circuit fault-clearing time remains almost unchanged compared with the Gramian matrix of the power system dynamic model near the initial operation point,a method,that combines coherency identification based on k-medoids method with the initial acceleration, is proposed to rapidly identify the transient instability mode (and thus to identify the critical cut-set)near the initial operation point.

· A graph model corresponding to the grid structure has been established. The spectral clustering algorithm has been used to determine network partitioning in coherent areas of static voltage stability (the non-weak node set and the weak node set; i.e., the critical cut-set), and thus to determine the CVSR of great concern.

Compliance with ethics guidelines

Yixin Yu,Yanli Liu,Chao Qin,and Tiankai Yang declare that they have no conflict of interest or financial conflicts to disclose.