A method for disturbance source identification driven by synchronous measurement data in distribution network

Abstract

The invention discloses a disturbance source identification method driven by synchronous measurement data of a distribution network, belongs to the technical field of synchronous measurement, and comprises the following steps: step 1, obtaining real-time disturbance data in the synchronous measurement data of the distribution network and performing data preprocessing; step 2, performing autoregressive modeling on the preprocessed synchronous measurement disturbance data; step 3, calculating the causal strength between the synchronous measurement data of each point and constructing a causal correlation matrix; step 4, performing conditional judgment on the causal correlation matrix according to the topology of the measurement points, and identifying the disturbance source; the method for identifying a disturbance source driven by synchronous measurement data of a distribution network provided by the invention can quickly distinguish between transmission side disturbance and distribution network disturbance in disturbance data, thereby extracting the real distribution network disturbance.

Description

Technical Field

[0001] The invention belongs to the technical field of synchronous measurement, and in particular relates to a disturbance source identification method driven by synchronous measurement data of a distribution network. Background Art

[0002] With the increase in system complexity and the access of a large number of distributed renewable energy sources, the safety and stability of distribution networks have received increasing attention. At the same time, the introduction of synchronous measurement technology on the distribution network side has also made the collection of high-granularity measurement data gradually feasible. With the help of synchronous measurement data obtained by devices such as distribution network PMU, real-time and accurate identification of distribution network disturbances will help the distribution network's event perception, equipment health diagnosis, state estimation, disturbance classification, and positioning, and establish a database foundation for data-driven research.

[0003] However, since the measured data of the distribution network is significantly affected by the transmission side, the identified disturbance data is often mixed with a large number of transmission side disturbances, which brings challenges to the identification of real distribution network disturbances. Specifically, when the real-time detection algorithm detects a disturbance event, it is often necessary to confirm with the site to recognize the source of the disturbance. This operating procedure is only applicable to subsequent analysis and lacks room for online decision-making. State perception requires that the disturbance occurring inside the feeder can be quickly identified after the disturbance is detected. Therefore, there is an urgent need for a disturbance source identification method driven by synchronous measurement data of the distribution network to quickly distinguish the transmission side disturbance and the distribution network disturbance in the disturbance data, so as to extract the real distribution network disturbance. Summary of the invention

[0004] The purpose of the present invention is to provide a disturbance source identification method driven by distribution network synchronous measurement data, so as to solve the problem in the prior art that it is difficult to quickly distinguish between transmission side disturbance and distribution network disturbance in disturbance data.

[0005] To achieve the above object, the present invention provides a method for identifying disturbance sources driven by synchronous measurement data of a distribution network, comprising the following steps:

[0006] Step 1: Obtain real-time disturbance data from the synchronous measurement data of the distribution network and perform data preprocessing;

[0007] Step 2, performing autoregressive modeling on the preprocessed synchronous measurement disturbance data;

[0008] Step 3, calculate the causal strength between the synchronous measurement data of each point and construct a causal correlation matrix;

[0009] Step 4: According to the topology of the measurement points, the causal correlation matrix is ​​conditionally determined to identify the disturbance source.

[0010] Preferably, the specific process of obtaining the real-time disturbance data in the distribution network synchronous measurement data and performing data preprocessing in step 1 is as follows:

[0011] Obtain the real-time disturbance data y from the synchronous measurement data of the distribution network, apply the data preprocessing method based on the BK filter to the real-time disturbance data y, and generate a new disturbance time series Specifically expressed as:

[0012]

[0013] In the formula, t represents time, K represents lag time, and y t-k represents the time series value at lag K, that is, the time series value at time tk, a k Represents the hysteresis coefficient, and the calculation expression of the hysteresis coefficient is as follows:

[0014]

[0015] Among them, j=0,±1,±2,...,±K represents the index range of the time series, a0 represents the constant term, θ is the normalization constant that makes the sum of weights zero, ω1 and ω2 are the low cutoff frequency and the high cutoff frequency, respectively.

[0016] Preferably, the autoregressive modeling of the preprocessed synchronous measurement disturbance data in step 2 is specifically as follows: the synchronous measurement data of the point number i is represented by a differential integrated moving average autoregressive ARIMA (p, d, q) model, where p is the autoregressive order, q is the moving average order, and d is the order of difference required to make the autoregressive order a stationary sequence; wherein, for the new disturbance time series y t * Perform a d-order difference step to stabilize the time series. The specific expression is as follows:

[0017]

