CN111289796A - Method for detecting subsynchronous oscillation of high-proportion renewable energy power system - Google Patents

Method for detecting subsynchronous oscillation of high-proportion renewable energy power system Download PDF

Info

Publication number
CN111289796A
CN111289796A CN202010200063.7A CN202010200063A CN111289796A CN 111289796 A CN111289796 A CN 111289796A CN 202010200063 A CN202010200063 A CN 202010200063A CN 111289796 A CN111289796 A CN 111289796A
Authority
CN
China
Prior art keywords
time
frequency
signal
subsynchronous oscillation
synchronous compression
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN202010200063.7A
Other languages
Chinese (zh)
Other versions
CN111289796B (en
Inventor
蔡东升
黄琦
马钺
宫彦
张锋
亢朋朋
王衡
樊国伟
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
University of Electronic Science and Technology of China
State Grid Xinjiang Electric Power Co Ltd
Original Assignee
University of Electronic Science and Technology of China
State Grid Xinjiang Electric Power Co Ltd
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by University of Electronic Science and Technology of China, State Grid Xinjiang Electric Power Co Ltd filed Critical University of Electronic Science and Technology of China
Priority to CN202010200063.7A priority Critical patent/CN111289796B/en
Publication of CN111289796A publication Critical patent/CN111289796A/en
Application granted granted Critical
Publication of CN111289796B publication Critical patent/CN111289796B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01RMEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R23/00Arrangements for measuring frequencies; Arrangements for analysing frequency spectra
    • G01R23/02Arrangements for measuring frequency, e.g. pulse repetition rate; Arrangements for measuring period of current or voltage
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01RMEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R31/00Arrangements for testing electric properties; Arrangements for locating electric faults; Arrangements for electrical testing characterised by what is being tested not provided for elsewhere
    • G01R31/08Locating faults in cables, transmission lines, or networks
    • G01R31/081Locating faults in cables, transmission lines, or networks according to type of conductors
    • G01R31/086Locating faults in cables, transmission lines, or networks according to type of conductors in power transmission or distribution networks, i.e. with interconnected conductors
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01RMEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R31/00Arrangements for testing electric properties; Arrangements for locating electric faults; Arrangements for electrical testing characterised by what is being tested not provided for elsewhere
    • G01R31/08Locating faults in cables, transmission lines, or networks
    • G01R31/088Aspects of digital computing
    • HELECTRICITY
    • H02GENERATION; CONVERSION OR DISTRIBUTION OF ELECTRIC POWER
    • H02JCIRCUIT ARRANGEMENTS OR SYSTEMS FOR SUPPLYING OR DISTRIBUTING ELECTRIC POWER; SYSTEMS FOR STORING ELECTRIC ENERGY
    • H02J3/00Circuit arrangements for ac mains or ac distribution networks

Landscapes

  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • Power Engineering (AREA)
  • Mathematical Physics (AREA)
  • Theoretical Computer Science (AREA)
  • Complex Calculations (AREA)

Abstract

The invention discloses a method for detecting subsynchronous oscillation of a high-proportion renewable energy power system, which comprises the steps of firstly collecting a voltage signal or a current signal at the output end of the high-proportion renewable energy power system, then calculating a time frequency spectrum of short-time Fourier transform of the voltage signal or the current signal, mapping the time frequency spectrum of the short-time Fourier transform to a new time frequency spectrum by using a multiple synchronous compression algorithm, then determining a mode to which a coefficient after the multiple synchronous compression transform is higher than a set threshold value by using a mean shift method, calculating subsynchronous oscillation frequency, and finally performing signal reconstruction and least square estimation on the subsynchronous oscillation mode to obtain parameters of each subsynchronous oscillation mode.

Description

Method for detecting subsynchronous oscillation of high-proportion renewable energy power system
Technical Field
The invention belongs to the technical field of power systems, and particularly relates to a method for detecting subsynchronous oscillation of a high-proportion renewable energy power system.
Background
Subsynchronous Oscillation (SSO) is an abnormal phenomenon caused by the interaction between the converter-based generator and the grid of the power system. With the continuing development of environmental pollution and energy crisis, the penetration and use of renewable energy sources (in particular wind and photovoltaic) is increasing, but at the same time, many SSO events are triggered. In these incidents, the SSO becomes so destructive as to threaten the safety of the equipment and the stable operation of the system.
Nowadays, synchronous measurement techniques including Wide Area Measurement Systems (WAMS) and Phasor Measurement Units (PMU) have covered most power transmission networks and power plants and have been widely applied to dynamic monitoring of large power systems, which makes it a promising platform SSO monitoring. When SSO occurs, the voltage and power of the generator or grid will contain harmonics and inter-harmonics above the low frequency and below the nominal frequency (i.e., 50/60 Hz). The signals measured by PMUs in actual electrical applications often contain some noise, which requires that the analysis method have high mode resolution and noise immunity for SSO detection.
Some efforts have been made to develop efficient methods for SSO, including Fast Fourier Transform (FFT), hilbert-yellow transform (HHT), and Prony. In practical applications of the FFT method, a fence effect and spectral leakage occur. The HHT decomposes the signal into a set of eigenmode type functions (IMTs) by Empirical Mode Decomposition (EMD), and then extracts the information of each IMT by Hilbert Transform (HT). However, spectral aliasing often occurs during EMD decomposition, and inter-harmonic harmonics cannot be accurately detected [6 ]. The noise immunity of the Prony method is poor. In 2014, Oberlin et al. Synchronous compression transforms (SSTs) based on short-time fourier transforms (STFTs) are proposed. By rearranging the time-frequency information obtained by Fourier transform, the problem of time-frequency spectrum energy divergence is improved, and the anti-aliasing and anti-noise performance is better. But is not ideal for processing that emphasizes frequency signals, d. Pham and s.meignen propose higher order SST levels to obtain better results, but an increase in order brings more computation. Thus, G. Yu proposes a multiple synchronous compressive transform (MSST) algorithm that inherits the advantages of SST, increases the readability of the time-frequency spectrum, and has good processing effect on signals with strong time-varying and non-stationary characteristics. Furthermore, the calculated amplitude is similar to STFT, and the signal can be reconstructed almost perfectly. Therefore, MSST has great potential to improve SSO detection.
Disclosure of Invention
The invention aims to overcome the defects of the prior art and provides a method for detecting subsynchronous oscillation of a high-proportion renewable energy power system.
To achieve the above objects, the present invention
A method for detecting subsynchronous oscillation of a high-proportion renewable energy power system is characterized by comprising the following steps of:
(1) collecting a voltage signal or a current signal at the output end of the high-proportion renewable energy power system;
(2) calculating a time-frequency spectrum of a short-time Fourier transform of the voltage signal or the current signal;
(2.1), setting a time-varying signal model of the voltage signal or the current signal as s (t):
s(t)=A(t)eiφ(t)(1)
wherein A (t) represents the amplitude of the voltage signal or the current signal, phi (t) represents the phase of the voltage signal or the current signal, i is an imaginary unit, and t is time;
(2.2) the formula (1) is subjected to a first-order Taylor expansion formula to obtain:
Figure BDA0002419053040000021
wherein u is an integral variable, phi' (t) is expressed as a derivative of phase with time, and local instantaneous frequency is obtained;
further, the signal s (t) can be expressed as:
s(u)=A(t)ei(φ(t)+φ'(t)(u-t))(3)
computing the time-frequency spectrum G (t, ω) of the signal s (t) short-time fourier transform:
Figure BDA0002419053040000031
wherein g (-) is a window function,
Figure BDA0002419053040000032
is the fourier transform of a window function, ω is the frequency;
(3) mapping the time frequency spectrum of the short-time Fourier transform to a new time frequency spectrum by adopting a multiple synchronous compression algorithm;
(3.1), calculating the derivative of G (t, ω) with respect to time as:
Figure BDA0002419053040000033
(3.2) calculating the instantaneous frequency of any point on the time scale domain
Figure BDA0002419053040000034
Figure BDA0002419053040000035
(3.3), performing a first synchronous compression transform on G (t, ω):
Figure BDA0002419053040000036
where δ (·) is the Dirichlet function, η is the frequency corresponding to a certain time t, Ts[1](t, η) represents the synchronous compressed transform coefficients of s (t) at the time bin (t, η) when the synchronous compressed transform is first performed;
then, synchronous compression transformation is performed again by using the time frequency spectrum obtained by the formula (7), and the obtained coefficient after secondary synchronous compression transformation is as follows:
Figure BDA0002419053040000037
in the same way, the time frequency spectrum obtained by the formula (8) is used for executing synchronous compression transformation again, and the analogy is repeated until m times of synchronous compression transformation are executed, and the coefficient after multiple synchronous compression transformation is obtained:
Figure BDA0002419053040000038
(3.4) calculating an instantaneous frequency estimation value after the multiple synchronous compression transformation on the basis of the formula (9);
Figure BDA0002419053040000041
(3.5) sequentially arranging the coefficients subjected to the multiple synchronous compression transformation according to time t to form a mapped time-frequency spectrum;
(4) determining the mode of the coefficient after the multiple synchronous compression transformation higher than the set threshold value by using a mean shift method, and calculating the subsynchronous oscillation frequency;
(4.1) setting a threshold value To(ii) a Comparing coefficients after multiple simultaneous compression transformations with ToWill be lower than ToThe coefficient after the multiple synchronous compression transform is abandoned and is higher than ToMapping the multiple synchronous compression transformed coefficients into Xs(t,η):Ts[m](t, η) → η in which Xs(t, η) is the coefficient Ts[m](t, η) a set of all time frequency bins (t, η) clustered at frequency η;
(4.2) applying the mean shift method to Xs(t, η) clustering to calculate modal frequencies for each class
Figure BDA0002419053040000042
Let XsThe number of modes of (t, η) is N, χnRepresents class N, N is 1,2, …, N;
calculating modal frequencies of class n
Figure BDA0002419053040000043
Figure BDA0002419053040000044
Wherein, | · | represents solving an absolute value;
(4.3) when the working frequency is h0In Hz, the subsynchronous oscillation frequency range is [ h1,h2]Modal frequency of Hz
Figure BDA0002419053040000045
Recording as subsynchronous oscillation frequency, wherein the mode is a subsynchronous oscillation mode;
(5) reconstructing a signal;
according to the inverse transform of the multiple simultaneous compression, each mode is reconstructed by the following formula:
Figure BDA0002419053040000046
where b is the bandwidth of the multiple simultaneous compression transform, g (0) is the value of the window function at 0, phin' (t) denotes the derivative of the phase of the nth mode over time;
reconstructed signal of signal s (t)
Figure BDA0002419053040000047
Reconstructing a superposition of signals for the respective modalities, namely:
Figure BDA0002419053040000048
(6) carrying out least square estimation on the subsynchronous oscillation mode reconstruction signals to obtain parameters of each subsynchronous oscillation mode;
(6.1) defining a function Gn(Annn);
Figure BDA0002419053040000051
Wherein the content of the first and second substances,
Figure BDA0002419053040000052
for each subsynchronous oscillation mode, AnTo a corresponding amplitude, thetanIs a phase angle and epsilonnTo be damping coefficient, TsThe sampling interval of the current/voltage signal is shown, and k is 0,1,2, … are different sampling points;
(6.2) solving function G based on least square estimation algorithmnTo AnnnThe parameter values when the partial differential is equal to 0 respectively, so as to obtain the subsynchronous oscillation mode parameters:
Figure BDA0002419053040000053
and finally acquiring the subsynchronous oscillation frequency, amplitude, damping coefficient and phase angle.
The invention aims to realize the following steps:
the invention discloses a method for detecting subsynchronous oscillation of a high-proportion renewable energy power system, which comprises the steps of firstly collecting a voltage signal or a current signal at the output end of the high-proportion renewable energy power system, then calculating a time frequency spectrum of short-time Fourier transform of the voltage signal or the current signal, mapping the time frequency spectrum of the short-time Fourier transform to a new time frequency spectrum by using a multiple synchronous compression algorithm, then determining a mode to which a coefficient after the multiple synchronous compression transform is higher than a set threshold value by using a mean shift method, calculating subsynchronous oscillation frequency, and finally performing signal reconstruction and least square estimation on a subsynchronous oscillation mode to obtain parameters of each subsynchronous oscillation mode.
Meanwhile, the method for detecting the subsynchronous oscillation of the high-proportion renewable energy power system has the following beneficial effects:
(1) the method can accurately realize the detection of the subsynchronous oscillation, has good robustness to noise interference, and can more accurately extract the information of the subsynchronous oscillation by improving the calculation precision of each mode through iteration.
(2) The method provides a good data base for researching the generation mechanism of the subsynchronous oscillation in the new energy convergence region, and has reference value for the research of the control strategy.
(3) The invention develops a subsynchronous oscillation detection scheme based on phasor measurement by utilizing the most advanced multi-mode signal analysis tool.
(4) The present invention has great potential in improving the detection of subsynchronous oscillations and the development of countermeasures against the occurring subsynchronous control interactions.
Drawings
FIG. 1 is a flow chart of a method for detecting sub-synchronous oscillation in a high-proportion renewable energy power system according to the present invention;
FIG. 2 is a diagram of multiple simultaneous compression transforms of a signal.
Detailed Description
The following description of the embodiments of the present invention is provided in order to better understand the present invention for those skilled in the art with reference to the accompanying drawings. It is to be expressly noted that in the following description, a detailed description of known functions and designs will be omitted when it may obscure the subject matter of the present invention.
Examples
Fig. 1 is a flow chart of a method for detecting sub-synchronous oscillation of a high-proportion renewable energy power system according to the invention.
In this embodiment, as shown in fig. 1, the method for detecting sub-synchronous oscillation of a high-proportion renewable energy power system of the present invention includes the following steps:
s1, signal acquisition;
collecting a voltage signal or a current signal at the output end of a high-proportion renewable energy power system;
s2, calculating a time frequency spectrum of the short-time Fourier transform of the voltage signal or the current signal;
s2.1, setting a time-varying signal model of the voltage signal or the current signal as S (t):
s(t)=A(t)eiφ(t)(1)
wherein A (t) represents the amplitude of the voltage signal or the current signal, phi (t) represents the phase of the voltage signal or the current signal, i is an imaginary unit, and t is time;
s2.2, performing a first-order Taylor expansion equation on the formula (1) to obtain:
Figure BDA0002419053040000061
wherein u is an integral variable, phi' (t) is expressed as a derivative of phase with time, and local instantaneous frequency is obtained;
further, the signal s (t) can be expressed as:
s(u)=A(t)ei(φ(t)+φ'(t)(u-t))(3)
computing the time-frequency spectrum G (t, ω) of the signal s (t) short-time fourier transform:
Figure BDA0002419053040000071
wherein g (-) is a window function,
Figure BDA0002419053040000072
is the fourier transform of a window function, ω is the frequency;
s3, mapping the time-frequency spectrum of the short-time Fourier transform to a new time-frequency spectrum by adopting a multiple synchronous compression algorithm;
s3.1, calculating the derivative of G (t, ω) with respect to time as:
Figure BDA0002419053040000073
s3.2, calculating the instantaneous frequency of any point on the time scale domain
Figure BDA0002419053040000074
Figure BDA0002419053040000075
S3.3, the energy of a time-frequency spectrum G (t, omega) obtained by short-time Fourier transform is seriously dispersed, so that the observation is unclear, and the first synchronous compression transform is performed on the G (t, omega):
Figure BDA0002419053040000076
where δ (·) is the Dirichlet function, η is the frequency corresponding to a certain time t, Ts[1](t, η) represents the synchronous compressed transform coefficients of s (t) at the time bin (t, η) when the synchronous compressed transform is first performed;
however, experiments show that the synchronous compression coefficient is not suitable for strong time-varying signal processing, so that when sub-synchronous oscillation occurs, the method of performing synchronous compression transformation on the time frequency spectrum is continuously executed, and an equation of multiple synchronous compression transformation coefficients can be obtained, that is, the time frequency spectrum obtained by the equation (7) is executed again for synchronous compression transformation, and the coefficients obtained after secondary synchronous compression transformation are substituted as follows:
Figure BDA0002419053040000081
in the same way, the time frequency spectrum obtained by the formula (8) is used for executing synchronous compression transformation again, and the analogy is repeated until m times of synchronous compression transformation are executed, and the coefficient after multiple synchronous compression transformation is obtained:
Figure BDA0002419053040000082
s3.4, calculating an instantaneous frequency estimation value after multiple synchronous compression transformation on the basis of the formula (9);
Figure BDA0002419053040000083
for each iteration, the multiple synchronous compression transformation constructs a new instantaneous frequency estimation value to redistribute the energy of the fuzzy short-time Fourier transformation, obviously, the frequency estimation value of the multiple synchronous compression transformation is closer to the actual instantaneous frequency of the signal through multiple iterations;
s3.5, mapping the signals from a time scale domain to a time-frequency domain based on frequency combination so as to obtain time-frequency information of the signals, and sequentially arranging coefficients subjected to multiple synchronous compression transformations according to time t to form a mapped time-frequency spectrum;
s4, determining the mode of the coefficient after multiple synchronous compression transformation higher than a set threshold value by using a mean shift method, and calculating the subsynchronous oscillation frequency;
s4.1, setting a threshold value ToIn this embodiment, the threshold ToSet to 0;
from the analysis results of the multiple simultaneous compressively transforming processes, it can be seen that the multiple simultaneous compressively transforming coefficients from the same pattern are concentrated around the same frequency, and thus, we can make the lower than threshold ToWill be ignored and will be above the threshold ToMapping the multiple synchronous compression transformed coefficients into Xs(t,η):Ts[m](t, η) → η in which Xs(t, η) is the coefficient Ts[m](t, η) a set of all time frequency bins (t, η) clustered at frequency η;
s4.2, using mean shift method to Xs(t, η) clustering to calculate modal frequencies for each class
Figure BDA0002419053040000084
Let XsThe number of modes of (t, η) is N, χnRepresents class N, N is 1,2, …, N;
calculating modal frequencies of class n
Figure BDA0002419053040000091
The modal frequency can be calculated as the centroid of the corresponding aggregated multiple synchronous compression transform coefficient, and the specific calculation formula is as follows:
Figure BDA0002419053040000092
wherein, | · | represents solving an absolute value;
s4.3, when the working frequency is 50Hz, the frequency range is [5,45 ]]Hz (when the working frequency is 60Hz, the subsynchronous oscillation frequency range is [5,55 ]]Hz) modal frequency
Figure BDA0002419053040000093
Recording as subsynchronous oscillation frequency, wherein the mode is a subsynchronous oscillation mode;
s5, signal reconstruction
Considering that the multiple simultaneous transforms redistribute the time-frequency coefficients only in the frequency direction and no information is lost, the multiple simultaneous transforms can perfectly reconstruct the signal, so each mode can be reconstructed by the following formula:
Figure BDA0002419053040000094
wherein b is the bandwidth of multiple synchronous compression transformation, in this embodiment, the bandwidth is 50 Hz; g (0) is the value of the window function at 0, phin' (t) denotes the derivative of the phase of the nth mode over time;
reconstructed signal of signal s (t)
Figure BDA0002419053040000095
Reconstructing a superposition of signals for the respective modalities, namely:
Figure BDA0002419053040000096
s6, carrying out least square estimation on the subsynchronous oscillation mode reconstruction signal to obtain parameters of each subsynchronous oscillation mode;
s6.1, defining function Gn(Annn);
According to the least square method, the best matching parameter is solved by minimizing the square sum of the modal reconstruction signal and the assumed parameter difference, i.e. the modal parameter can be formulated as the following problem:
Figure BDA0002419053040000097
wherein the content of the first and second substances,
Figure BDA0002419053040000101
the oscillation frequency for each subsynchronous oscillation mode; a. thenTo a corresponding amplitude, thetanIs a phase angle and epsilonnIs a damping coefficient; t issFor the sampling interval of the current/voltage signal, in the present embodiment, the sampling interval is set to 8.33 milliseconds, that is, the sampling frequency is 120 Hz; k is 0,1,2, … is a different sample point;
s6.2, solving a function G based on a least square estimation algorithmnTo AnnnThe parameter values when the partial differential is equal to 0 respectively, so as to obtain the subsynchronous oscillation mode parameters:
Figure BDA0002419053040000102
and finally acquiring the subsynchronous oscillation frequency, amplitude, damping coefficient and phase angle.
In order to facilitate the understanding of the technical contents of the invention by those skilled in the art, the following further explains the technical contents of the invention with reference to the drawings. The invention performs parameter identification of subsynchronous oscillation on an example signal of the IEEE second standard model (operating frequency is 60Hz) as shown in table 1.
Figure BDA0002419053040000103
TABLE 1
The signal waveform is as shown in fig. 2(a), similar to the set of subsynchronous oscillation components, which are caused by interference occurring at 0.5s, and then cleared after 0.1 s. Then, STFT, SST, and MSST are performed on the signals.
Fig. 2(b), fig. 2(c) and fig. 2(d) show results after the STFT, SST and MSST processing, and fig. 2(b) and fig. 2(c) are graphs showing SST and MSST coefficients. It can be seen that the perturbed instantaneous frequency components are concentrated at three frequencies. It can be clearly seen that the "pseudo" time-frequency representation consisting of STFT and SST, as shown in fig. 2(b) and 2(c), is a blurred image due to "interference" from STFT and SST. The primary frequency was 51 Hz. To avoid interference from other frequencies, filters are used to decompose the signal into several sub-bands. Compared with other methods, the result of MSST signal processing has no pseudo time-frequency phenomenon, and the image is clearest. Thus, the process does not require those pretreatment measures.
From the analysis results in fig. 2(d) of MSST processing, the number of patterns can be easily determined by visual observation. However, for more efficient recognition, an algorithm that automatically obtains the number of patterns is needed. As can be seen from fig. 2(d), MSST coefficients from the same pattern are all concentrated around the same frequency. Therefore, we can define the mapping as a feature where the MSST coefficients are above a set threshold; all other MSST coefficients below this threshold will be ignored. Then, clustering analysis is carried out by using the characteristics, and the number of the modes to which the nonzero MSST coefficient belongs and the frequency corresponding to the modes are determined.
As shown in fig. 2(d), it can also reduce computational complexity by using a linear kernel function when applying the mean-shift algorithm, considering that all clusters can be clearly distinguished. The mode number is determined by a mean shift method, and the mode frequency can be calculated as the centroid of the corresponding MSST coefficient cluster. To evaluate the detection frequency
Figure BDA0002419053040000111
The accuracy of the result, we adopt the estimation error EE, consisting of:
Figure BDA0002419053040000112
wherein f isnIs the actual frequency;
for each detected modality, calculating a frequency
Figure BDA0002419053040000113
And EE and recorded in table 2. As can be seen from table 2, the modal frequency estimation is accurate, and the estimation errors of the three modes are all less than 1%.
Figure BDA0002419053040000114
TABLE 2
Reconstructing each modal signal, introducing a signal-to-noise ratio (SNR) for evaluating the accuracy of the reconstructed signal, and calculating the following formula:
Figure BDA0002419053040000115
wherein s isnFor the original signal of each of the modalities,
Figure BDA0002419053040000116
the reconstructed signal for each mode found for equation (12) is the reconstructed signal for each detected subsynchronous oscillation mode in the example and is given in table 3. In order to verify the noise immunity of the method, 20dB of white Gaussian noise is added into the signal. The reconstruction accuracy is shown in the third row of the table, indicating that the method can accurately reconstruct the individual components of the signal.
Figure BDA0002419053040000121
TABLE 3
Performing parameter identification on the reconstructed subsynchronous oscillation mode by adopting a least square method, wherein the obtained result is shown in table 4;
Figure BDA0002419053040000122
TABLE 4
Although illustrative embodiments of the present invention have been described above to facilitate the understanding of the present invention by those skilled in the art, it should be understood that the present invention is not limited to the scope of the embodiments, and various changes may be made apparent to those skilled in the art as long as they are within the spirit and scope of the present invention as defined and defined by the appended claims, and all matters of the invention which utilize the inventive concepts are protected.

Claims (2)

1. A method for detecting subsynchronous oscillation of a high-proportion renewable energy power system is characterized by comprising the following steps of:
(1) collecting a voltage signal or a current signal at the output end of the high-proportion renewable energy power system;
(2) calculating a time-frequency spectrum of a short-time Fourier transform of the voltage signal or the current signal;
(2.1), setting a time-varying signal model of the voltage signal or the current signal as s (t):
s(t)=A(t)eiφ(t)(1)
wherein A (t) represents the amplitude of the voltage signal or the current signal, phi (t) represents the phase of the voltage signal or the current signal, i is an imaginary unit, and t is time;
(2.2) the formula (1) is subjected to a first-order Taylor expansion formula to obtain:
Figure FDA0002419053030000011
wherein u is an integral variable, phi' (t) is expressed as a derivative of phase with time, and local instantaneous frequency is obtained;
further, the signal s (t) can be expressed as:
s(u)=A(t)ei(φ(t)+φ'(t)(u-t))
(3)
computing the time-frequency spectrum G (t, ω) of the signal s (t) short-time Fourier transform:
Figure FDA0002419053030000012
wherein g (-) is a window function,
Figure FDA0002419053030000013
is the fourier transform of the window function;
(3) mapping the time frequency spectrum of the short-time Fourier transform to a new time frequency spectrum by adopting a multiple synchronous compression algorithm;
(3.1), calculating the derivative of G (t, ω) with respect to time as:
Figure FDA0002419053030000014
(3.2) calculating the instantaneous frequency of any point on the time scale domain
Figure FDA0002419053030000021
Figure FDA0002419053030000022
(3.3), performing a first synchronous compression transform on G (t, ω):
Figure FDA0002419053030000023
where δ (·) is the Dirichlet function, η is the frequency corresponding to a certain time t, Ts[1](t, η) represents the synchronous compressed transform coefficients of s (t) at the time bin (t, η) when the synchronous compressed transform is first performed;
then, synchronous compression transformation is performed again by using the time frequency spectrum obtained by the formula (7), and the obtained coefficient after secondary synchronous compression transformation is as follows:
Figure FDA0002419053030000024
in the same way, the time frequency spectrum obtained by the formula (8) is used for executing synchronous compression transformation again, and the analogy is repeated until m times of synchronous compression transformation are executed, and the coefficient after multiple synchronous compression transformation is obtained:
Figure FDA0002419053030000025
(3.4) calculating an instantaneous frequency estimation value after the multiple synchronous compression transformation on the basis of the formula (9);
Figure FDA0002419053030000026
(3.5) sequentially arranging the coefficients subjected to the multiple synchronous compression transformation according to time t to form a mapped time-frequency spectrum;
(4) determining the mode of the coefficient after the multiple synchronous compression transformation higher than the set threshold value by using a mean shift method, and calculating the subsynchronous oscillation frequency;
(4.1) setting a threshold value To(ii) a Comparing coefficients after multiple simultaneous compression transformations with ToWill be lower than ToThe coefficient after the multiple synchronous compression transform is abandoned and is higher than ToMapping the multiple synchronous compression transformed coefficients into Xs(t,η):Ts[m](t, η) → η in which Xs(t, η) is the coefficient Ts[m](t, η) a set of all time frequency bins (t, η) clustered at frequency η;
(4.2) applying the mean shift method to Xs(t, η) clustering to calculate modal frequencies for each class
Figure FDA0002419053030000027
Let XsThe number of modes of (t, η) is N, χnDenotes the nth class, n is 1,2, …,N;
Calculating modal frequencies of class n
Figure FDA0002419053030000031
Figure FDA0002419053030000032
Wherein, | · | represents solving an absolute value;
(4.3) when the working frequency is h0In Hz, the subsynchronous oscillation frequency range is [ h1,h2]Modal frequency of Hz
Figure FDA0002419053030000033
Recording as subsynchronous oscillation frequency, wherein the mode is a subsynchronous oscillation mode;
(5) reconstructing a signal;
according to the inverse transform of the multiple simultaneous compression, each mode is reconstructed by the following formula:
Figure FDA0002419053030000034
where b is the bandwidth of the multiple synchronous compression transform, g (0) is the value of the window function at 0, phi'n(t) represents the derivative of the phase of the nth mode over time;
reconstructed signal of signal s (t)
Figure FDA0002419053030000035
Reconstructing a superposition of signals for the respective modalities, namely:
Figure FDA0002419053030000036
(6) carrying out least square estimation on the subsynchronous oscillation mode reconstruction signals to obtain parameters of each subsynchronous oscillation mode;
(6.1) defining a function Gn(Annn);
Figure FDA0002419053030000037
Wherein the content of the first and second substances,
Figure FDA0002419053030000038
for each subsynchronous oscillation mode, AnTo a corresponding amplitude, thetanIs a phase angle and epsilonnTo be damping coefficient, TsThe sampling interval of the current/voltage signal is shown, and k is 0,1,2, … are different sampling points;
(6.2) solving function G based on least square estimation algorithmnTo AnnnThe parameter values when the partial differential is equal to 0 respectively, so as to obtain the subsynchronous oscillation mode parameters:
Figure FDA0002419053030000041
and finally acquiring the subsynchronous oscillation frequency, amplitude, damping coefficient and phase angle.
2. The method according to claim 1, wherein the operating frequency h is a frequency of sub-synchronous oscillation of the high-proportion renewable energy power system0Sub-synchronous oscillation frequency range [ h ] at 50Hz1,h2]=[5,45]Hz; the working frequency h0Sub-synchronous oscillation frequency range [ h ] at 60Hz1,h2]=[5,55]Hz。
CN202010200063.7A 2020-03-20 2020-03-20 Method for detecting subsynchronous oscillation of high-proportion renewable energy power system Active CN111289796B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010200063.7A CN111289796B (en) 2020-03-20 2020-03-20 Method for detecting subsynchronous oscillation of high-proportion renewable energy power system

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010200063.7A CN111289796B (en) 2020-03-20 2020-03-20 Method for detecting subsynchronous oscillation of high-proportion renewable energy power system

