STATE ESTIMATION METHOD FOR MULTI-STAGE VOLTAGE SAG

Abstract

A state estimation method for a multi-stage voltage sag comprises analyzing characteristics of multi-stage voltage sags due to different causes. A method for determining the cause of sudden change time of an amplitude of the multi-stage voltage sag is based on a jump characteristic of the sag amplitude. Calculation methods of a relay protection action matrix and a faulty line set are utilized to preliminarily obtain a faulty line set based on the fault clearing time of a relay protection device and other characteristics to effectively reduce a calculation amount of sag state estimation. Based on a substantive characteristic that different events cause a change of system impedance, a method for inferring a cause of an event in each stage of the multi-stage voltage sag is used to improve the accuracy of the state estimation. The state estimation method reduces the difficulty in applying to the multi-stage voltage sag.

Claims

1. A state estimation method for a multi-stage voltage sag, comprising the following steps: step 1: analyzing causes of a multi-stage voltage sag, determining characteristics of multi-stage voltage sags due to different causes, and identifying causes of voltage sags at different sudden change points accordingly; step 2: building an action matrix of a relay protection device based on a specified fault clearing time at different positions of a power grid; step 3: obtaining one or more possible faulty lines through solving to form a faulty line set to reduce a calculation amount of state estimation of the multi-stage voltage sag; step 4: calculating system impedance matrices before and after a fault, an impedance matrix in the case of a relay protection action on one side of the faulty line, and an impedance matrix when a distributed generation (DG) device is disconnected from the power grid; step 5: using an optimization model to infer a specific faulty line from the faulty line set obtained in step 3, and inferring a specific event of each stage of the multi-stage voltage sag based on a voltage sag type; and step 6: estimating duration of each sag stage based on the voltage sag type and a time length between adjacent sudden voltage sag change points, and estimating a sag amplitude based on an event inference result in step 5.

2. The state estimation method for a multi-stage voltage sag according to claim 1, wherein the causes and corresponding characteristics of the multi-stage voltage sag are as follows: cause I: during a short circuit fault, one or more DG devices in a power system trip, resulting in a loss of a power supply of the power grid; and a corresponding characteristic is that a voltage amplitude on a right side of a corresponding sudden voltage sag change point when the voltage sag occurs is smaller than that on a left side; and cause II: during the short circuit fault, relay protection devices on two sides of the faulty line trip at different time, resulting in a change of a topology of the power grid; and a corresponding characteristic is that the voltage amplitude on the right side of the corresponding sudden voltage sag change point when the voltage sag occurs is greater than that on the left side.

3. The state estimation method for a multi-stage voltage sag according to claim 2, wherein the identifying causes of voltage sags at different sudden change points specifically comprises: assuming that for a waveform segment containing the voltage sag, a total of s + 1 sudden voltage sag change points are detected by using a waveform point detection method, which are respectively expressed as .sup.MP.sub.0,.sup.MP.sub.1,...,.sup.MP.sub.x,...,.sup.MP.sub.s and correspond to s + 1 sudden change time points t.sub.0,t.sub.1,...,t.sub.x,...,t.sub.s, wherein these sudden voltage sag change points divide a waveform of the voltage sag into .sup.s segments; and voltages of the different waveform segments are respectively expressed as u.sub.1,u.sub.2,...,u.sub.x,...,u.sub.s, and an .sup.m th voltage .sup.u.sub.m is calculated as follows: $u_m={\displaystyle {.Math.}_{j=ff*\left(t_{m-1}+1/f\right)}^{j=ff*t_m}{x}_j}/\left(\left(t_m-t_{m-1}+1/f\right)*ff\right);j\in \left[1,k\right]$ wherein ƒƒ represents a sampling rate in units of piece/second; ƒ represents a power-frequency current frequency in units of Hz; .sup.x.sub.j represents a j .sup.th value in an effective value sequence; m ∈ [1,s] ; and k represents a total quantity of sampling points; and forming a binary vector W according to the above formula: $\left\{\begin{array}{c}W=\left[w_1,w_2,\text{L},w_x,\text{L},w_s\right],x\in \left[1,s\right]\\ w_x=\left\{\begin{array}{l}0,\text{if}{u}_x<u_{x-1}\text{or}x=1\\ 1,\text{if}{u}_x<u_{x-1}\end{array}\right)\end{array}\right)$ wherein a value of the element .sup.w.sub.x in the vector W is 0 or 1; when the value of the .sup.w.sub.x is 0, it indicates that a voltage on a right side of .sup.MP.sub.x-1 corresponding to a time point .sup.t.sub.x-1 is less than that on a left side, which means that the multi-stage voltage sag is formed due to the cause I at this time point, wherein x > 1; on the contrary, when the value of the .sup.w.sub.x is 1, it means that the multi-stage voltage sag is formed due to the cause II at the time point .sup.t.sub.x-1, wherein x > 1.

4. The state estimation method for a multi-stage voltage sag according to claim 3, wherein the building an action matrix of a relay protection device in step 2 specifically comprises: step 2.1: building a basic action matrix, in other words, a universal matrix describing a fault removal behavior of the relay protection device: $\Gamma =\left(\begin{array}{ccccc}0& {\gamma }_{12}& {\gamma }_{13}& \text{L}& {\gamma }_{1n}\\ {\gamma }_{21}& 0& {\gamma }_{23}& \text{L}& {\gamma }_{2n}\\ {\gamma }_{31}& {\gamma }_{32}& 0& \text{O}& \text{M}\\ \text{M}& \text{M}& \text{O}& 0& {\gamma }_{\left(n-1\right)n}\\ {\gamma }_{n1}& {\gamma }_{n2}& \text{L}& {\gamma }_{n\left(n-1\right)}& 0\end{array}\right)$ wherein n represents a quantity of buses in the system, .sup.γ.sub.ij and .sup.γ.sub.ji respectively represent fault clearing time of protection devices close to i -side and j -side buses on a line .sup.l.sub.ij, i, j ∈ [1, n];i ≠ j, and .sup.γ.sub.ij,.sup.γ.sub.ji = 0 represents that there is no physical connection between the i -side bus and the j -side bus; step 2.2: improving the universal matrix and decoupling an improved universal matrix into two matrices: $\Gamma \widetilde{{\rm N}}=\left(\begin{array}{ccccc}0& 0& 0& 0& 0\\ {\gamma }_{21}& 0& 0& 0& 0\\ {\gamma }_{31}& {\gamma }_{32}& 0& 0& 0\\ \text{M}& \text{M}& \text{O}& 0& 0\\ {\gamma }_{n1}& {\gamma }_{n2}& \text{L}& {\gamma }_{n\left(n-1\right)}& 0\end{array}\right)$ $\Gamma \Delta =\left(\begin{array}{ccccc}0& {\gamma }_{12}& {\gamma }_{13}& \text{L}& {\gamma }_{1n}\\ 0& 0& {\gamma }_{23}& \text{L}& {\gamma }_{2n}\\ 0& 0& 0& \text{O}& \text{M}\\ 0& 0& 0& 0& {\gamma }_{\left(n-1\right)n}\\ 0& 0& 0& 0& 0\end{array}\right)$ wherein the upper and lower triangular matrices respectively represent fault clearing time of protective devices of a same type on two terminal buses of the line; a lower triangular element in the lower triangular matrix represents the fault clearing time of the protection device of the i -side bus of the line .sup.l.sub.ij, and .sup.γ.sub.ij.sup.=.sup.0 and .sup.i.sup.>.sup.j represent that there is no physical connection between the i -side bus and the .sup.j -side bus; and an upper triangular element in the upper triangular matrix represents a parameter on the other side of the line; and to represent two-stage protection, two other similar matrices are constructed to represent parameters of backup protection: $\Lambda \widetilde{{\rm N}}=\left(\begin{array}{ccccc}0& 0& 0& 0& 0\\ {\lambda }_{21}& 0& 0& 0& 0\\ {\lambda }_{31}& {\lambda }_{32}& 0& 0& 0\\ \text{M}& \text{M}& \text{O}& 0& 0\\ {\lambda }_{n1}& {\lambda }_{n2}& \text{L}& {\lambda }_{n\left(n-1\right)}& 0\end{array}\right)$ $\Lambda \Delta =\left(\begin{array}{ccccc}0& {\lambda }_{12}& {\lambda }_{13}& \text{L}& {\lambda }_{1n}\\ 0& 0& {\lambda }_{23}& \text{L}& {\lambda }_{2n}\\ 0& 0& 0& \text{O}& \text{M}\\ 0& 0& 0& 0& {\lambda }_{\left(n-1\right)n}\\ 0& 0& 0& 0& 0\end{array}\right)$ step 2.3: determining a cooperation relationship of the protection devices, wherein for the two-stage protection, there are four cooperation relationships between main protection and the backup protection in the line, and the improved protection action matrix is expressed as the following cooperation relationships: $\Gamma \widetilde{{\rm N}}+\Gamma \Delta ={\Theta }_{\text{I}}$ $\Gamma \widetilde{{\rm N}}+\Lambda \Delta ={\Theta }_{\text{II}}$ $\Lambda \widetilde{{\rm N}}+\Lambda \Delta ={\Theta }_{\text{III}}$ $\Lambda \widetilde{{\rm N}}+\Gamma \Delta ={\Theta }_{\text{IV}}$ wherein Θ.sub.I represents that the main protection on the two sides of the faulty line cooperates with each other to remove the fault; Θ.sub.II represents that the backup protection on the two sides of the faulty line cooperates with each other to remove the fault; Θ.sub.III and Θ.sub.IV represent that the main protection on one side of the faulty line cooperates with the backup protection on the other side of the faulty line to remove the fault; and step 2.4: correcting fault clearing time in the action matrix, wherein the fault clearing time is corrected as follows: $\left\{\begin{array}{l}{\gamma }_{ij}={\gamma }_{ij,set}+{\delta }_{1,ij}\\ {\lambda }_{ij}={\lambda }_{ij,set}+{\delta }_{2,ij}\end{array}\right);i,j\in \left[1,n\right];i\ne j$ wherein .sup.γij,set and .sup.ηij,set respectively represent specified values of the main protection and the backup protection; .sup.δ.sub.1,ij and .sup.δ.sub.2,.sub.ij respectively represent deviations between actual fault clearing time and the specified values; and the deviations are random numbers within [0, .sup.δ], wherein δ represents a maximum error value during testing or historical operation of relay protection devices of a same model.

5. The state estimation method for a multi-stage voltage sag according to claim 4, wherein step 3 specifically comprises: step 3.1: defining the following four voltage sag types based on the sag cause identification in step 1: type I: single-stage rectangular sag; type II: multi-stage voltage sag due to the cause I; type III: multi-stage voltage sag due to the cause II; and type IV: multi-stage voltage sag due to the cause I and the cause II; and step 3.2: obtaining the faulty line set for different types of sags, wherein for the type I and the type II: assuming that time points at which first and last sudden voltage sag change points of the waveform of the voltage sag are detected in step 1 are .sup.t.sub.0 and .sup.t.sub.s respectively, a time length from occurrence of the fault to removal of the fault in the system is .sup.t.sub.s-.sup.t.sub.0 ; and if two elements of symmetric positions of a main diagonal line of a matrix in the four matrices in step 2.3 are equal to t.sub.s-t.sub.0 within an error threshold, lines corresponding to these two elements are possible faulty lines, and therefore, a solution model of the faulty line set LF is as follows: $\left\{\begin{array}{l}LF=\left\{\text{L,}{l}_{ij},\text{L}\right\}\\ \text{s}\text{.t}\text{.}LF.Math.LN\\ \sqrt{\left|{\theta }_{ij}-\left(t_s-t_0\right)\right|+\left|{\theta }_{ji}-\left(t_s-t_0\right)\right|}\le \\ \sqrt{\max \left({\delta }_{1,ij},{\delta }_{1,ji}\right)+\max \left({\delta }_{2,ij},{\delta }_{2,ji}\right)}\le \sqrt{2\delta }\\ i,j\in \left[1,n\right];i\ne j\end{array}\right)$ wherein .sup.θ.sub.ij represents an element in an i .sup.th row and a j .sup.th column in the matrices ${\Theta }_{\text{I}}$ to ${\Theta }_{\text{IV}},$ and LN represents a set of lines in an intersection of sag domains of a bus of a monitoring device; and for the type III and the type IV: assuming that the time point at which the first sudden voltage sag change point of the waveform of the voltage sag detected in step 1 is t.sub.0, and time points of sudden voltage sag change points of two corresponding relay protection actions are t.sub.x-1 and t.sub.s respectively, a time length from the occurrence of the fault to the protection actions of the relay protection devices on the two sides of the faulty line are t.sub.x-1-t.sub.0 and t.sub.s-t.sub.0 respectively; and in this case, the solution model of the faulty line set is as follows: ${\left\{\begin{array}{l}LF=\left\{\text{L,}{l}_{ij},\text{L}\right\}\\ \text{s}\text{.t}\text{.}LF.Math.LN\\ \sqrt{\left|{\theta }_{ij}-\left(t_s-t_0\right)\right|+\left|{\theta }_{ji}-\left(t_s-t_0\right)\right|}\le \\ \sqrt{\max \left({\delta }_{1,ij},{\delta }_{1,ji}\right)+\max \left({\delta }_{2,ij},{\delta }_{2,ji}\right)}\le \sqrt{2\delta }\\ i,j\in \left[1,n\right];i\ne j;s\ge 2;w_x=1\end{array}\right)}_{.}$ .