[0018] In the formula, represents the d-order difference operator, which performs d-time difference on the time series. Represents the time series value after the first-order difference, y t * 'Indicates calculation Abbreviation of y * t-1 represents the value at the moment before t;

[0019] The new perturbation time series The regression formula is expressed as:

[0020]

[0021] In the formula, θi is the sliding average coefficient, r i is the autoregressive coefficient, ε t is the residual value, represents the value at time i before t, μ represents the zero-mean Gaussian white noise process;

[0022] Define the regression error ε in the resulting ARIMA model t is the autoresidual error of the measurement point i, and the autoregression coefficient r is used to evaluate the regression effect. i The calculation method is as follows:

[0023]

[0024] In the formula, r i It is expressed as the cumulative sum of the residual values ​​at all time points when the time series at the i-th position is used as an exogenous variable in the ARIMAX model.

[0025] Preferably, step 3 specifically comprises: calculating the causal correlation matrix based on the quantitative verification method of the causal relationship between the predicted time series by using the synchronous measurement data of all measurement points installed in the same substation (feeder) and the synchronous measurement data of all measurement points installed in another substation (feeder) adjacent thereto, including the following steps:

[0026] S31. Assume that the time series vectors x and y represent the synchronous measurement data of the points with sequence numbers i and j respectively, and construct the corresponding regression formula expression, which is as follows:

[0027]

[0028] In the formula, y tx represents a regression that uses only the past values ​​of y to predict the current value of y, y ty represents a regression that uses the past values ​​of y and x to predict the current value of y. j With b j represents the coefficient in the model, and d represents the number of differences;

[0029] S32, judge y tx With y ty The prediction accuracy of tx y ty If the predicted value is accurate, it is determined that the time series x has a Granger causal relationship with y. Based on the Granger causal relationship and the regression model in step 2, a new ARIMA regression model is constructed as follows:

[0030]

[0031] Among them, m jrepresents the jth correlation coefficient in the model, x t-j Represents the value of the time series x at time tj, and defines the regression error ε in this formula t ′ is the mutual residual between measurement points i and j, and the parameter r used to evaluate the regression effect of the new ARIMA model is i,j It is expressed as:

[0032] S33, according to r i,j and r i And combined with the information theory formula, the Granger causality strength from measurement point j to measurement point i is obtained as follows:

[0033]

[0034] Among them, j→i represents the direction of Granger causality, that is, j causes i; g j→i The larger the value, the more accurate the regression considering point j is than the regression considering only point i, and the stronger the Granger causality from j to i is.

[0035] S34, using elements g between different points j→i The matrix G composed of is the causal correlation matrix between the measurements at different locations.

[0036] Preferably, in step 4, the causal correlation matrix is ​​conditionally determined according to the topology of the measurement points, and the disturbance source is specifically identified as follows: assuming that there are M points in a substation, using the diagonal and off-diagonal submatrices G m,m With G m,n , m,n=1,…,M,m≠n, respectively representing the Granger causality strength within the feeder and the Granger causality strength between feeders. If the causal correlation matrix satisfies the following conditions, the disturbance is determined to come from the transmission side; otherwise, it is considered that the disturbance comes from the local distribution network or the frequency estimation algorithm misjudges the disturbance. The specific expression is as follows:

[0037]

[0038] Therefore, the present invention adopts the above-mentioned distribution network synchronous measurement data-driven disturbance source identification method, which can identify disturbances in the synchronous measurement data, distinguish the transmission side disturbances and distribution network disturbances in the data, thereby identifying the real distribution network disturbances, and providing a basis for data-driven distribution network disturbance analysis and intervention.

[0039] The technical solution of the present invention is further described in detail below through the accompanying drawings and embodiments. BRIEF DESCRIPTION OF THE DRAWINGS

[0040] Figure 1It is an overall flow chart of a disturbance source identification method driven by synchronous measurement data of a distribution network of the present invention;

[0041] Figure 2 A synchronous measurement data diagram of a distribution network measured when a transmission side disturbance occurs provided by an embodiment of the present invention;

[0042] Figure 3 A synchronous measurement data diagram of a distribution network measured when no transmission side disturbance occurs provided by an embodiment of the present invention;

[0043] Figure 4 A schematic diagram of the calculation results of the causal correlation matrix when a transmission side disturbance occurs provided by an embodiment of the present invention;

[0044] Figure 5 A schematic diagram of the calculation results of the causal correlation matrix when no transmission side disturbance occurs provided by an embodiment of the present invention;