Publications (2)

Publication Number Publication Date
CN111289796A true CN111289796A (en) 2020-06-16
CN111289796B CN111289796B (en) 2021-03-30

Family

ID=71025856

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010200063.7A Active CN111289796B (en) 2020-03-20 2020-03-20 Method for detecting subsynchronous oscillation of high-proportion renewable energy power system

Country Status (1)

Country Link
CN (1) CN111289796B (en)

Cited By (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111856562A (en) * 2020-07-30 2020-10-30 成都理工大学 Generalized high-order synchronous extrusion seismic signal time-frequency decomposition and reconstruction method
CN112186759A (en) * 2020-09-28 2021-01-05 西安热工研究院有限公司 Doubly-fed wind power plant subsynchronous oscillation suppression method for adaptively capturing frequency points
CN112446290A (en) * 2020-10-19 2021-03-05 电子科技大学 Key parameter identification method for wind power subsynchronous oscillation
CN112505452A (en) * 2020-11-25 2021-03-16 东南大学 Wide-area system broadband oscillation monitoring method
CN113489000A (en) * 2021-07-14 2021-10-08 河海大学 Method for inhibiting low-frequency oscillation of water-light complementary energy base

Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
DE3604996A1 (en) * 1986-02-17 1987-08-20 Siemens Ag Method for detecting sub-synchronous oscillations on AC voltage cables, and devices for carrying out the method
US8564379B2 (en) * 2008-10-15 2013-10-22 Robert Bosch Gmbh Device and method for testing a frequency-modulated clock generator
CN104821579A (en) * 2015-04-24 2015-08-05 中国电力科学研究院 Convertor station electrical signals-based subsynchronous oscillation monitoring analysis method
CN106199183A (en) * 2016-08-16 2016-12-07 国电南瑞科技股份有限公司 A kind of PMU realizing sub-synchronous oscillation on-line identification alarm and method
CN106383270A (en) * 2016-08-26 2017-02-08 清华大学 Wide-area measurement information based electric power system sub-synchronous oscillation monitoring method and system
CN109683057A (en) * 2017-10-18 2019-04-26 中国电力科学研究院 A kind of online disturbance source locating method and system of subsynchronous oscillation of electrical power system

Patent Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
DE3604996A1 (en) * 1986-02-17 1987-08-20 Siemens Ag Method for detecting sub-synchronous oscillations on AC voltage cables, and devices for carrying out the method
US8564379B2 (en) * 2008-10-15 2013-10-22 Robert Bosch Gmbh Device and method for testing a frequency-modulated clock generator
CN104821579A (en) * 2015-04-24 2015-08-05 中国电力科学研究院 Convertor station electrical signals-based subsynchronous oscillation monitoring analysis method
CN106199183A (en) * 2016-08-16 2016-12-07 国电南瑞科技股份有限公司 A kind of PMU realizing sub-synchronous oscillation on-line identification alarm and method
CN106383270A (en) * 2016-08-26 2017-02-08 清华大学 Wide-area measurement information based electric power system sub-synchronous oscillation monitoring method and system
CN109683057A (en) * 2017-10-18 2019-04-26 中国电力科学研究院 A kind of online disturbance source locating method and system of subsynchronous oscillation of electrical power system

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
YAN GONG等: "Analysis on Oscillation Propagation Characteristics based on Impedance Model", 《IEEE》 *
喻敏: "基于同步挤压小波变换的电能质量扰动分析研究", 《中国博士学位论文全文数据库工程科技Ⅱ辑》 *

Cited By (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111856562A (en) * 2020-07-30 2020-10-30 成都理工大学 Generalized high-order synchronous extrusion seismic signal time-frequency decomposition and reconstruction method
CN111856562B (en) * 2020-07-30 2022-07-26 成都理工大学 Generalized high-order synchronous extrusion seismic signal time-frequency decomposition and reconstruction method
CN112186759A (en) * 2020-09-28 2021-01-05 西安热工研究院有限公司 Doubly-fed wind power plant subsynchronous oscillation suppression method for adaptively capturing frequency points
CN112186759B (en) * 2020-09-28 2022-11-15 西安热工研究院有限公司 Doubly-fed wind power plant subsynchronous oscillation suppression method for adaptively capturing frequency points
CN112446290A (en) * 2020-10-19 2021-03-05 电子科技大学 Key parameter identification method for wind power subsynchronous oscillation
CN112505452A (en) * 2020-11-25 2021-03-16 东南大学 Wide-area system broadband oscillation monitoring method
CN112505452B (en) * 2020-11-25 2024-05-24 东南大学 Wide-band oscillation monitoring method for wide-area system
CN113489000A (en) * 2021-07-14 2021-10-08 河海大学 Method for inhibiting low-frequency oscillation of water-light complementary energy base
CN113489000B (en) * 2021-07-14 2023-08-08 河海大学 Method for restraining low-frequency oscillation of water-light complementary energy base

Also Published As

Publication number Publication date
CN111289796B (en) 2021-03-30

Similar Documents

Publication Publication Date Title
CN111289796B (en) Method for detecting subsynchronous oscillation of high-proportion renewable energy power system
Gaouda et al. Application of multiresolution signal decomposition for monitoring short-duration variations in distribution systems
CN109638862B (en) CEEMDAN algorithm-based electric power system low-frequency oscillation mode identification method
CN110648088B (en) Electric energy quality disturbance source judgment method based on bird swarm algorithm and SVM
CN110068759A (en) A kind of fault type preparation method and device
CN105572501A (en) Power quality disturbance identification method based on SST conversion and LS-SVM
CN109659957B (en) APIT-MEMD-based power system low-frequency oscillation mode identification method
Li et al. Period-assisted adaptive parameterized wavelet dictionary and its sparse representation for periodic transient features of rolling bearing faults
Bruna et al. Selection of the most suitable decomposition filter for the measurement of fluctuating harmonics
CN115618213A (en) Charger voltage disturbance analysis method, system, equipment and storage medium
Mei et al. Wavelet packet transform and improved complete ensemble empirical mode decomposition with adaptive noise based power quality disturbance detection
Lv et al. Longitudinal synchroextracting transform: A useful tool for characterizing signals with strong frequency modulation and application to machine fault diagnosis
CN109635430A (en) Grid power transmission route transient signal monitoring method and system
CN105606892A (en) Power grid harmonic and inter-harmonic analysis method based on SST transformation
CN110147637A (en) Based on the small impact-rub malfunction diagnostic method for involving the greedy sparse identification of harmonic components
CN111650436B (en) Subsynchronous oscillation identification method for high-proportion renewable energy power system
Thomas et al. Machine learning based detection and classification of power system events
Li et al. A Supraharmonic measurement method based on matrix pencil with high-order mixed mean cumulant
Lin et al. Denoising of wind speed data by wavelet thresholding
Irhamsyah et al. Parameters Quality Performance of Signal Fluctuation Based on Data Grouping
Xingang et al. Supraharmonics measurement algorithm based on CS-SAMP
He et al. Fault diagnosis for electronic transformer with wavelet energy entropy and SVM
Lu et al. Subsynchronous oscillation detection using synchrosqueezing wavelet transforms and K-means clustering
Yang et al. New Measurement Algorithm for Supraharmonic Real-time monitoring Based on Dynamic Compressed Sensing
Kapisch et al. Novelty Detection in Power Quality Signals with Surrogates: a Time-Frequency Technique

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant