Forced power oscillation disturbance source positioning method based on KDE transfer entropy causal analysis

Abstract

The invention discloses a forced power oscillation disturbance source positioning method based on KDE transfer entropy causal analysis, and belongs to the technical field of power system operation and maintenance. The forced power disturbance source positioning is realized by analyzing the causal relationship between disturbance source and non-disturbance source unit measurement signals in the forced power oscillation process of the system. According to the method, the causal analysis method is used for solving the forced power oscillation disturbance source positioning problem, dependence on a system model and parameters is eliminated, rapid online positioning of the forced oscillation disturbance source can be achieved only according to measurement data, and the positioning precision and positioning efficiency of the forced oscillation disturbance source are effectively improved.

Description

Technical Field

[0001] The invention belongs to the technical field of power system operation and maintenance, and in particular relates to a method for locating a forced power oscillation disturbance source based on KDE transfer entropy causal analysis. Background Art

[0002] Forced oscillation of power system is an oscillation phenomenon caused by continuous periodic disturbance. When the external disturbance frequency is close to the natural oscillation frequency of the system, it may cause resonance, seriously affecting the safety and stability of the power system. The characteristics of forced power oscillation include rapid onset, strong randomness, and rapid calming after the disturbance source is lost. Therefore, in the context of the expansion of the power grid and the widespread application of new energy, quickly and accurately locating the oscillation source is the key to suppressing forced oscillation.

[0003] There are two main types of forced power oscillation source positioning methods that are currently widely recognized: positioning methods based on physical mechanisms and positioning methods based on data-driven. The accuracy of methods based on physical mechanisms depends on the detailed models and parameters of the system, but in actual complex systems, it is often difficult to obtain these models and parameters, which limits the application of the methods. Data-driven methods have high requirements on the quality of input data. If there is noise, missing or outliers in the data, it will significantly affect the accuracy of disturbance source positioning.

[0004] In recent years, the scale of power system measurement data has continued to expand, making causal analysis show great potential in revealing the causal relationship between variables. This analysis method can explore the potential driving factors in time series, reveal the mechanism of event occurrence, and provide guidance for intervention measures. Using causal analysis methods to solve the problem of locating the source of forced power oscillation disturbances will get rid of the dependence on system models and parameters. However, at this stage, there are few studies on the application of causal analysis methods to the location of forced power oscillation disturbance sources.

[0005] Therefore, a new technical solution is urgently needed in the prior art to solve this problem. Summary of the invention

[0006] The technical problem to be solved by the present invention is to provide a forced power oscillation disturbance source locating method based on KDE transfer entropy causal analysis to solve the technical problem that causal analysis methods are rarely applied to forced power oscillation disturbance source locating research at this stage.

[0007] The method for locating the source of forced power oscillation disturbance based on KDE transfer entropy causal analysis includes the following steps, and the following steps are performed in sequence:

[0008] Step 1: Obtain the generator G based on time domain simulation or vector measurement unit PMU measurement A , G BActive power response signal, low-pass filtering and band-pass filtering pre-process the active power input signal, filter out the frequency range of no concern, and normalize the oscillation data to eliminate the influence of the magnitude difference of the measured data on the accuracy of the positioning result; the processed data is used to establish the generator G A , G B Active power time series: and Among them, P GAn Represents the generator G A The active power at the nth sampling point, P GBn Indicates the generator G B Active power at the nth sampling point, time series P GA Satisfies the k-order Markov process, the time series P GB If the i-th order Markov process is satisfied, the active power time series P is defined as GA With P GB Entropy is transferred between them;

[0009] Step 2: Based on Bayesian theorem, the active power time series P GA With P GB The transfer entropy is converted into joint probability density;

[0010] Step 3: Establish a multivariate kernel density estimation model, and select the Gaussian function as the multivariate kernel function in the model to improve the accuracy and applicability of transfer entropy calculation;

[0011] Step 4: Use the multivariate kernel density estimation model established in step 3 to improve the conditional probability of the joint probability density in step 2 to obtain the improved conditional probability; according to the active power time series P in step 1 GA With P GB The entropy transfer between them and the improved conditional probability calculation of the kernel probability density KDE transfer entropy of two active power time series;

[0012] Step 5: To avoid the influence of the randomness of the single causal analysis calculation results on the accuracy of the disturbance source location results, a sliding window analysis method is used to repeat the process from step 1 to step 4 to calculate and obtain multiple groups of generators G A , G B The transfer entropy data group between active powers is calculated, and the mean value of the transfer entropy data group is obtained. According to the size of the mean value, the generator G is realized A , G B Accurate quantification of causal relationships between active powers;