6. The state estimation method for a multi-stage voltage sag according to claim 5, wherein step 4 specifically comprises: step 4.1: calculating the system impedance matrix before the fault, wherein a system admittance matrix Y.sup.SE is expressed as a sum of a line admittance matrix $Y_L^{SE}$ and a generator admittance matrix, $Y_G^{SE}$ as shown below: $Y^{SE}=Y_G^{SE}+Y_L^{SE}$ assuming there are .sup.n buses in the power system, the line admittance matrix $Y_L^{SE}$ is calculated as follows based on a line topology relationship and an impedance parameter: $Y_L^{SE}=\left(\begin{array}{cccc}{{\displaystyle \alpha }}_{11}^{se}& {{\displaystyle \alpha }}_{12}^{se}& \text{L}& {{\displaystyle \alpha }}_{1n}^{se}\\ {{\displaystyle \alpha }}_{21}^{se}& {{\displaystyle \alpha }}_{22}^{se}& \text{L}& {{\displaystyle \alpha }}_{2n}^{se}\\ \text{M}& \text{M}& \text{O}& \text{M}\\ {{\displaystyle \alpha }}_{n1}^{se}& {{\displaystyle \alpha }}_{n2}^{se}& \text{L}& {{\displaystyle \alpha }}_{nm}^{se}\end{array}\right)$ wherein se = 1,2,0 represents a positive sequence, a negative sequence, and a zero sequence, .sup.α.sub.ij represents mutual admittance of nodes i and j, α.sub.ii represents self-admittance of the node .sup.i, and i ≠ j; i, j ∈ [1,n]; the matrix $Y_G^{SE}$ is a diagonal matrix, and an element value on a diagonal line is equal to self-admittance of a corresponding generator, as shown below: $Y_G^{SE}=\left(\begin{array}{cccc}{{\displaystyle \beta }}_{11}^{se}& 0& \text{L}& 0\\ 0& {{\displaystyle \beta }}_{22}^{se}& \text{L}& 0\\ \text{M}& \text{M}& \text{O}& \text{M}\\ 0& 0& \text{L}& {{\displaystyle \beta }}_{nm}^{se}\end{array}\right)$ wherein .sup.β.sub.ii = 0 represents that there is no generator for the bus; and the system impedance matrix is calculated as follows: $Z_N^{SE}={\left(Y_L^{SE}\right)}^{-1}=\left(\begin{array}{cccc}{{\displaystyle Z}}_{11}^{se}& {{\displaystyle Z}}_{12}^{se}& \text{L}& {{\displaystyle Z}}_{1n}^{se}\\ {{\displaystyle Z}}_{21}^{se}& {{\displaystyle Z}}_{22}^{se}& \text{L}& {{\displaystyle Z}}_{2n}^{se}\\ \text{M}& \text{M}& \text{O}& \text{M}\\ {{\displaystyle Z}}_{n1}^{se}& {{\displaystyle Z}}_{n2}^{se}& \text{L}& {{\displaystyle Z}}_{nm}^{se}\end{array}\right)$ step 4.2: calculating the system impedance matrix after the short circuit fault, wherein mutual impedance $Z_{mf}^{se}$ between a fault position ƒ .sub.l and a to-be-solved node .sup.m after the circuit short fault, and self-impedance $Z_{ff}^{se}$ of the fault position ƒ .sub.l are respectively obtained through solving according to the following two formulas: $Z_{mf}^{se}=\left(1-p\right)Z_{mi}^{se}+p{Z}_{mj}^{se}=g_{mf1}\left(i,j,p\right)$ $Z_{ff}^{se}=\left\{\begin{array}{c}{\left(1-p\right)}^2{Z}_{ii}^{se}+p^2{Z}_{jj}^{se}+\\ 2p\left(1-p\right)Z_{ij}^{se}+p\left(1-p\right)z_{ij}^{se}\end{array}\right\}=g_{jj1}\left(i,j,p\right)$ wherein $Z_{mi}^{se}$ , $Z_{mj}^{se}$ , $Z_{ij}^{se}$ , $Z_{ii}^{se}$ , and $Z_{jj}^{se}$ are elements in the matrix $Z_N^{SE}$ ; and $z_{ij}^{se}$ represents impedance of the line .sup.l.sub.ij, wherein the impedance is further expressed as functions g.sub.mf1(i,j,p) and g.sub.ff1(i,j,p), and the two functions respectively represent mutual impedance between the fault position and the target node .sup.m and self-impedance of the fault position when the short circuit fault occurs at the position ƒ.sub.l having a distance P away from the i -side of the line l.sub.ij; step 4.3: calculating the impedance matrix in the case of the relay protection action on one side of the faulty line, wherein assuming that after the short circuit fault occurs on the line .sup.l.sub.ij, a j -side protection device acts to cut off the line on the corresponding side, and the system impedance is calculated as follows: first, a branch with impedance of .sup.-z.sub.ij is added between the i -side bus and the j -side bus in the original system, and in this case, the system impedance matrix $Z_{AT}^{SE}$ is corrected according to the following formula: $\begin{array}{c}{Z}_{AT}^{SE}=Z_N^{SE}-\frac{\text{Δ}Z\text{Δ}{Z}^T}{-z_{ij}^{se}+Z_{ii}^{se}+Z_{jj}^{se}-2Z_{ij}^{se}}\\ =\left(\begin{array}{cccc}{Z}_{11,\text{AT}}^{se}& Z_{12,\text{AT}}^{se}& \text{L}& Z_{1n,\text{AT}}^{se}\\ Z_{21,\text{AT}}^{se}& Z_{22,\text{AT}}^{se}& \text{L}& Z_{2n,\text{AT}}^{se}\\ \text{M}& \text{M}& \text{O}& \text{M}\\ Z_{n1,\text{AT}}^{se}& Z_{n2,\text{AT}}^{se}& \text{L}& Z_{nn,\text{AT}}^{se}\end{array}\right)\end{array}$ wherein ΔZ represents a process quantity, which is calculated according to the following formula: $\text{Δ}Z={\left[\begin{array}{cccc}{Z}_{1j}^{se}-Z_{1i}^{se}& Z_{2j}^{se}-Z_{2i}^{se}& \text{L}& Z_{nj}^{se}-Z_{ni}^{se}\end{array}\right]}^{\text{T}}$ then a branch with impedance of .sup.pz.sub.ij is added to the i -side bus, and in this case, the system impedance matrix is further corrected as $Z_{ARP}^{SE}$ according to the following formula: $Z_{ARP}^{SE}=\left(\begin{array}{ccccc}{Z}_{11,AT}^{se}& Z_{12,AT}^{se}& \text{L}& Z_{1n,AT}^{se}& Z_{1i,AT}^{se}\\ Z_{21,AT}^{se}& Z_{22,AT}^{se}& \text{L}& Z_{2n,AT}^{se}& Z_{2i,AT}^{se}\\ \text{M}& \text{M}& \text{O}& \text{M}& \text{M}\\ Z_{n1,AT}^{se}& Z_{n2,AT}^{se}& \text{L}& Z_{nn,AT}^{se}& Z_{ni,AT}^{se}\\ Z_{i1,AT}^{se}& Z_{i2,AT}^{se}& \text{L}& Z_{in,AT}^{se}& Z_{ii,AT}^{se}+p{z}_{ij}\end{array}\right)$ compared with the $Z_{AT}^{SE}$ , the $Z_{ARP}^{SE}$ adds one line and one column to represent the mutual impedance $Z_{mf}^{se}$ between the fault position ƒ .sub.l and each target bus .sup.m or the self-impedance $Z_{ff}^{se}$ of the fault position; and the impedance is further expressed as functions g .sub.mf 2 (i,j,p,d) and g .sub.ff 2 (i, j, p, d), wherein the two functions respectively represent mutual impedance between the fault position and the target node .sup.m and self-impedance of the fault position when a d -side protection device on the line acts to cut off a part of the line after the short circuit fault occurs at the position ƒ.sub.l far from the terminal P of the node i on the line .sup.l.sub.ij, as shown in the following formulas: $Z_{mf}^{se}=Z_{im,\text{AT}}^{se}=g_{mf2}\left(i,j,p,d\right);d=i\text{or}j$ $Z_{ff}^{se}=Z_{ii,\text{AT}}^{se}+p{z}_{ij}=g_{ff2}\left(i,j,p,d\right);d=i\text{or}j$ step 4.4: calculating the impedance matrix when the DG device is disconnected from the power grid, wherein assuming that DG devices of all buses in a bus set h are disconnected from the power grid, the diagonal matrix first needs to be corrected as follows: $Y_G^{SE}=\left(\begin{array}{cccc}{\beta }_{11}^{se}& 0& \text{L}& 0\\ 0& {\beta }_{22}^{se}& \text{L}& 0\\ \text{M}& \text{M}& \text{O}& \text{M}\\ 0& 0& \text{L}& {\beta }_{nn}^{se}\end{array}\right);{\beta }_{ii}^{se}=0;\forall i\in h$ the system impedance matrix before the fault is calculated according to the above formula, and finally system impedance when the DG device is disconnected from the power grid in the case of the short circuit fault is calculated; and the impedance is further expressed as functions g.sub.mf3 (i, j, p, h) and g.sub.ff3 (i, j, p, h), wherein the two functions respectively represent mutual impedance between the fault position and the target node .sup.m and self-impedance of the fault position after the short circuit fault occurs at the position ƒ.sub.l far from the terminal P of the node i on the line l.sub.ij and the DG devices of all the buses in the bus set h are disconnected from the power grid, as shown in the following two formulas: $Z_{mf}^{se}=g_{mf3}\left(i,j,p,h\right)$ $Z_{ff}^{se}=g_{ff3}\left(i,j,p,h\right).$ .

7. The state estimation method for a multi-stage voltage sag according to claim 6, wherein in step 5, the following four optimization models are available based on the sag types: (1) for the type I, the following optimization model is used to infer the faulty line and its short circuit condition, in other words, i, j, and P ; $\left\{\begin{array}{l}\max \text{-}\sqrt{\left|f\left(g_{mf1}\left(i,j,p\right),g_{ff1}\left(i,j,p\right)\right)\right|-u_1}\hfill \\ \text{s}\text{.t}\text{.}i,j\in \left[1,n\right];i\ne j\hfill \end{array}\right),$ wherein ƒ (.Math.) represents a short-circuit calculation function; (2) for the type II, the following optimization model is used to infer the faulty line and its short circuit condition, a time point and a sequence of disconnecting the DG device from the power grid, in other words, i, j, p, and .sup.h.sub.q ; $\left\{\begin{array}{l}\max \text{-}\left(\begin{array}{c}\sqrt{\left|f\left(g_{mf1}\left(i,j,p\right),g_{ff1}\left(i,j,p\right)\right)\right|-u_1}+\\ {\displaystyle \underset{x=2}{\overset{s}{.Math.}}{\displaystyle \underset{q=1}{\overset{s-1}{.Math.}}\sqrt{\left|f\left(g_{mf3}\left(i,j,p,h_q\right),g_{ff3}\left(i,j,p,h_q\right)\right)\right|-u_x}}}\end{array}\right)\hfill \\ \text{s}\text{.t}\text{.}i,j\in \left[1,n\right];i\ne j;h_1.Math.h_2.Math..Math..Math.h_{s\text{-1}}\hfill \end{array}\right)$ wherein .sup.s represents a quantity of stages of the multi-stage voltage sag; and .sup.h.sub.q represents a set of DG devices disconnected from the power grid during a q-stage voltage sag; (3) for the type III, the following optimization model is used to infer the faulty line and its short circuit condition, and a tripping time point and sequence of the relay protection device, in other words, i, j, P, and d ; $\left\{\begin{array}{l}\max \text{-}\left(\begin{array}{c}\sqrt{\left|f\left(g_{mf1}\left(i,j,p\right),g_{ff1}\left(i,j,p\right)\right)\right|-u_1}+\\ \sqrt{\left|f\left(g_{mf2}\left(i,j,p,d\right),g_{ff2}\left(i,j,p,d\right)\right)\right|-u_2}\end{array}\right)\hfill \\ \text{s}\text{.t}\text{.}i,j\in \left[1,n\right];i\ne j\hfill \end{array}\right)$ (4) for the type IV, the following optimization model is used to infer the faulty line and its short circuit condition, the tripping time point and sequence of the relay protection device, and the time point and the sequence of disconnecting the DG device from the power grid, in other words, i, j, p, d, and .sup.h.sub.q; $\left\{\begin{array}{l}\max \text{-}\left(\begin{array}{c}\begin{array}{c}\sqrt{\left|f\left(g_{mf1}\left(i,j,p\right),g_{ff1}\left(i,j,p\right)\right)\right|-u_1}+\\ {\displaystyle \underset{x=2}{\overset{o}{.Math.}}{\displaystyle \underset{q=1}{\overset{o-1}{.Math.}}\sqrt{\left|f\left(g_{mf3}\left(i,j,p,h_q\right),g_{ff3}\left(i,j,p,h_q\right)\right)\right|-u_x}}}\end{array}\\ \begin{array}{c}+\left|f\left(g_{mf2}\left(i,j,p,d\right),g_{ff2}\left(i,j,p,d\right)\right)\right|-u_{o+1}+\\ {\displaystyle \underset{x=o+2}{\overset{s}{.Math.}}{\displaystyle \underset{q=o}{\overset{s-2}{.Math.}}\sqrt{\left|f\left(g_{mf3}\left(i,j,p,h_q\right),g_{ff3}\left(i,j,p,h_q\right)\right)\right|-u_x}}}\end{array}\end{array}\right)\hfill \\ \text{s}\text{.t}\text{.}i,j\in \left[1,n\right];i\ne j;h_1.Math.h_2.Math.\text{L}.Math.h_{s\text{-2}}\hfill \end{array}\right)$ wherein .sup.o represents a time point at which a relay protection device on a side of the faulty line first acts to trip; and .sup.u.sub.o+1 represents an amplitude of an o+1-stage voltage sag.

8. The state estimation method for a multi-stage voltage sag according to claim 7, wherein in step 6, the following four voltage sag state estimation methods are available based on the sag types: 1) for the type I, a sag amplitude of any unmonitored bus .sup.m is estimated according to the following formula: $f\left(g_{mf1}\left(i,j,p\right),g_{ff1}\left(i,j,p\right)\right)$ 2) for the type II, a sag amplitude of any unmonitored bus .sup.m is estimated according to the following formula: $\left\{\begin{array}{l}{u}_1=f\left(g_{mf1}\left(i,j,p\right),g_{ff1}\left(i,j,p\right)\right)\hfill \\ u_x=f\left(g_{mf3}\left(i,j,p,h_{x-1}\right),g_{ff3}\left(i,j,p,h_{x-1}\right)\right)\hfill \end{array}\right)$ 3) for the type III, a sag amplitude of any unmonitored bus .sup.m is estimated according to the following formula: $\left\{\begin{array}{l}{u}_1=f\left(g_{mf1}\left(i,j,p\right),g_{ff1}\left(i,j,p\right)\right)\hfill \\ u_2=f\left(g_{mf2}\left(i,j,p,d\right),g_{ff2}\left(i,j,p,d\right)\right)\hfill \end{array}\right)$ 4) for the type IV, a sag amplitude of any unmonitored bus .sup.m is estimated according to the following formula: $\left\{\begin{array}{l}{u}_1=f\left(g_{mf1}\left(i,j,p\right),g_{ff1}\left(i,j,p\right)\right)\hfill \\ u_x=f\left(g_{mf3}\left(i,j,p,h_{x-1}\right),g_{ff3}\left(i,j,p,h_{x-1}\right)\right)\hfill \\ u_y=f\left(g_{mf2}\left(i,j,p,d\right),g_{ff2}\left(i,j,p,d\right)\right)\hfill \\ u_z=f\left(g_{mf3}\left(i,j,p,h_{z-1}\right),g_{ff3}\left(i,j,p,h_{z-1}\right)\right)\hfill \\ x\in \left[2,y-1\right];z\in \left[y+1,s-2\right]\hfill \end{array}\right)$ wherein .sup.y and z respectively represent y.sup.th and z.sup.th stages of the multi-stage voltage sag.

Description

BRIEF DESCRIPTION OF THE DRAWINGS

[0094] FIG. 1 shows a short-circuit calculation model of a power system;

[0095] FIG. 2 shows a short-circuit calculation model of a power system after a faulty line is cut off on one side according to the present disclosure; and

[0096] FIG. 3 is a process of building a state estimation model of a multi-stage voltage sag according to the present disclosure.

DETAILED DESCRIPTION OF THE EMBODIMENTS

[0097] The present disclosure will be further described below in conjunction with the accompanying drawings and specific embodiments.

[0098] The technical solution of the present disclosure is mainly divided into six steps: identifying causes of a multi-stage voltage sag, building an action matrix of a relay protection device, obtaining a faulty line set by solving equations, obtaining an impedance matrix by solving equations, and inferring a specific event at each stage of the multi-stage voltage sag and a state estimation method for the multi-stage voltage sag. Details of each step and its substeps are as follows:

[0099] 1. The causes of the multi-stage voltage sag are identified.

There Are Two Causes of the Multi-Stage Voltage Sag

[0100] Cause 1: During a short circuit fault, one or more DG device are in a power system trip, resulting in a loss of a power supply from a power grid.

[0101] Cause 2: During the short circuit fault, relay protection devices on two sides of a faulty line trip at different times, resulting in a change of the topology of the power grid.

Characteristics of Multi-Stage Voltage Sags Due to Different Causes

[0102] Cause 1 means that after the fault, the system further loses the power supply, so a monitored voltage sag amplitude is further reduced. The corresponding characteristic is that a voltage amplitude on a right side of the corresponding sudden voltage sag change point when cause 1 occurs is smaller than that on a left side.

[0103] On the contrary, cause 2 means that after the fault occurs, a protection device on one side of the line partially isolates the fault from the system, so the monitored voltage sag amplitude increases. The corresponding characteristic is that a voltage amplitude on the right side of the corresponding sudden voltage sag change point when the cause 2 occurs is greater than that on the left side.

Cause Identification of Voltage Sags at Different Sudden Change Points

[0104] It is assumed that for a waveform segment containing the voltage sag, a total of s+1 sudden voltage sag change points are detected by using a waveform point detection method, which are respectively expressed as MP.sub.0,MP.sub.1,...,MP.sub.x,...,MP.sub.s and correspond to s +1 sudden change time points t.sub.0, t.sub.1,..., t.sub.x,..., t.sub.s . These sudden voltage sag change points divide a waveform of the voltage sag into .sup.S segments. Voltages of different waveform segments are respectively expressed as .sup.u.sub.1, .sup.u.sub.2, .. .,.sup.u.sub.x,.. ., u.sub.s, and an m (m ∈ [1, s]) .sup.th voltage u.sub.m can be calculated as follows:

[00059]$\begin{array}{cc}{u}_m={\displaystyle {.Math.}_{j=ff*\left(t_{m-1}+1/f\right)}^{j=ff*t_m}{x}_j}/\left(\left(t_m-t_{m-1}+1/f\right)*ff\right);j\in \left[1,k\right]& \text{­­­(7)}\end{array}$

[0105] In the above formula, ff represents a sampling rate in units of piece/second, f represents a power-frequency current frequency in units of Hz, and .sup.x.sub.j represents a j .sup.th value in an effective value sequence. Binary vector W shown in formula (8) can be formed according to the formula (7):

[00060]$\begin{array}{cc}\left\{\begin{array}{l}W=\left[w_1,w_2,\text{L}\text{,}{w}_x,\text{L,}{w}_s\right],x\in \left[1,s\right]\hfill \\ w_x=\left\{\begin{array}{l}0,\text{if}{u}_x<u_{x-1}\text{or}x=1\hfill \\ 1,\text{if}{u}_x>u_{x-1}\hfill \end{array}\right)\hfill \end{array}\right)& \text{­­­(8)}\end{array}$

[0106] The value of the element w.sub.x in the vector W is 0 or 1. When the value is 0, it indicates that a voltage on a right side of MP.sub.x-1 corresponding to time point t.sub.x.sub.-1 (x > 1) is less than that on a left side, which means that the multi-stage voltage sag is formed due to cause 1 at this time point. On the contrary, when the value is 1, it means that the multi-stage voltage sag is formed due to cause 2 at the time point t.sub.x-1 (x > 1).

[0107] 2. Action matrix of the relay protection device.

[0108] Actual tripping time of the relay protection device in a sag event can be roughly obtained based on multi-stage voltage sag causes at time points corresponding to different sudden voltage sag change points. In order to use the information, a method for building the action matrix of the relay protection device based on specified fault clearing time at different positions of the power grid is proposed in step 2.

Basic Action Matrix

[0109] There are different types of relay protection, including low-voltage protection, over-current protection, and distance protection. Each type of protection has its own trigger rule and fault clearing time to clear the fault. The present disclosure provides a universal matrix for describing a fault removal behavior of the relay protection device, as shown in the following formula:

[00061]$\begin{array}{cc}\text{Γ=}\left(\begin{array}{ccccc}0& {\gamma }_{12}& {\gamma }_{13}& \text{L}& {\gamma }_{1n}\\ {\gamma }_{21}& 0& {\gamma }_{23}& \text{L}& {\gamma }_{2n}\\ {\gamma }_{31}& {\gamma }_{32}& 0& \text{O}& \text{M}\\ \text{M}& \text{M}& \text{O}& 0& {\gamma }_{\left(n-1\right)n}\\ {\gamma }_{n1}& {\gamma }_{n2}& \text{L}& {\gamma }_{n\left(n-1\right)}& 0\end{array}\right)& \text{­­­(9)}\end{array}$

[0110] Where n represents a quantity of buses in the system, .sup.γ.sub.ij and .sup.γ.sub.ji (i, j ∈ [1, n]; i ≠ j) respectively represent fault clearing time of protection devices close to i -side and j -side buses on line l.sub.ij in units of ms, and .sup.γ.sub.ij, .sup.γ.sub.ji= 0 represents that there is no physical connection between the i -side bus and the j -side bus. The fault clearing time is a specified value of the relay protection device recorded in the software of a power distribution automation system and other power enterprises.

The Action Matrix is Improved

[0111] In practical application, staged protection is widely used in a power distribution system, such as two-stage protection. The two-stage protection means that two types of protection are provided on one line, namely, main protection and backup protection. The two types of protection have different fault clearing time, and the main protection has shorter fault clearing time. When the short circuit fault occurs, the main protection is triggered first, and the backup protection works when the main protection cannot be triggered. In general, the fault can be removed by the main protection independently, and sometimes may be removed by the backup protection independently or by both the main protection and the backup protection. Considering that it is difficult for the basic action matrix defined in the formula (9) to present a cooperation relationship of the staged protection, the present disclosure improves the basic action matrix and decouples an improved basic action matrix into two matrices, as shown in equations (10) and (11):

[00062]$\begin{array}{cc}\text{Γ}\tilde{\text{N}}\text{=}\left(\begin{array}{ccccc}0& 0& 0& 0& 0\\ {\gamma }_{21}& 0& 0& 0& 0\\ {\gamma }_{31}& {\gamma }_{32}& 0& 0& 0\\ \text{M}& \text{M}& \text{O}& 0& 0\\ {\gamma }_{n1}& {\gamma }_{n2}& \text{L}& {\gamma }_{n\left(n-1\right)}& 0\end{array}\right)& \text{­­­(10)}\end{array}$

[00063]$\begin{array}{cc}\text{ΓΔ=}\left(\begin{array}{ccccc}0& {\gamma }_{12}& {\gamma }_{13}& \text{L}& {\gamma }_{1n}\\ 0& 0& {\gamma }_{23}& \text{L}& {\gamma }_{2n}\\ 0& 0& 0& \text{O}& \text{M}\\ 0& 0& 0& 0& {\gamma }_{\left(n-1\right)n}\\ 0& 0& 0& 0& 0\end{array}\right)& \text{­­­(11)}\end{array}$

[0112] In the above formulas, the upper and lower triangular matrices respectively represent fault clearing time of protection devices of the same type on two terminal buses of the line. The lower triangular element in the lower triangular matrix represents the fault clearing time of the protection device of the i -side bus of the line l.sub.ij, and .sup.γ.sub.ij = 0(i > j) represents that there is no physical connection between the i -side bus and the j -side bus. The upper triangular element in the upper triangular matrix represents a parameter on the other side of the line, which has a similar meaning. The formulas (10) and (11) represent parameters of the main protection. To represent the two-stage protection, two other similar matrices are constructed to represent parameters of the backup protection, as shown in formulas (12) and (13):

[00064]$\begin{array}{cc}\text{Λ}\tilde{\text{N}}=\left(\begin{array}{ccccc}0& 0& 0& 0& 0\\ {\lambda }_{21}& 0& 0& 0& 0\\ {\lambda }_{31}& {\lambda }_{32}& 0& 0& 0\\ \text{M}& \text{M}& \text{O}& 0& 0\\ {\lambda }_{n1}& {\lambda }_{n2}& \text{L}& {\lambda }_{n\left(n-1\right)}& 0\end{array}\right)& \text{­­­(12)}\end{array}$

[00065]$\begin{array}{cc}\text{ΛΔ}=\left(\begin{array}{ccccc}0& {\lambda }_{12}& {\lambda }_{13}& \text{L}& {\lambda }_{1n}\\ 0& 0& {\lambda }_{23}& \text{L}& {\lambda }_{2n}\\ 0& 0& 0& \text{O}& \text{M}\\ 0& 0& 0& 0& {\lambda }_{\left(n-1\right)n}\\ 0& 0& 0& 0& 0\end{array}\right)& \text{­­­(13)}\end{array}$

Cooperation Relationship of the Protection Devices

[0113] For the widely used two-stage protection, there are four cooperation relationships between the main protection and the backup protection in the line. The improved protection action matrix can conveniently express these cooperation relationships as follows:

[00066]$\begin{array}{cc}\text{Γ}\tilde{\text{N}}\text{+}\text{ΓΔ=}{\text{Θ}}_{\text{I}}& \text{­­­(14)}\end{array}$

[00067]$\begin{array}{cc}\text{Γ}\tilde{\text{N}}\text{+}\text{ΛΔ=}{\text{Θ}}_{\text{II}}& \text{­­­(15)}\end{array}$

[00068]$\begin{array}{cc}\text{Λ}\tilde{\text{N}}\text{+}\text{ΛΔ=}{\text{Θ}}_{\text{III}}& \text{­­­(16)}\end{array}$

[00069]$\begin{array}{cc}\text{Λ}\tilde{\text{N}}\text{+}\text{ΓΔ=}{\text{Θ}}_{\text{IV}}& \text{­­­(17)}\end{array}$

[0114] Relationship 1: As shown in the formula (14), Θ.sub.I represents that the main protection on two sides of the faulty line cooperates with each other to remove the fault.

[0115] Relationship 2: As shown in the formula (15), Θ.sub.II represents that the backup protection on the two sides of the faulty line cooperates with each other to remove the fault.

[0116] Relationship ¾: As shown in the formulas (16) and (17), Θ.sub.III and Θ.sub.IV represent that the main protection on one side of the faulty line cooperates with the backup protection on the other side of the faulty line to remove the fault.

Fault Clearing Time in the Action Matrix is Corrected

[0117] The fault clearing time used in the formulas (9) to (17) is a specified value obtained from a power company. Theoretically, the fault clearing time of the main protection or the backup protection on the two sides of the line should be the same and equal to the specified value. However, in actual operation, due to error, the elements in these matrices may deviate from specified values, so necessary corrections need to be made, as shown in formula (18):

[00070]$\begin{array}{cc}\left\{\begin{array}{l}{\gamma }_{ij}={\gamma }_{ij,set}+{\delta }_{1,ij}\hfill \\ {\lambda }_{ij}={\lambda }_{ij,set}+{\delta }_{2,ij}\hfill \end{array}\right);i,j\in \left[1,n\right];i\ne j& \text{­­­(18)}\end{array}$

[0118] In the above formula, .sup.γ.sub.ij,.sub.set and .sup.λ.sub.ij,.sub.set respectively represent specified values of the main protection and the backup protection. .sup.δ.sub.1,ij and .sup.δ.sub.2,.sub.ij respectively represent deviations between actual fault clearing time and the specified values. The deviations are random numbers within [0, δ ], where δ represents a maximum error value during the testing or historical operation of relay protection devices of the same model.

[0119] 3. The faulty line set is obtained by solving equations.

[0120] One or more possible faulty lines can be obtained by solving based on the sag cause identification in step 1 and the content of the protection action matrix in step 2 to reduce the calculation amount of the state estimation of the multi-stage voltage sag.

Sag Type

[0121] The following four voltage sag types are defined based on the sag cause identification in step 1: [0122] Type 1: single-stage rectangular sag, [0123] Type 2: multi-stage voltage sag due to the cause 1, [0124] Type 3: multi-stage voltage sag due to the cause 2, and [0125] Type 4: multi-stage voltage sag due to the cause 1 and the cause 2.

The Faulty Line Set is Obtained for Different Types of Sags by Solving Equations

[0126] Type 1 and 2: In these two types of voltage sag events, the relay protection devices on the two sides of the faulty line act simultaneously at a time point of the last sudden voltage sag change to remove the fault. Assuming that time points at which first and last sudden voltage sag change points of the waveform of the voltage sag are detected in step 1 are t.sub.0 and t.sub.s respectively, the time length from the occurrence of the fault to removal of the fault in the system is t.sub.s-t.sub.0 . If two elements of symmetric positions of the main diagonal line of a matrix in formulas (14) to (17) are equal to t.sub.s-t.sub.0 within an error threshold, lines corresponding to these two elements are possible faulty lines. Therefore, a solution model of the faulty line set LF is as follows:

[00071]$\begin{array}{cc}\left\{\begin{array}{l}LF=\left\{\text{L},l_{ij},\text{L}\right\}\\ \text{s}\text{.t}\text{.}LF.Math.LN\\ \sqrt{\left|{\theta }_{ij}-\left(t_s-t_0\right)\right|+\left|{\theta }_{ji}-\left(t_s-t_0\right)\right|}\le \\ \sqrt{\max \left({\delta }_{1,ij},{\delta }_{1,ji}\right)+\max \left({\delta }_{2,ij},{\delta }_{2,ji}\right)}\le \sqrt{2\delta }\\ i,j\in \left[1,n\right];i\ne j\end{array}\right)& \text{­­­(19)}\end{array}$

[0127] Where .sup.θ.sub.ij represents an element in an i .sup.th row and a j .sup.th column in the matrices Θ.sub.I to Θ.sub.IV and LN represents a set of lines in the intersection of sag domains of a bus of a monitoring device.

[0128] Type 3 and 4: In these two types of voltage sag events, the relay protection devices on the two sides of the faulty line cannot act simultaneously to remove the fault. It is assumed that the time point at which the first sudden voltage sag change point of the waveform of the voltage sag is detected in step 1 and time points of sudden voltage sag change points of two corresponding relay protection actions are t.sub.0, t.sub.x-1, and t.sub.s respectively. The time length from the occurrence of the fault to the protection actions of the relay protection devices on the two sides of the faulty line are t.sub.x-1-t.sub.0 and t.sub.s -t.sub.0, respectively. In this case, the solution model of the faulty line set is as follows:

[00072]$\begin{array}{cc}\left\{\begin{array}{l}LF=\left\{\text{L},l_{ij},\text{L}\right\}\\ \text{s}\text{.t}\text{.}LF.Math.LN\\ \sqrt{\left|{\theta }_{ij}-\left(t_s-t_0\right)\right|+\left|{\theta }_{ji}-\left(t_s-t_0\right)\right|}\le \\ \sqrt{\max \left({\delta }_{1,ij},{\delta }_{1,ji}\right)+\max \left({\delta }_{2,ij},{\delta }_{2,ji}\right)}\le \sqrt{2\delta }\\ i,j\in \left[1,n\right];i\ne j\end{array}\right)& \text{­­­(20)}\end{array}$

[0129] 4. The impedance matrix is obtained by solving equations.

[0130] Whether caused by cause 1 or cause 2, the multi-stage voltage sag is essentially a change of a system impedance matrix due to an event during the fault.

System Impedance Matrix Before the Fault

[0131] System admittance matrix Y.sup.SE may be expressed as the sum of line admittance matrix

[00073]$Y_L^{SE}$

and generator admittance matrix

[00074]$Y_G^{SE}$

, as shown below:

[00075]$\begin{array}{cc}{Y}^{SE}=Y_G^{SE}+Y_L^{SE}& \text{­­­(21)}\end{array}$

[0132] Assuming there are .sup.n buses in the power system, the line admittance matrix

[00076]$Y_L^{SE}$

can be calculated as follows based on a line topology relationship and an impedance parameter:

[00077]$\begin{array}{cc}{Y}_L^{SE}=\left(\begin{array}{cccc}{\alpha }_{11}^{se}& {\alpha }_{12}^{se}& \text{L}& {\alpha }_{1n}^{se}\\ {\alpha }_{21}^{se}& {\alpha }_{22}^{se}& \text{L}& {\alpha }_{2n}^{se}\\ \text{M}& \text{M}& \text{O}& \text{M}\\ {\alpha }_{n1}^{se}& {\alpha }_{n2}^{se}& \text{L}& {\alpha }_{nn}^{se}\end{array}\right)& \text{­­­(22)}\end{array}$

[0133] In the above formula, se=1,2,0 represents a positive sequence, a negative sequence, and a zero sequence, α .sub.ij (i ≠ j; i, j ∈ [1,n]) represents mutual admittance of nodes i and j, and α.sub.ii (i ∈ [1, n]) represents self-admittance of the node i.

[0134] The matrix

[00078]$Y_G^{SE}$

is a diagonal matrix, and an element value on a diagonal line is equal to self-admittance of a corresponding generator, as shown below:

[00079]$\begin{array}{cc}{Y}_G^{SE}=\left(\begin{array}{cccc}{\beta }_{11}^{se}& 0& \text{L}& 0\\ 0& {\beta }_{22}^{se}& \text{L}& 0\\ \text{M}& \text{M}& \text{O}& \text{M}\\ 0& 0& \text{L}& {\beta }_{nn}^{se}\end{array}\right)& \text{­­­(23)}\end{array}$

where β.sub.ii = 0(i ∈ [1, n]) represents that there is no generator for the bus.

[0135] The system impedance matrix can be calculated as follows according to the formula (21):

[00080]$\begin{array}{cc}{Z}_N^{SE}={\left(Y_G^{SE}\right)}^{-1}=\left(\begin{array}{cccc}{Z}_{11}^{se}& Z_{12}^{se}& \text{L}& Z_{1n}^{se}\\ Z_{21}^{se}& Z_{22}^{se}& \text{L}& Z_{2n}^{se}\\ \text{M}& \text{M}& \text{O}& \text{M}\\ Z_{n1}^{se}& Z_{n2}^{se}& \text{L}& Z_{nn}^{se}\end{array}\right)& \text{­­­(24)}\end{array}$

System Impedance Matrix After the Short Circuit Fault

[0136] When the short circuit fault occurs in the system, system impedance can be calculated using a short-circuit calculation model of a complex power system shown in FIG. 1. Mutual impedance

[00081]$Z_{mf}^{se}$

between fault position f.sub.l and to-be-solved node .sup.m after the circuit short fault and self-impedance

[00082]$Z_{ff}^{se}$

of the fault position f.sub.l can be respectively obtained by solving the following two formulas:

[00083]$\begin{array}{cc}{Z}_{mf}^{se}=\left(1-p\right)Z_{mi}^{se}+p{Z}_{mj}^{se}=g_{mf1}\left(i,j,p\right)& \text{­­­(25)}\end{array}$

[00084]$\begin{array}{cc}{Z}_{ff}^{se}=\left\{\begin{array}{c}{\left(1-p\right)}^2{Z}_{ii}^{se}+p^2{Z}_{jj}^{se}+\\ 2p\left(1-p\right)Z_{ij}^{se}+p\left(1-p\right)z_{ij}^{se}\end{array}\right\}=g_{ff1}\left(i,j,p\right)& \text{­­­(26)}\end{array}$

[0137] Where

[00085]$Z_{mi}^{se},Z_{mj}^{se},Z_{ij}^{se},Z_{ii}^{se},\text{and}{Z}_{jj}^{se}$

are elements in the matrix

[00086]$Z_N^{SE};\text{and}{z}_{ij}^{se}$

represents impedance of the line .sup.l.sub.ij, where the impedance is further expressed as functions g.sub.mf1 (i, j,p) and g .sub.ff1 (i, j,p), and the two functions respectively represent mutual impedance between the fault position and the target node .sup.m and self-impedance of the fault position when the short circuit fault occurs at the position .sup.f.sub.l far from terminal P of the node i on the line .sup.l.sub.ij .

Impedance Matrix in the Case of a Relay Protection Action on One Side of the Faulty Line

[0138] Assuming that after the short circuit fault occurs on the line l.sub.ij, a j -side protection device acts to cut off the line on the corresponding side, and the system impedance can be calculated by using the model shown in FIG. 2.

[0139] First, a branch with impedance of -z.sub.ij is added between the i -side bus and the j -side bus in the original system. In this case, the system impedance matrix

[00087]$Z_{AT}^{SE}$

can be corrected according to the following formula:

[00088]$\begin{array}{cc}\begin{array}{c}{Z}_{AT}^{SE}=Z_N^{SE}-\frac{\text{Δ}Z\text{Δ}{Z}^T}{-z_{ij}^{se}+Z_{ii}^{se}+Z_{jj}^{se}-2Z_{ij}^{se}}\\ =\left(\begin{array}{cccc}{Z}_{11,\text{AT}}^{se}& Z_{12,\text{AT}}^{se}& \text{L}& Z_{1n,\text{AT}}^{se}\\ Z_{21,\text{AT}}^{se}& Z_{22,\text{AT}}^{se}& \text{L}& Z_{2n,\text{AT}}^{se}\\ \text{M}& \text{M}& \text{O}& \text{M}\\ Z_{n1,\text{AT}}^{se}& Z_{n2,\text{AT}}^{se}& \text{L}& Z_{nn,\text{AT}}^{se}\end{array}\right)\end{array}& \text{­­­(27)}\end{array}$

[0140] Where ΔZ can be calculated according to the following formula:

[00089]$\begin{array}{cc}\text{Δ}Z={\left[\begin{array}{cccc}{Z}_{1j}^{se}-Z_{1i}^{se}& Z_{2j}^{se}-Z_{2i}^{se}& .Math.& Z_{nj}^{se}-Z_{ni}^{se}\end{array}\right]}^T& \text{­­­(28)}\end{array}$

[0141] Then, a branch with impedance of .sup.pz.sub.ij is added to the i -side bus. In this case, the system impedance matrix is further corrected as

[00090]$Z_{ARP}^{SE}$

according to the following formula:

[00091]$\begin{array}{cc}{Z}_{ARP}^{SE}=\left(\begin{array}{ccccc}{Z}_{11,AT}^{se}& Z_{12,AT}^{se}& \text{L}& Z_{1n,AT}^{se}& Z_{1i,AT}^{se}\\ Z_{21,AT}^{se}& Z_{22,AT}^{se}& \text{L}& Z_{2n,AT}^{se}& Z_{2i,AT}^{se}\\ \text{M}& \text{M}& \text{O}& \text{M}& \text{M}\\ Z_{n1,AT}^{se}& Z_{n2,AT}^{se}& \text{L}& Z_{nn,AT}^{se}& Z_{ni,AT}^{se}\\ Z_{i1,AT}^{se}& Z_{i2,AT}^{se}& \text{L}& Z_{in,AT}^{se}& Z_{ii,AT}^{se}+p{z}_{ij}\end{array}\right)& \text{­­­(29)}\end{array}$

[0142] Compared with the

[00092]$Z_{AT}^{SE},Z_{ARP}^{SE}$

adds one line and one column to represent the mutual impedance

[00093]$Z_{mf}^{se}$

between the fault position f and each target bus .sup.m or the self-impedance

[00094]$Z_{ff}^{se}$

of the fault position. The impedance is further expressed as functions g .sub.mf 2 (i, j, p, d) and g.sub.ff 2 (i, j, p, d), where the two functions respectively represent mutual impedance between the fault position and the target node .sup.m and self-impedance of the fault position when a d -side protection device on the line acts to cut off a part of the line after the short circuit fault occurs at the position f far from the terminal P of the node i on the line l.sub.ij, as shown in the following formulas:

[00095]$\begin{array}{cc}{Z}_{mf}^{se}=Z_{im,\text{AT}}^{se}=g_{mf2}\left(i,j,p,d\right);d=i\text{or}j& \text{­­­(30)}\end{array}$

and

[00096]$\begin{array}{cc}{Z}_{ff}^{se}=Z_{ii,\text{AT}}^{se}+p{z}_{ij}=g_{ff2}\left(i,j,p,d\right);d=i\text{or}j& \text{­­­(31)}\end{array}$

Impedance Matrix When the DG Device is Disconnected From the Power Grid

[0143] When the short circuit fault occurs, DG devices on some nodes may be disconnected from the power grid due to insufficient low-voltage ride through capabilities. Assuming that DG devices of all buses in bus set h are disconnected from the power grid, the matrix shown in the formula (23) first needs to be corrected as follows:

[00097]$\begin{array}{cc}{Y}_G^{SE}=\left(\begin{array}{cccc}{\beta }_{11}^{se}& 0& \text{L}& 0\\ 0& {\beta }_{22}^{se}& \text{L}& 0\\ \text{M}& \text{M}& \text{O}& \text{M}\\ 0& 0& \text{L}& {\beta }_{nn}^{se}\end{array}\right);{\beta }_{ii}^{se}=0;\forall i\in h& \text{­­­(32)}\end{array}$

[0144] Then, the system impedance matrix before the fault is calculated by substituting the above formula into the formula (24). Finally, system impedance when the DG device is disconnected from the power grid in the case of the short circuit fault is calculated according to the formulas (1) to (26). The impedance is further expressed as functions g .sub.mf 3 (i, j, p, h) and g .sub.ff.sub.3 (i, j, p,h), where the two functions respectively represent mutual impedance between the fault position and the target node .sup.m and self-impedance of the fault position after the short circuit fault occurs at the position f far from the terminal P of the node i on the line l.sub.ij, and the DG devices of all the buses in the bus set h are disconnected from the power grid, as shown in the following formulas:

[00098]$\begin{array}{cc}{Z}_{mf}^{se}=g_{mf3}\left(i,j,p,h\right)& \text{­­­(33)}\end{array}$

and

[00099]$\begin{array}{cc}{Z}_{ff}^{se}=g_{ff3}\left(i,j,p,h\right)& \text{­­­(34)}\end{array}$

[0145] 5. The specific event at each stage of the multi-stage voltage sag is inferred.

[0146] After the system impedance matrices under different situations are obtained by solving equations in step 4, an optimization model in step 5 may be used to infer a specific faulty line from the faulty line set obtained in step 3, and the specific event of each stage of the multi-stage voltage sag is inferred. There are also the following four inference models based on the sag types:

Type 1

[0147] For this type of sag, the following optimization model may be used to infer the faulty line and its short circuit condition (in other words, i, j, and p ):

[00100]$\begin{array}{cc}\left\{\begin{array}{l}\max \text{-}\sqrt{\left|f\left(g_{mf1}\left(i,j,p\right),g_{ff1}\left(i,j,p\right)\right)\right|-u_1}\\ \text{s}\text{.t}\text{.}i,j\in \left[1,n\right];i\ne j\end{array}\right)& \text{­­­(35)}\end{array}$

[0148] Where f ( ) represents a short-circuit calculation function, and calculation formulas of different types of short circuit faults are shown in the formulas (2) to (6).

Type 2

[0149] For this type of sag, the following optimization model may be used to infer the faulty line and its short circuit condition, and a time point and a sequence of disconnecting the DG device from the power grid (in other words, i, j, p, and h.sub.q).sub.:

[00101]$\begin{array}{cc}\left\{\begin{array}{l}\max \text{-}\left(\begin{array}{c}\sqrt{\left|f\left(g_{mf1}\left(i,j,p\right),g_{ff1}\left(i,j,p\right)\right)\right|-u_1}+\\ {\displaystyle \underset{x=2}{\overset{s}{.Math.}}{\displaystyle \underset{q=1}{\overset{s-1}{.Math.}}\sqrt{\left|f\left(g_{mf3}\left(i,j,p,h_q\right),g_{ff3}\left(i,j,p,h_q\right)\right)\right|-u_x}}}\end{array}\right)\\ \text{s}\text{.t}\text{.}i,j\in \left[1,n\right];i\ne j;h_1.Math.h_2.Math..Math..Math.h_{s-1}\end{array}\right)& \text{­­­(36)}\end{array}$

Type 3

[0150] For this type of sag, the following optimization model may be used to infer the faulty line and its short circuit condition, and a tripping time point and sequence of the relay protection device (in other words, i, j, p, and d ):

[00102]$\begin{array}{cc}\left\{\begin{array}{l}\max \text{-}\left(\begin{array}{c}\sqrt{\left|f\left(g_{mf1}\left(i,j,p\right),g_{ff1}\left(i,j,p\right)\right)\right|-u_1}+\\ \sqrt{\left|f\left(g_{mf2}\left(i,j,p,d\right),g_{ff2}\left(i,j,p,d\right)\right)\right|-u_2}\end{array}\right)\\ \text{s}\text{.t}\text{.}i,j\in \left[1,n\right];i\ne j\end{array}\right)& \text{­­­(37)}\end{array}$

Type 4

[0151] It is assumed that an o .sup.th sudden voltage sag change point detected in a waveform segment of the voltage sag corresponds to a single-side protection trip. For this type of sag, the following optimization model may be used to infer the faulty line and its short circuit condition, the tripping time point and sequence of the relay protection device, and the time point and the sequence of the DG device disconnected from the power grid (in other words, i, j, p, d, and .sup.hq ):

[00103]$\begin{array}{cc}\left\{\begin{array}{l}\max \text{-}\left(\begin{array}{c}\sqrt{\left|f\left(g_{mf1}\left(i,j,p\right),g_{ff1}\left(i,j,p\right)\right)\right|-u_1}+\\ {\displaystyle \underset{x=2}{\overset{o}{.Math.}}{\displaystyle \underset{q=1}{\overset{o-1}{.Math.}}\sqrt{\left|f\left(g_{mf3}\left(i,j,p,h_q\right),g_{ff3}\left(i,j,p,h_q\right)\right)\right|-u_x}}}\\ +\left|f\left(g_{mf2}\left(i,j,p,d\right),g_{ff2}\left(i,j,p,d\right)\right)\right|-u_{o+1}+\\ {\displaystyle \underset{x=o+2}{\overset{s}{.Math.}}{\displaystyle \underset{h=0}{\overset{s-2}{.Math.}}\sqrt{\left|f\left(g_{mf3}\left(i,j,p,h_q\right),g_{ff3}\left(i,j,p,h_q\right)\right)\right|-u_x}}}\end{array}\right)\\ \text{s}\text{.t}\text{.}i,j\in \left[1,n\right];i\ne j;h_1.Math.h_2.Math..Math..Math.h_{s-2}\end{array}\right)& \text{­­­(38)}\end{array}$

[0152] It should be noted that the formulas (35) to (38) are optimization models for a single monitoring device, which can be obtained through solving by using an intelligent algorithm such as a quantum genetic algorithm. When there are a plurality of monitoring devices, an objective function in the optimization model is the sum of respective objectives functions of the monitoring devices.

[0153] 6. The state estimation method for the multi-stage voltage sag is inferred.

[0154] In voltage sag state estimation, a sag amplitude and duration of an unmonitored bus must be estimated. Duration of each sag stage can be estimated based on the time length between adjacent sudden voltage sag change points. The sag amplitude is estimated based on an event inference result in step 5. Similarly, there may be the following four voltage sag state estimation methods based on the voltage sag types:

Type 1

[0155] For this type of sag, a sag amplitude of any unmonitored bus .sup.m may be estimated according to the following formula:

[00104]$\begin{array}{cc}f\left(g_{mf1}\left(i,j,p\right),g_{ff1}\left(i,j,p\right)\right)& \text{­­­(39)}\end{array}$

Type 2

[0156] For this type of sag, a sag amplitude of any unmonitored bus .sup.m may be estimated according to the following formula:

[00105]$\begin{array}{cc}\left\{\begin{array}{l}{u}_1=f\left(g_{mf1}\left(i,j,p\right),g_{ff1}\left(i,j,p\right)\right)\\ u_x=f\left(g_{mf3}\left(i,j,p,h_{x-1}\right),g_{ff3}\left(i,j,p,h_{x-1}\right)\right)\end{array}\right)& \text{­­­(40)}\end{array}$

Type 3

[0157] For this type of sag, a sag amplitude of any unmonitored bus .sup.m may be estimated according to the following formula:

[00106]$\begin{array}{cc}\left\{\begin{array}{l}{u}_1=f\left(g_{mf1}\left(i,j,p\right),g_{ff1}\left(i,j,p\right)\right)\\ u_2=f\left(g_{mf2}\left(i,j,p,d\right),g_{ff2}\left(i,j,p,d\right)\right)\end{array}\right)& \text{­­­(41)}\end{array}$

Type 4

[0158] For this type of sag, a sag amplitude of any unmonitored bus .sup.m may be estimated according to the following formula:

[00107]$\begin{array}{cc}\left\{\begin{array}{l}{u}_1=f\left(g_{mf1}\left(i,j,p\right),g_{ff1}\left(i,j,p\right)\right)\\ u_x=f\left(g_{mf3}\left(i,j,p,h_{x-1}\right),g_{ff3}\left(i,j,p,h_{x-1}\right)\right)\\ u_y=f\left(g_{mf2}\left(i,j,p,d\right),g_{ff2}\left(i,j,p,d\right)\right)\\ u_z=f\left(g_{mf3}\left(i,j,p,h_{z-1}\right),g_{ff3}\left(i,j,p,h_{z-1}\right)\right)\\ x\in \left[2,y-1\right];z\in \left[y+1,s-2\right]\end{array}\right)& \text{­­­(42)}\end{array}$

[0159] The method in the present disclosure analyzes characteristics of multi-stage voltage sags due to different causes, considers sudden change time of each stage and fault clearing time of a relay protection device in the monitoring data of the multi-stage voltage sag, infers a specific event that may occur in each stage of the multi-stage voltage sag, and finally obtains a multi-stage state estimation method that meets project requirements.