[0045] Figure 6 A schematic diagram of the causality comparison between feeders and within a feeder provided in an embodiment of the present invention. DETAILED DESCRIPTION

[0046] The following detailed description of the embodiments of the present invention provided in the accompanying drawings is not intended to limit the scope of the invention claimed for protection, but merely represents selected embodiments of the present invention. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without creative work are within the scope of protection of the present invention.

[0047] See also Figure 1-6 , a method for identifying disturbance sources driven by synchronous measurement data of a distribution network, comprising the following steps:

[0048] Step 1: Obtain real-time disturbance data from the synchronous measurement data of the distribution network and perform data preprocessing; the specific process is as follows:

[0049] Obtain the real-time disturbance data y from the synchronous measurement data of the distribution network, apply the data preprocessing method based on the BK filter to the real-time disturbance data y, and generate a new disturbance time series Specifically expressed as:

[0050]

[0051] In the formula, t represents time, K represents lag time, and y t-k represents the time series value at lag K, that is, the time series value at time tk, a k Represents the hysteresis coefficient, and the calculation expression of the hysteresis coefficient is as follows:

[0052]

[0053] Among them, j = 0, ±1, ±2, …, ±K represents the index range of the time series, a0 represents the constant term, θ is the normalization constant that makes the sum of weights zero, ω1 and ω2 are the low cutoff frequency and the high cutoff frequency, respectively.

[0054] Step 2: Autoregressive modeling is performed on the preprocessed synchronous measurement disturbance data. Specifically, the synchronous measurement data of point i is represented by the differential integrated moving average autoregressive ARIMA (p, d, q) model, where p is the autoregressive order, q is the moving average order, and d is the order of difference required to make the autoregressive order a stationary series; Perform a d-order difference step to stabilize the time series. The specific expression is as follows:

[0055]

[0056] In the formula, represents the d-order difference operator, which performs d-time difference on the time series. Represents the time series value after the first-order difference, y t * 'Indicates calculation Abbreviation of y * t-1 represents the value at the moment before t;

[0057] The new perturbation time series The regression formula is expressed as:

[0058]

[0059] In the formula, θ i is the sliding average coefficient, r i is the autoregressive coefficient, ε t is the residual value, y * t-i represents the value at time i before t, μ represents the zero-mean Gaussian white noise process, and if the regression process does not add additional noise signals, μ is 0.

[0060] Define the regression error ε in the resulting ARIMA model t is the autoresidual error of the measurement point i, and the autoregression coefficient r is used to evaluate the regression effect. i The calculation method is as follows:

[0061]

[0062] In the formula, r i It is expressed as the cumulative sum of the residual values ​​at all time points when the time series at the i-th position is used as an exogenous variable in the ARIMAX model.

[0063] Step 3, calculate the causal strength between the synchronous measurement data of each point and construct a causal correlation matrix; specifically, the synchronous measurement data of all measurement points installed in the same substation (feeder) and the synchronous measurement data of all measurement points installed in another adjacent substation (feeder) are used to jointly calculate the causal correlation matrix based on the quantitative test method of the causal relationship between the predicted time series, including the following steps:

[0064] S31. Assume that the time series vectors x and y represent the synchronous measurement data of the points with sequence numbers i and j respectively, and construct the corresponding regression formula expression, which is as follows:

[0065]

[0066] In the formula, y tx represents a regression that uses only the past values ​​of y to predict the current value of y, y ty represents a regression that uses the past values ​​of y and x to predict the current value of y. j With b j represents the coefficient in the model, and d represents the number of differences;

[0067] S32, judge y tx With y ty The prediction accuracy of tx y ty If the predicted value is accurate, it is determined that the time series x has a Granger causal relationship with y. Based on the Granger causal relationship and the regression model in step 2, a new ARIMA regression model is constructed as follows:

[0068]

[0069] Among them, m j represents the jth correlation coefficient in the model, x t-j Represents the value of the time series x at time tj, and defines the regression error ε in this formula t ′ is the mutual residual between measurement points i and j, and the parameter r used to evaluate the regression effect of the new ARIMA model is i,j It is expressed as:

[0070] S33, according to r i,j and r i And combined with the information theory formula, the Granger causality strength from measurement point j to measurement point i is obtained as follows:

[0071]

[0072] Among them, j→i represents the direction of Granger causality, that is, j causes i; gj→i The larger the value, the more accurate the regression considering point j is than the regression considering only point i, and the stronger the Granger causality from j to i is.

[0073] S34, using elements g between different points j→i The matrix G composed of is the causal correlation matrix between the measurements at different locations.

[0074] Step 4: According to the topology of the measurement points, the causal correlation matrix is ​​conditionally determined to identify the disturbance source; specifically: suppose there are M points in a substation, and use the diagonal and off-diagonal submatrices G m,m With G m,n , m,n=1,…,M,m≠n, respectively representing the Granger causality strength within the feeder and the Granger causality strength between feeders. If the causal correlation matrix satisfies the following conditions, the disturbance is determined to come from the transmission side; otherwise, it is considered that the disturbance comes from the local distribution network or the frequency estimation algorithm misjudges the disturbance. The specific expression is as follows:

[0075]

[0076] Example

[0077] A continuous trip occurred in the Jibei power grid. The transmission side disturbance caused the synchronous measurement data in the two feeders A and B under the same substation to trigger the disturbance detection alarm. A1, A2, B1, and B2 represent the four distribution network PMUs on the two feeders under the substation. The synchronous measurement data when the disturbance occurred is as follows: Figure 2 In order to compare with the situation when the disturbance occurs, the normal data without disturbance is also analyzed based on Granger causality strength. The synchronous measurement data without transmission side disturbance is shown as follows: Figure 3 In this case, measurements from the same feeder, such as A1 and A2, are considered intra-feeder measurements, and measurements from different feeders connected to the same substation, such as A1 and B1, are considered inter-feeder measurements.

[0078] After calculation using the embodiment of the invention, when a disturbance occurs on the transmission side, the intra-feeder and inter-feeder Granger causality matrices composed of A1, A2 and B1, B2 are summarized as follows: Figure 4 When there is no disturbance on the transmission side, the intra-feeder and inter-feeder Granger causality matrices composed of A1, A2 and B1, B2 are as follows: Figure 5 shown. Figure 4In the above figure, the Granger causality between two PMUs on the same feeder is considered as intra-feeder causality, and the matrix is ​​marked with a dotted line. The Granger causality between different feeders is considered as inter-feeder causality, and the matrix is ​​marked with a solid line. The average value of the intra-feeder causality strength is used as a benchmark. When the average inter-feeder causality strength is higher than the benchmark, it is identified as a transmission side disturbance. Figure 4 and Figure 5 The numbers in the dotted lines are the summed mean values ​​of the intra-feeder causal calculation results (the upper left corner is the intra-feeder causal mean value of feeder A, and the lower right corner is the intra-feeder causal mean value of feeder B); the numbers in the solid lines are the summed mean values ​​of the inter-feeder causal calculation results (the upper right corner is the inter-feeder causal mean value of feeder A to feeder B, and the lower left corner is the inter-feeder causal mean value of feeder B to feeder A).

[0079] The comparison of the Granger causality strength with and without transmission side disturbance is as follows: Figure 6 As shown in the figure, it can be seen that in the data with transmission side disturbance, the average intra-feeder Granger causality strength is smaller than the inter-feeder causality strength, while the opposite is true when there is no transmission side disturbance.

[0080] The above-mentioned field cases show that there are significant differences in the calculated values ​​of the causal strength matrix inside and outside the feeder when a disturbance occurs on the transmission side and when no disturbance occurs, which verifies the effectiveness and feasibility of analyzing and tracing the source of field disturbance data based on the Granger causal strength calculation method in the embodiment of the invention.

[0081] Therefore, the present invention adopts the above-mentioned distribution network synchronous measurement data-driven disturbance source identification method, which can identify disturbances in the synchronous measurement data, distinguish the transmission side disturbances and distribution network disturbances in the data, thereby identifying the real distribution network disturbances, and providing a basis for data-driven distribution network disturbance analysis and intervention.

[0082] Finally, it should be noted that the above embodiments are only used to illustrate the technical solution of the present invention rather than to limit it. Although the present invention has been described in detail with reference to the preferred embodiments, those skilled in the art should understand that they can still modify or replace the technical solution of the present invention with equivalents, and these modifications or equivalent replacements cannot cause the modified technical solution to deviate from the spirit and scope of the technical solution of the present invention.

Claims

1. A method for identifying disturbance sources driven by synchronous measurement data of a distribution network, characterized in that: The following steps are involved: Step 1: Obtain real-time disturbance data from the synchronous measurement data of the distribution network and perform data preprocessing; Step 2, performing autoregressive modeling on the preprocessed synchronous measurement disturbance data; Step 3, calculate the causal strength between the synchronous measurement data of each point and construct a causal correlation matrix; Step 4: Condition the causal correlation matrix according to the topology of the measurement points to identify the disturbance source; Step 3 specifically includes: calculating the causal correlation matrix based on the quantitative test method of the causal relationship between the predicted time series of the synchronous measurement data of all measurement points installed on the same feeder and the synchronous measurement data of all measurement points installed on another adjacent feeder, including the following steps: S31. Assume that the time series vectors x and y represent the synchronous measurement data of the points with sequence numbers i and j respectively, and construct the corresponding regression formula expression, which is as follows: In the formula, y tx represents a regression that uses only the past values ​​of y to predict the current value of y, y ty represents a regression that uses the past values ​​of y and x to predict the current value of y. j With b j represents the coefficient in the model, and d represents the number of differences; S32, judge y tx With y ty The prediction accuracy of tx y ty If the predicted value is accurate, it is determined that the time series x has a Granger causal relationship with y. Based on the Granger causal relationship and the regression model in step 2, a new ARIMA regression model is constructed as follows: Among them, m j represents the jth correlation coefficient in the model, x t-j Represents the value of the time series x at time tj, and defines the regression error ε in this formula t ′ is the mutual residual between measurement points i and j, and the parameter r used to evaluate the regression effect of the new ARIMA model is i,j It is expressed as: S33, according to r i,j and r i And combined with the information theory formula, the Granger causality strength from measurement point j to measurement point i is obtained as follows: Among them, j→i represents the direction of Granger causality, that is, j causes i; g j→i The larger the value, the more accurate the regression considering point j is than the regression considering only point i, and the stronger the Granger causality from j to i is. S34, using elements g between different points j→i The composed matrix G obtains the causal correlation matrix between the measurement values ​​at different points; In step 4, the causal correlation matrix is ​​conditionally determined based on the topology of the measurement points. The specific identification of the disturbance source is as follows: Assume that there are M points in a substation, and use the diagonal and off-diagonal submatrices G m,m With G m,n , m,n=1,…,M,m≠n, respectively representing the Granger causality strength within the feeder and the Granger causality strength between feeders. If the causal correlation matrix satisfies the following conditions, the disturbance is determined to come from the transmission side; otherwise, it is considered that the disturbance comes from the local distribution network or the frequency estimation algorithm misjudges the disturbance. The specific expression is as follows:

2. The method for identifying disturbance sources driven by synchronous measurement data of a distribution network according to claim 1, characterized in that: The specific process of obtaining the real-time disturbance data in the distribution network synchronous measurement data and performing data preprocessing in step 1 is as follows: Obtain the real-time disturbance data y from the synchronous measurement data of the distribution network, apply the data preprocessing method based on the BK filter to the real-time disturbance data y, and generate a new disturbance time series Specifically expressed as: In the formula, t represents time, K represents lag time, and y t-k represents the time series value at lag K, that is, the time series value at time tk, a k Represents the hysteresis coefficient, and the calculation expression of the hysteresis coefficient is as follows: Among them, j=0,±1,±2,...,±K represents the index range of the time series, a0 represents the constant term, θ is the normalization constant that makes the sum of weights zero, ω1 and ω2 are the low cutoff frequency and the high cutoff frequency, respectively.

3. The method for identifying disturbance sources driven by synchronous measurement data of a distribution network according to claim 2, characterized in that: In step 2, the autoregressive modeling of the preprocessed synchronous measurement disturbance data is specifically as follows: the synchronous measurement data of the point number i is represented by the differential integrated moving average autoregressive ARIMA (p, d, q) model, where p is the autoregressive order, q is the moving average order, and d is the order of difference required to make the autoregressive order a stationary series; Perform a d-order difference step to stabilize the time series. The specific expression is as follows: and t * =▽ d and t * ,▽and t * =and t * '=and t * -and * t-1 In the formula, ▽ d represents the d-order difference operator, which performs d-time difference on the time series. Represents the time series value after the first-order difference, y t * 'Indicates calculation Abbreviation of y * t-1 represents the value at the moment before t; The new perturbation time series The regression formula is expressed as: In the formula, θ i is the sliding average coefficient, r i is the autoregressive coefficient, ε t is the residual value, y * t-i represents the value at time i before t, μ represents the zero-mean Gaussian white noise process; Define the regression error ε in the resulting ARIMA model t is the autoresidual error of the measurement point i, and the autoregression coefficient r is used to evaluate the regression effect. i The calculation method is as follows: In the formula, r i It is expressed as the cumulative sum of the residual values ​​at all time points when the time series at the i-th position is used as an exogenous variable in the ARIMAX model.

Images