[0013] Step 6: By comparing the generator G A , G B The size of KDE transfer entropy of active power time series Identify the generator G A , G B The causal relationship between active power and generator G A , G B The direction of active power transmission is as follows:

[0014] 1)P GA P GB The transfer entropy is greater than P GB P GA The transfer entropy is Then put P GA is called the dependent variable, P GB It is called the result variable, indicating that the power transfer direction in this forced oscillation mode is P GA →P GB , then the forced power oscillation disturbance source is determined to be the generator G A ;

[0015] 2) P GA P GB The transfer entropy is less than P GB P GA The transfer entropy is Then put P GB is called the dependent variable, P GA It is called the result variable, indicating that the power transfer direction in this forced oscillation mode is P GB →P GA , then the forced power oscillation disturbance source is determined to be the generator G B ;

[0016] 3) P GA P GB The transfer entropy and P GB P GA The difference in the transfer entropy of is within the set approximate relationship threshold, then Determine P GA With P GB There is no obvious causal relationship between them.

[0017] In step 1, the active power time series P GA With P GB The entropy transferred between them is:

[0018]

[0019] Where: express and The joint probability density of express and The joint probability density of and represents conditional probability;

[0020] in, Indicates the generator G A The active power at the n+1th sampling point, Indicates the generator G B The active power at the n+1th sampling point, Represents a time series P that satisfies a k-order Markov process GA , It represents the time series P that satisfies the i-th order Markov process GB ; P GB Take the (n+1)th sampling point, and P GA Take the conditional probability of the nth sampling point; P GA Take the (n+1)th sampling point, and P GB Take the conditional probability of the nth sampling point.

[0021] The joint probability density in step 2 is shown in formulas (2) and (3):

[0022]

[0023]

[0024] The multivariate kernel density estimation model in step 3 is shown in formulas (4) and (5):

[0025]

[0026] Where: H is the coefficient matrix, which is an i×i-dimensional symmetric positive definite matrix, and its coefficients are selected to meet the requirement of maximizing the accuracy of the joint probability density; K(.) is the multivariate kernel function;

[0027] K(.) is represented by the function K(x) and meets the following requirements to minimize its impact on the accuracy of probability density estimation:

[0028]

[0029] Selecting Gaussian function as the multivariate kernel function K(x) further improves the accuracy of transfer entropy calculation:

[0030]

[0031] In formula (6) and formula (7), R i represents the set of real numbers, I i Represents the set of imaginary numbers, x i represents the time series, x represents the time series x i The sampling value of θ represents the time series xi The variance of

[0032] The improved conditional probability in step 4 is expressed as:

[0033]

[0034] Where: H 1 , H 2 and H 3 are all symmetric positive definite coefficient matrices; K(.) is a multivariate kernel function, where the Gaussian function is selected; N is the corresponding discrete number in the univariate kernel probability density estimation model, and m represents the generator G B The mth item of the active power time series, d represents the generator G A The dth term of the active power time series.

[0035] Through the above design scheme, the present invention can bring the following beneficial effects:

[0036] The present invention realizes the positioning of the forced power disturbance source by analyzing the causal relationship between the disturbance source and the non-disturbance source unit measurement signal during the system forced power oscillation process. The present invention uses the causal analysis method to solve the problem of locating the forced power oscillation disturbance source, gets rid of the dependence on the system model and parameters, and can realize the rapid online positioning of the forced oscillation disturbance source based only on the measurement data, effectively improving the positioning accuracy and positioning efficiency of the forced oscillation disturbance source. BRIEF DESCRIPTION OF THE DRAWINGS

[0037] The present invention is further described below with reference to the accompanying drawings and specific embodiments:

[0038] Figure 1 This is a flow chart of the method for locating the source of forced power oscillation disturbance based on KDE transfer entropy causal analysis of the present invention.

[0039] Figure 2 This is a wiring diagram of an IEEE 4-machine 2-area structure in an embodiment of a method for locating a forced power oscillation disturbance source based on KDE transfer entropy causal analysis of the present invention.

[0040] Figure 3 It is a relationship diagram of the directivity of the active power transfer entropy of the generator in the embodiment of the method for locating the source of forced power oscillation disturbance based on KDE transfer entropy causal analysis of the present invention.

[0041] Figure 4 It is the calculation result of the change of the transfer entropy difference over time in the embodiment of the forced power oscillation disturbance source positioning method based on KDE transfer entropy causal analysis of the present invention. DETAILED DESCRIPTION

[0042] A method for locating the source of forced power oscillation disturbance based on KDE transfer entropy causal analysis, such as Figure 1 As shown, the following steps are included and performed sequentially.

[0043] Step 1: Obtain the generator G based on time domain simulation or vector measurement unit PMU measurement A , G B Active power response signal, low-pass filtering and band-pass filtering pre-process the active power input signal, filter out the frequency range of no concern, and normalize the oscillation data to eliminate the influence of the magnitude difference of the measured data on the accuracy of the positioning result; the processed data is used to establish the generator G A , G B Active power time series: and Among them, P GAn Indicates P GA The active power at the nth sampling point, P GBn Indicates P GB Active power at the nth sampling point, time series P GA Satisfies the k-order Markov process, the time series P GB If the i-th order Markov process is satisfied, the active power time series P is defined as GA With P GB Entropy is transferred between them;

[0044] Step 2: Based on Bayesian theorem, the active power time series P GA With P GB The transfer entropy is converted into joint probability density;

[0045] Step 3: Establish a multivariate kernel density estimation model, and select the Gaussian function as the multivariate kernel function in the model to improve the accuracy and applicability of transfer entropy calculation;

[0046] Step 4: Use the multivariate kernel density estimation model established in step 3 to improve the conditional probability of the joint probability density in step 2 to obtain the improved conditional probability; according to the active power time series P in step 1 GA With P GB The transfer entropy between them and the improved conditional probability calculation of the KDE (kernel probability density) transfer entropy of two active power time series;

[0047] Step 5: To avoid the influence of the contingency of the single causal analysis calculation result on the accuracy of the disturbance source location result, the present invention adopts a sliding window analysis method, repeats the process from step 1 to step 4, and calculates multiple groups of generators G A , G B The transfer entropy data group between active powers is further calculated and the mean value of the transfer entropy data group is obtained. Thus, the generator G is realized according to the size of the mean value. A , GB Accurate quantification of cause-effect relationships between active powers.

[0048] Step 6: By comparing the generator G A , G B The size of KDE transfer entropy of active power time series Identify the generator G A , G B The causal relationship between active power and generator G A , G B The direction of active power transmission. Specifically:

[0049] 1) If P GA P GB The transfer entropy is greater than P GB P GA The transfer entropy is Just put P GA is called the “dependent” variable, P GB It is called the “result” variable, indicating that the power transfer direction in this forced oscillation mode is P GA (cause) → P GB (Result), then it is determined that the source of the forced power oscillation disturbance is the generator G A .

[0050] 2) If P GA P GB The transfer entropy is less than P GB P GA The transfer entropy is Then put P GB is called the “dependent” variable, P GA It is called the “result” variable; it indicates that the power transfer direction in this forced oscillation mode is P GB (cause) → P GA (Result), then it is determined that the source of the forced power oscillation disturbance is the generator G B .

[0051] 3) If P GA P GB The transfer entropy and P GB P GA The difference in the transfer entropy of is within the set approximate relationship threshold, then Then determine P GA With P GB There is no obvious causal relationship between them.

[0052] In step 1, the active power time series P GA With P GB The entropy transferred between them is:

[0053]

[0054] Where: express and The joint probability density of express and The joint probability density of and represents conditional probability;

[0055] in, Indicates the generator G A The active power at the n+1th sampling point, Indicates the generator G B The active power at the n+1th sampling point, Represents a time series P that satisfies a k-order Markov process GA , It represents the time series P that satisfies the i-th order Markov process GB ; P GB Take the (n+1)th sampling point, and P GA Take the conditional probability of the nth sampling point; P GA Take the (n+1)th sampling point, and P GB Take the conditional probability of the nth sampling point.

[0056] The joint probability density in step 2 is shown in formulas (2) and (3):

[0057]

[0058] The multivariate kernel density estimation model in step 3 is shown in formulas (4) and (5):

[0059]

[0060] Where: H is the coefficient matrix, which is an i×i-dimensional symmetric positive definite matrix, and its coefficients are selected to meet the requirement of maximizing the accuracy of the joint probability density; K(.) is the multivariate kernel function;

[0061] K(.) is represented by the function K(x) and meets the following requirements to minimize its impact on the accuracy of probability density estimation:

[0062]

[0063] Selecting Gaussian function as the multivariate kernel function K(x) further improves the accuracy of transfer entropy calculation:

[0064]

[0065] In formula (6) and formula (7), R i represents the set of real numbers, I i Represents the set of imaginary numbers, x i represents the time series, x represents the time series x i The sampling value of θ represents the time series x i The variance of

[0066] The improved conditional probability in step 4 is expressed as:

[0067]

[0068] Where: H 1 , H 2 and H 3 are all symmetric positive definite coefficient matrices; K(.) is a multivariate kernel function, where the Gaussian function is selected; N is the corresponding discrete number in the univariate kernel probability density estimation model, and m represents the generator G B The mth item of the active power time series, d represents the generator G A The dth term of the active power time series.

[0069] The present invention is further described in detail below in conjunction with embodiments:

[0070] like Figure 2 The typical IEEE 4-machine 2-region system is shown. Under the basic operation state, the small disturbance calculation results show that the system has 3 electromechanical oscillation modes, and the damping characteristics of each oscillation mode are good. The present invention takes the inter-regional oscillation mode as an example to verify the accuracy and effectiveness of the proposed forced power oscillation disturbance source location method.

[0071] In the simulation, a sinusoidal disturbance signal with an amplitude of 0.015 pu and a frequency of 0.6 Hz is continuously injected into the excitation system of generator G1 to resonate with the 0.61 Hz inter-regional mode of the system. At the same time, in order to simulate environmental excitation, all loads in the system fluctuate randomly at 3% of the base value, and the simulation lasts for 200 seconds.

[0072] The KDE transfer entropy method of the present invention is used to calculate the transfer entropy between the active powers of the four generators in the system. The single transfer entropy calculation result (at the 50th second) is shown in Table 1.

[0073] Table 1 Single transfer entropy calculation results

[0074]

[0075] Here, if X is G1 and Y is G2, then it represents the transfer entropy T of the active power of generator G1 to the active power of generator G2. G1→G2; If X is G2 and Y is G1, then the transfer entropy T of the active power of generator G2 to generator G1 is G2→G1 The calculation method of other transfer entropies is similar. The directional relationship of the active power transfer entropy of the four generators is as follows: Figure 3 as shown in .

[0076] From Table 1 and Figure 3 From the analysis, it can be seen that the KDE transfer entropy (T) of the active power of generator G1 to the other three generators (G2, G3, G4) G1→G2 、T G1→G3 、T G1→G4 ) are 1.004, 1.131, and 1.333, respectively, which are all higher than the active power KDE transfer entropy (T G2→G1 、T G3→G1 、T G4→G1 ) are 0.859, 1.010, and 1.037, respectively. According to the causal analysis disturbance source location method of KDE transfer entropy, it can be seen that the active power oscillation of generator G1 is the main cause of the active power oscillation of generators G2, G3, and G4. The location result is consistent with the disturbance source setting, which further verifies the accuracy of the method proposed in the present invention.

[0077] In the actual operation of the system, due to the different noise intensities at different moments, the measured data will have certain differences in different time sections. In order to reduce the impact of the randomness of a single calculation result on the accuracy of the positioning result, the present invention uses a sliding window analysis method to calculate the transfer entropy between the active power of each generator. Similarly, a 50-second analysis window is selected, the active power of the generator is input, and the window is slid once every 2 seconds. The transfer entropy between the active power of each generator is calculated within 200 seconds, and the mean statistical results are listed in Table 2. The difference results of different transfer entropies are shown in Table 2. Figure 4 shown.

[0078] Table 2 System transfer entropy mean statistics

[0079]

[0080] from Figure 4 It can be seen that as the calculation window slides, T G1→G2、G3、G4 With T G2、G3、G4→G1 The difference between them also fluctuates to a certain extent, but the difference in transfer entropy is always greater than 0. This means that T G1→G2、G3、G4 Always higher than T G2、G3、G4→G1 , indicating that the oscillation transmission direction of the active power of the generator is always from the generator G1 to other generators in the system. Therefore, the generator G1 is identified as the disturbance source unit, which further verifies the accuracy of the disturbance source positioning method proposed in the present invention.

Claims

1. A method for locating the source of forced power oscillation disturbance based on KDE transfer entropy causal analysis, characterized by: The process includes the following steps, which are performed in sequence: Step 1: Obtain the generator G based on time domain simulation or vector measurement unit PMU measurement A , G B Active power response signal, low-pass filtering and band-pass filtering pre-process the active power input signal, filter out the frequency range of no concern, and normalize the oscillation data to eliminate the influence of the magnitude difference of the measured data on the accuracy of the positioning result; the processed data is used to establish the generator G A , G B Active power time series: P GA :{P GA1 ,P GA2 ,…,P GAn } and P GB :{P GB1 ,P GB2 ,…,P GBn }, where P GAn Indicates the generator G A The active power at the nth sampling point, P GBn Indicates the generator G B Active power at the nth sampling point, time series P GA Satisfies the k-order Markov process, the time series P GB If the i-th order Markov process is satisfied, the active power time series P is defined as GA With P GB Entropy is transferred between them; Step 2: Based on Bayesian theorem, the active power time series P GA With P GB The transfer entropy is converted into joint probability density; Step 3: Establish a multivariate kernel density estimation model, and select the Gaussian function as the multivariate kernel function in the model to improve the accuracy and applicability of transfer entropy calculation; Step 4: Use the multivariate kernel density estimation model established in step 3 to improve the conditional probability of the joint probability density in step 2 to obtain the improved conditional probability; according to the active power time series P in step 1 GA With P GB The entropy transfer between them and the improved conditional probability calculation of the kernel probability density KDE transfer entropy of two active power time series; Step 5: To avoid the influence of the randomness of the single causal analysis calculation results on the accuracy of the disturbance source location results, a sliding window analysis method is used to repeat the process from step 1 to step 4 to calculate and obtain multiple groups of generators G A , G B The transfer entropy data group between active powers is calculated, and the mean value of the transfer entropy data group is obtained. According to the size of the mean value, the generator G is realized A , G B Accurate quantification of causal relationships between active powers; Step 6: By comparing the generator G A , G B The size of KDE transfer entropy of active power time series (T PGA→PGB , T PGB→PGA ) Identify the generator G A , G B The causal relationship between active power and generator G A , G B The direction of active power transmission is: 1)P GA P GB The transfer entropy is greater than P GB P GA The transfer entropy, T PGA→PGB >T PGB→PGA , then P GA is called the dependent variable, P GB It is called the result variable, indicating that the power transfer direction in this forced oscillation mode is P GA →P GB , then the forced power oscillation disturbance source is determined to be the generator G A ; 2) P GA P GB The transfer entropy is less than P GB P GA The transfer entropy, T PGA→PGB <T PGB→PGA , then P GB is called the dependent variable, P GA It is called the result variable, indicating that the power transfer direction in this forced oscillation mode is P GB →P GA , then the forced power oscillation disturbance source is determined to be the generator G B ; 3) P GA P GB The transfer entropy and P GB P GA If the difference of the transfer entropy is within the set approximate relationship threshold, then TP GA →P GB ≈T PGB→PGA , determine P GA With P GB There is no obvious causal relationship between them.

2. The method for locating a forced power oscillation disturbance source based on KDE transfer entropy causal analysis according to claim 1 is characterized in that: In step 1, the active power time series P GA With P GB The entropy transferred between them is: Where: express and The joint probability density of express and The joint probability density of and represents conditional probability; in, Indicates the generator G A The active power at the n+1th sampling point, Indicates the generator G B The active power at the n+1th sampling point, Represents a time series P that satisfies a k-order Markov process GA , It represents the time series P that satisfies the i-th order Markov process GB ; P GB Take the (n+1)th sampling point, and P GA Take the conditional probability of the nth sampling point; P GA Take the (n+1)th sampling point, and P GB Take the conditional probability of the nth sampling point.

3. The method for locating the forced power oscillation disturbance source based on KDE transfer entropy causal analysis according to claim 1 is characterized in that: The joint probability density in step 2 is shown in formulas (2) and (3):

4. The method for locating the forced power oscillation disturbance source based on KDE transfer entropy causal analysis according to claim 3 is characterized in that: The multivariate kernel density estimation model in step 3 is shown in formulas (4) and (5): Where: H is the coefficient matrix, which is an i×i-dimensional symmetric positive definite matrix, and its coefficients are selected to meet the requirement of maximizing the accuracy of the joint probability density; K(.) is the multivariate kernel function; K(.) is represented by the function K(x) and meets the following requirements to minimize its impact on the accuracy of probability density estimation: Selecting Gaussian function as the multivariate kernel function K(x) further improves the accuracy of transfer entropy calculation: In formula (6) and formula (7), R i represents the set of real numbers, I i Represents the set of imaginary numbers, x i represents the time series, x represents the time series x i The sampling value of θ represents the time series x i The variance of .

5. The method for locating the forced power oscillation disturbance source based on KDE transfer entropy causal analysis according to claim 4 is characterized in that: The improved conditional probability in step 4 is expressed as: Where: H1, H2 and H3 are all symmetric positive definite coefficient matrices; K(.) is the multivariate kernel function, and the Gaussian function is selected here; N is the corresponding discrete number in the univariate kernel probability density estimation model, and m represents the generator G B The mth item of the active power time series, d represents the generator G A The dth term of the active power time series.

Images