Subsynchronous oscillation mode identification method and device for power system containing constructed network type converter based on matrix pencil and LM, and storage medium

Abstract

The invention discloses a subsynchronous oscillation mode identification method and device for a power system containing a constructed network converter based on a matrix pencil and an LM, and a storage medium, and belongs to the technical field of subsynchronous oscillation mode identification in the power system. The method comprises the following steps: performing mode parameter identification by using a matrix pencil method to obtain an initial value of each mode parameter, and performing iteration by using an LM method to solve an optimal value of each mode parameter. A trust region is introduced into the LM algorithm. And an optimization problem is changed from an unconstrained least square problem to a constrained least square problem, so that the identification precision is improved. Compared with a traditional method, the method provided by the invention has the advantages that the calculation is simpler, the accuracy is higher, and the method has a better development prospect in subsynchronous oscillation early warning.

Description

Technical Field

[0001] The present invention belongs to the technical field of subsynchronous oscillation mode identification in power systems, and in particular relates to a method, device and computer-readable storage medium for identifying subsynchronous oscillation mode parameters by combining a matrix bundle method and a Levenberg-Marquardt (hereinafter referred to as LM) algorithm. Background Art

[0002] At present, vigorously developing new energy is a major trend in the development of power grids, but large-scale new energy grid connection is bound to cause system stability problems. In addition, the grid connection of new energy units leads to the connection of a large number of power electronic equipment, resulting in a decrease in grid inertia. Traditional static and transient stability problems, as well as new problems such as parameter coupling problems, have become problems that need to be solved urgently.

[0003] The connection of a large number of new energy sources to the grid has led to the connection of a large number of power converter devices to the grid. Due to its low inertia characteristics, it is difficult to maintain the inertia of the synchronous generator in the original power system. At the same time, its control links are coupled with each other, resulting in more new stability problems. The widely used converters are mainly divided into two types: grid-following converters and grid-building converters. The grid-following converter mainly relies on the phase-locked loop to achieve synchronization with the grid, while the grid-building converter mainly adopts power synchronization and has a certain voltage and frequency support capability for the grid. However, the converter is ultimately a power electronic device. Unlike the inertia of the synchronous generator, its large-scale grid connection is bound to cause various stability problems such as subsynchronous oscillation.

[0004] The nonlinear and multi-coupled characteristics of power electronic equipment make the eigenvalue method suffer from the "curse of dimensionality". With the development of science and technology, the wide area measurement system (WAMS) has been widely used in power systems. Many experts and scholars use the grid data measured by WAMS to identify subsynchronous oscillation parameters, and related methods have also been further studied.

[0005] The matrix bundle method is a modal parameter estimation method proposed in the late 1980s and early 1990s. It has been widely used in aircraft, automobiles, bridges and power systems. When there is noise in the identification signal, the statistical variance of the system identification parameters obtained by the matrix bundle method is always better than the Prony method.

[0006] LM method is an estimation method for least squares estimation of regression parameters in nonlinear regression. It was proposed by DW Marquardt in 1963, and he developed it based on a paper by K. Levenbevg in 1944. This method is a combination of the steepest descent method and the linearization method (Taylor series). Because the steepest descent method is applicable to the situation where the parameter estimate is far from the optimal value at the beginning of the iteration, while the linearization method, that is, the Gauss-Newton method, is applicable to the later stage of the iteration, when the parameter estimate is close to the optimal value.

[0007] At present, there is no report on combining the matrix bundle and LM methods for subsynchronous oscillation mode identification of power systems containing grid-type converters to improve identification accuracy. Therefore, a new technical solution is urgently needed in the prior art to solve this problem. Summary of the invention

[0008] The technical problem to be solved by the present invention is: to provide a method, device and storage medium for identifying subsynchronous oscillation modes of a power system containing a grid-type converter based on matrix bundle and LM, so as to solve the technical problem of related reports that the matrix bundle and LM methods are not combined for subsynchronous oscillation mode identification of a power system containing a grid-type converter to improve the identification accuracy.

[0009] The method for identifying subsynchronous oscillation modes of a power system containing a grid-type converter based on a matrix bundle and LM comprises the following steps, and the following steps are performed in sequence:

[0010] Step 1: extracting eigenvalues ​​of oscillation information identification parameters of a given signal by using a matrix bundle method, and taking a vector composed of the extracted eigenvalues ​​as an initial value, wherein the oscillation information identification parameters include amplitude, attenuation coefficient, frequency and phase;

[0011] Step 2: Use the Levenberg-Marquardt algorithm to iteratively optimize the vector composed of each eigenvalue starting from the initial value, and the optimal solution obtained is the subsynchronous oscillation mode identification result of the power system.

[0012] The specific method of extracting the characteristic value of the oscillation information identification parameter of a given signal by using the matrix bundle method is:

[0013] S101. The matrix bundle method fits the function of a given signal by constructing a fitting function, which is expressed as:

[0014]

[0015] Where: A i represents the amplitude, σ i represents the attenuation factor, ω i Indicates frequency, Indicates phase;

[0016] The function of the given signal is expressed as:

[0017] y(t k )=y(k),k=0,1,…,N-1 (2);

[0018] In the formula, y(k) represents the kth sampling value among the sampling values ​​with equal time intervals; N is the total number of sampling values;

[0019] Rewrite equation (1) into the following complex exponential form:

[0020]

[0021] Where: λ i The characteristic value representing the identification parameter of the oscillation information;

[0022] Simplify formula (3) and express it as:

[0023]

[0024] In formula (4),

[0025]

[0026] S102, perform extreme point z i Solution:

[0027] S103, characteristic value λ of oscillation information identification parameter i Obtained by the following formula:

[0028]

[0029] S104, Solve B i , so that it satisfies that for all k

[0030] From formula (4), we can get:

[0031]

[0032] Formula (17) can be simplified as:

[0033] ZB=Y(18);

[0034] The matrix B is N × n, so equation (18) is solved by the least squares method;

[0035] S105, solve according to formula (3) and obtain the fitted That is, the extracted eigenvalues.

[0036] The pole z iThe solution process is as follows:

[0037] 1) The matrix bundle is given by:

[0038] [Y2]-λ[Y1]=[Z1][B]{[Z0]-λ[I]}[Z2] (6); where [B] represents the residual matrix; [I] represents the n×n unit matrix;

[0039] [Z0] = diag[z1,z2,…,z n ] (7);

[0040]

[0041]

[0042] Where n represents the desired number of eigenvalues; L represents the beam parameter;

[0043] 2) The steps of the matrix bundle method are as follows:

[0044] ①Select L, satisfying n≤L≤Nn;

[0045] ②Construct matrix [Y]:

[0046]

[0047] Wherein, the matrix [Y] is the matrix composed of the observation values ​​at each sampling time;

[0048] ③Perform singular value decomposition on [Y]:

[0049] [Y]=[U][S][V] T (11); where [U] and [V] are unitary matrices, containing [Y][Y] respectively. T and [Y] T The eigenvector of [Y];

[0050] ④Construct matrices [V1] and [V2] and satisfy:

[0051] V1=[v1 v2 v3 … v n-1 ] (12);

[0052] V2=[v2 v3 v4 … v n ] (13);

[0053] In the formula, v n-1 is the n-1th right singular eigenvector; v n is the nth right singular eigenvector;

[0054] ⑤Construct matrices [Y1] and [Y2]:

[0055] [Y1]=[V1] T [V1] (14);

[0056] [Y2]=[V2] T [V1] (15);

[0057] ⑥ The pole zi to be solved is calculated as the eigenvalue of the matrix pair {[Y1]; [Y2]}.

[0058] The Levenberg-Marquardt algorithm introduces a trust region, which changes the optimization problem from an unconstrained least squares problem to a constrained least squares problem, which is a nonlinear least squares problem. The specific steps are as follows:

[0059] 1) Minimizing the sum of squares of errors is used as the expression of nonlinear least squares problem:

[0060]

[0061] In the formula, y i Is the system at t i The output measurement value at time, x is the vector composed of the identification parameters of each oscillation information obtained by the fitting function in the matrix bundle method;

[0062] 2) Solve equation (19):

[0063] The Taylor series expansion of the function f(x) at x0 is:

[0064] min f(x0+Δx)=||f(x0)+J T Δx|| 2 (20);

[0065] Where J is the Jacobian matrix, Δx is the small perturbation of the function f(x) at x0;

[0066]

[0067] In formula (21), x n ,y n are the identification parameter value and output measurement value at the nth sampling time respectively;

[0068] Derivate the minimization function of equation (20), set the derivative to 0, and find the extreme value of the error function:

[0069] JJ T Δx=-Jf(x0) (23);

[0070] The iteration formula is:

[0071] Δx=-(J T J)-1 J T f(x0) (24);

[0072] Therefore, the matrix H is introduced to ensure the continued feasibility of iterative operations:

[0073] H≈J T J+μI(25);

[0074] In the formula, μ is the trust region, I is the unit matrix;

[0075] The update iteration formula is:

[0076] x k+1 =x k -(J T J+μI) -1 g k (26);

[0077] In the formula, g k is the negative gradient direction, x k is the iteration value of the kth iteration;

[0078] 3) Perform iterative operations:

[0079] ① Given the initial value x0 and the initial optimization radius μ1:

[0080] ②For the kth iteration, add the trust region to the Gauss-Newton method and solve:

[0081]

[0082] In formula (27), f(x k ) is the fitting function value of the kth iteration; Δx k is the step length;

[0083] ③ Calculate and obtain the quantitative index ρ to measure the degree of approximation:

[0084]

[0085] ④. Then set the next optimization radius μ1′=2μ1, and let x k+1 =x k +Δx k Iterate

[0086] Then set the next optimization radius μ1′=0.5μ1, let x k+1 =x k +Δx k Iterate

[0087] ⑤. Determine whether the above algorithm converges. If not, return to step ②. If converged, end. The optimal solution obtained is the subsynchronous oscillation mode identification result of the power system.

[0088] The subsynchronous oscillation mode identification device of a grid-type converter power system based on a matrix bundle and LM, using the subsynchronous oscillation mode identification method of a grid-type converter power system based on a matrix bundle and LM, comprises:

[0089] An eigenvalue extraction module is used to extract the eigenvalues ​​of the oscillation information identification parameters of a given signal by using a matrix bundle method, and use the vector composed of the extracted eigenvalues ​​as an initial value;

[0090] The Levenberg-Marquardt algorithm module is used to iteratively optimize the vector composed of each eigenvalue starting from the initial value to obtain the optimal solution.

[0091] The eigenvalue extraction module includes a fitting unit, a pole z i Solving for unit, eigenvalue λ i Solving unit and variable B with identification parameters i Solving unit;

[0092] The fitting unit is used to realize the following functions:

[0093] Construct a fitting function to fit the function of the given signal. The fitting function is expressed as:

[0094]

[0095] Where: A i represents the amplitude, σ i represents the attenuation factor, ω i Indicates frequency, Indicates phase;

[0096] The function of the given signal is expressed as:

[0097] y(t k )=y(k),k=0,1,…,N-1 (2);

[0098] In the formula, y(k) represents the kth sample value among the samples with equal time intervals; N is the total number of samples;

[0099] Rewrite equation (1) into the following complex exponential form:

[0100]

[0101] Where: λ i The characteristic value representing the identification parameter of the oscillation information;

[0102] Simplify formula (3) and express it as:

[0103]

[0104] In formula (4),

[0105]

[0106] In formula (5), z i represents the system pole;

[0107] The pole z i Solving for unit, eigenvalue λ i Solving unit and variable B with identification parameters i Substitute the solution result of the solving unit into the fitting function to obtain the fitted That is, the extracted eigenvalues.

[0108] The pole z i Solver element for pole z i The solution is as follows:

[0109] 1) The matrix bundle is given by:

[0110] [Y2]-λ[Y1]=[Z1][B]{[Z0]-λ[I]}[Z2] (6);

[0111] Where [B] represents the residual matrix; [I] represents the n×n identity matrix;

[0112] [Z0] = diag[z1,z2,…,z n ] (7);

[0113]

[0114] n represents the desired number of eigenvalues; L represents the beam parameter;

[0115] 2) The steps of the matrix bundle method are as follows:

[0116] ①Select L, satisfying n≤L≤Nn;

[0117] ②Construct matrix [Y]:

[0118]

[0119] The matrix [Y] is composed of the observation values ​​at each sampling time;

[0120] ③Perform singular value decomposition on [Y]:

[0121] [Y]=[U][S][V] T(11); where [U] and [V] are unitary matrices, containing [Y][Y] respectively. T and [Y] T The eigenvector of [Y];

[0122] ④Construct matrices [V1] and [V2] and satisfy:

[0123] V1=[v1 v2 v3 … v n-1 ] (12);

[0124] V2=[v2 v3 v4 … v n ] (13);

[0125] In the formula, v n-1 is the n-1th right singular eigenvector; v n is the nth right singular eigenvector;

[0126] ⑤Construct matrices [Y1] and [Y2]:

[0127] [Y1]=[V1] T [V1] (14);

[0128] [Y2]=[V2] T [V1] (15);

[0129] ⑥ The pole zi to be solved is calculated as the eigenvalue of the matrix pair {[Y1]; [Y2]}.

[0130] The eigenvalue λ i The solving unit is used to obtain the characteristic value λ of the power system oscillation information identification parameter through the following formula: i :

[0131]

[0132] The variable B containing the identification parameters i Solving unit, used to solve B i , so that it satisfies that for all k The details are as follows:

[0133] From formula (4), we can get:

[0134]

[0135] Formula (17) can be simplified as:

[0136] ZB=Y(18);

[0137] The matrix B is N×n, so equation (18) is solved by the least squares method.

[0138] The Levenberg-Marquardt algorithm module is used to implement the following functions:

[0139] 1) Minimizing the sum of squares of errors is used as the expression of nonlinear least squares problem:

[0140]

[0141] In the formula, y i Is the system at t i The output measurement value at time, x is the vector composed of the identification parameters of each oscillation information obtained by the fitting function in the matrix bundle method;

[0142] 2) Solve equation (19):

[0143] The Taylor series expansion of the function f(x) at x0 is:

[0144] min f(x0+Δx)=||f(x0)+J T Δx|| 2 (20);

[0145] Where J is the Jacobian matrix;

[0146]

[0147] Derivate the minimization function of equation (20), set the derivative to 0, and find the extreme value of the error function:

[0148] JJ T Δx=-Jf(x0) (23);

[0149] The iteration formula is:

[0150] Δx=-(J T J) -1 J T f(x0) (24);

[0151] Therefore, the unit matrix is ​​introduced to ensure the continued feasibility of iterative operations:

[0152] H≈J T J+μI(25);

[0153] In the formula, μ is the trust region, I is the unit matrix;

[0154] The update iteration formula is:

[0155] x k+1 =x k -(J T J+μI) -1 g k (26);

[0156] In the formula, g k is the negative gradient direction, x k is the iteration value of the kth iteration;

[0157] 3) Perform iterative operations:

[0158] ①Give the initial value x0 and the initial optimization radius μ1;

[0159] ②For the kth iteration, add the trust region to the Gauss-Newton method and solve:

[0160]

[0161] ③ Calculate and obtain the quantitative index ρ to measure the degree of approximation:

[0162]

[0163] ④. Then set the next optimization radius μ1′=2μ1, and let x k+1 =x k +Δx k Iterate

[0164] Then set the next optimization radius μ1′=0.5μ1, let x k+1 =x k +Δx k Iterate; (If equal is missing, please change one of them to be greater than or equal to or less than or equal to)

[0165] ⑤. Determine whether the above algorithm converges. If not, return to step ②. If converged, end. The optimal solution obtained is the subsynchronous oscillation mode identification result of the power system.

[0166] A storage medium, which is a computer-readable storage medium, stores a computer program on the computer-readable storage medium, and when the computer program is executed by a processor, the subsynchronous oscillation mode identification method of a power system containing a grid-type converter based on a matrix bundle and LM is implemented.

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

[0168] The present invention first uses the matrix bundle method to identify the mode parameters, obtains the initial values ​​of each mode parameter, and then uses the LM method to iterate and solve the optimal value of each modal parameter. The LM algorithm introduces the trust region. The optimization problem is changed from an unconstrained least squares problem to a constrained least squares problem, which improves the identification accuracy. Compared with the traditional method, the method proposed by the present invention is simpler in calculation and more accurate, and can have a good development prospect in subsynchronous oscillation early warning. BRIEF DESCRIPTION OF THE DRAWINGS

[0169] The present invention is further described below in conjunction with the attached drawings, tables and specific implementations:

[0170] Figure 1 It is a simulation diagram of subsynchronous oscillation of a single-machine infinite system including a grid-forming (GFM) converter in an embodiment of the present invention;

[0171] Figure 2 is a graph of active power delivered by a grid-type converter in an embodiment of the present invention;

[0172] Figure 3 It is a spectrum diagram after performing FFT transformation on the power curve in an embodiment of the present invention. DETAILED DESCRIPTION

[0173] The method for identifying subsynchronous oscillation modes of a power system containing a grid-type converter based on a matrix bundle and LM in the present invention comprises the following steps:

[0174] 1. Matrix bundle method

[0175] The matrix bundle method is a method for extracting system modal parameters based on signal poles. It can be used to extract oscillation information such as amplitude, attenuation coefficient, frequency, and phase of a given signal. The Prony method solves the characteristic roots of the system by solving the roots of the characteristic polynomial of the system, while the matrix bundle method calculates by constructing a matrix, and the eigenvalues ​​of the system can be directly solved through the eigenvalues ​​of the matrix.

[0176] The matrix bundle method constructs a fitting function to fit the function of a given signal. The fitting function can be expressed as

[0177]

[0178] The fitted signal is expressed as:

[0179] y(t k )=y(k),k=0,1,…,N-1 (2);

[0180] N is the number of samples, and the N sample values ​​are obtained by sampling at equal intervals, and the sampling time interval is Δt.

[0181] Formula (1) can be rewritten into the following complex exponential form:

[0182]

[0183] Convert the above formula into

[0184]

[0185] In the formula,

[0186]

[0187] The following describes how to solve z i method.

[0188] The matrix bundle is given by:

[0189] [Y2]-λ[Y1]=[Z1][B]{[Z0]-λ[I]}[Z2] (6); where,

[0190] [Z0] = diag[z1,z2,…,z n ] (7);

[0191]

[0192] [B] represents the residual matrix; [I] represents the n×n identity matrix; n represents the desired number of eigenvalues; L represents the bundle parameter, satisfying n≤L≤Nn.

[0193] Matrix bundle method steps:

[0194] 1). Select L such that n≤L≤Nn.

[0195] 2). Construct matrix [Y].

[0196]

[0197] Here, the matrix [Y] is a matrix composed of the observation values ​​at each sampling time.

[0198] 3). Perform singular value decomposition on [Y]:

[0199] [Y]=[U][S][V] T (11); where [U] and [V] are unitary matrices, respectively containing [Y][Y] T and [Y] T The feature vector of [Y].

[0200] 4). Construct matrices [V1] and [V2] to satisfy:

[0201] V1=[v1 v2 v3 … v n-1 ] (12);

[0202] V2=[v2 v3 v4 … v n ] (13);

[0203] In the formula, v i is the i-th right singular eigenvector of V.

[0204] 5) Construct [Y1] and [Y2]:

[0205] [Y1]=[V1] T [V1] (14);

[0206] [Y2]=[V2] T [V1] (15);

[0207] 6). The required pole zi can be found as the eigenvalue of the matrix pair {[Y1]; [Y2]}.

[0208] Eigenvalue λ of the oscillation information identification parameter i It can be obtained by the following formula:

[0209]

[0210] Once z is solved i , the characteristic value λ of the oscillation information identification parameter i It can be solved according to formula (16). The next step is to solve B i , so that it satisfies that for all k

[0211] From formula (4), we can get:

[0212]

[0213] The above formula can be simplified as:

[0214] ZB=Y(18);

[0215] The matrix B is N×n, so the above equation must also be solved by the least squares method. Then, according to equation (3), the estimated waveform can be solved:

[0216] 2. LM method

[0217] LM algorithm is one of the optimization algorithms. Optimization is to find the parameter vector that minimizes the function value. It has a wide range of applications, such as economics, management optimization, network analysis, optimal design, mechanical or electronic design, etc.

[0218] The LM algorithm introduces a trust region, which changes the optimization problem from an unconstrained least squares problem to a constrained least squares problem. It can be simply understood as: the gradient descent method is used in the early iteration, and the Gauss-Newton method is used in the later iteration.

[0219] The nonlinear least squares problem is formulated by minimizing the sum of squared errors:

[0220] min f(x0+Δx)=||f(x0+Δx)|| 2(19);

[0221] In the formula, y i Is the system at t i The output measurement value at time , x is a vector composed of the identification parameters such as amplitude, attenuation coefficient, frequency and phase in formula (1), which depends on the eigenvalue of the system state matrix.

[0222] To solve equation (19), we can use the same steps as the Newton-Raphson method. The function f(x) is expanded in a Taylor series around a certain x0:

[0223] min f(x0+Δx)=||f(x0)+J T Δx|| 2 (20);

[0224] Where J is the Jacobian matrix, Δx is the small perturbation of the function f(x) at x0;

[0225]

[0226] Derivative the minimization function of equation (20), set the derivative to 0, and find the extreme value of the error function:

[0227]

[0228] Solving for this yields:

[0229] JJ T Δx=-Jf(x0) (23);

[0230] The iteration formula is:

[0231] Δx=-(J T J) -1 J T f(x0) (24);

[0232] J T J is not necessarily reversible, so we need to introduce the unit matrix:

[0233] H≈J T J+μI(25);

[0234] The update iteration formula is:

[0235] x k+1 =x k -(J T J+μI) -1 g k (26);

[0236] In the formula, g k is the negative gradient direction, x kis the iteration value of the kth iteration;

[0237] The algorithm flow is:

[0238] ① Given the initial value x0 and the initial optimization radius μ1:

[0239] ②For the kth iteration, add the trust region to the Gauss-Newton method and solve:

[0240]

[0241] In formula (27), f(x k ) is the fitting function value of the kth iteration; Δx k is the step length;

[0242] ③ Calculate and obtain the quantitative index ρ to measure the degree of approximation:

[0243]

[0244] ④. Then set the next optimization radius μ1′=2μ1, and let x k+1 =x k +Δx k Iterate

[0245] Then set the next optimization radius μ1′=0.5μ1, let x k+1 =x k +Δx k Iterate

[0246] ⑤. Determine whether the above algorithm converges. If not, return to step ②. If converged, end. The optimal solution obtained is the subsynchronous oscillation mode identification result of the power system.

[0247] Like all nonlinear iterative methods, whether the LM method can converge depends on the choice of initial values. Choosing a suitable initial value is beneficial to the iterative process.

[0248] Example:

[0249] In order to verify the accuracy of the matrix bundle-LM method in the identification of subsynchronous oscillation parameters of power systems with grid-type converters, the present invention builds a subsynchronous oscillation simulation diagram of a single-machine infinite system with a grid-type converter in Simulink, as shown in Figure 1 As shown in the figure. In the example, the grid-type converter adopts the classic droop control. At the second second, the short-circuit ratio is changed by changing the line impedance, which stimulates the subsynchronous oscillation phenomenon. The active power curve sent by the grid-type converter is shown in Figure 2The points on the active power curve of the grid-connected converter are taken to the MATLAB workspace at equal intervals, and the matrix bundle-LM program is run for offline identification. The identification results of each mode (including amplitude, attenuation coefficient, frequency, and phase) are shown in Table 1.

[0250] Table 1 Identification results obtained by matrix bundle-LM method

[0251]

[0252] In order to verify the feasibility of the matrix bundle-LM method, the active power curve sent by the grid-type converter in the simulation example is transformed by fast Fourier transform FFT to obtain the spectrum diagram of the power signal as shown in the figure. Figure 3 As shown. Observing the spectrum, it is found that the oscillation frequency of the system is around 19 Hz, and the data identified by the matrix beam-LM method also contains a component with a frequency of 19 Hz. Therefore, the resonant frequencies obtained by the two methods are similar, which also verifies that this method has certain advantages and is simple to calculate, so it can have a good development prospect in subsynchronous oscillation early warning.

Claims

1. A method for identifying subsynchronous oscillation modes of a power system with a grid-connected converter based on matrix bundles and LM is characterized by: The method comprises the following steps, and the following steps are performed in sequence: Step 1: extracting eigenvalues ​​of oscillation information identification parameters of a given signal by using a matrix bundle method, and taking a vector composed of the extracted eigenvalues ​​as an initial value, wherein the oscillation information identification parameters include amplitude, attenuation coefficient, frequency and phase; Step 2: Use the Levenberg-Marquardt algorithm to iteratively optimize the vector composed of each eigenvalue starting from the initial value, and the optimal solution obtained is the subsynchronous oscillation mode identification result of the power system.

2. The method for identifying subsynchronous oscillation modes of a power system containing a grid-type converter based on matrix bundle and LM according to claim 1 is characterized by: The specific method of extracting the characteristic value of the oscillation information identification parameter of a given signal by using the matrix bundle method is: S101. The matrix bundle method fits the function of a given signal by constructing a fitting function, which is expressed as: Where: A i represents the amplitude, σ i represents the attenuation factor, ω i Indicates frequency, Indicates phase; The function of the given signal is expressed as: y(t k )=y(k),k=0,1,…,N-1 (2); In the formula, y(k) represents the kth sampling value among the sampling values ​​with equal time intervals; N is the total number of sampling values; Rewrite equation (1) into the following complex exponential form: Where: λ i The characteristic value representing the identification parameter of the oscillation information; Simplify formula (3) and express it as: In formula (4), S102, perform extreme point z i Solution: S103, characteristic value λ of oscillation information identification parameter i It is obtained by the following formula: S104, Solve B i , so that it satisfies that for all k From formula (4), we can get: Formula (17) can be simplified as: ZB=Y(18); The matrix B is N × n, so equation (18) is solved by the least squares method; S105, solve according to formula (3) and obtain the fitted That is, the extracted eigenvalues.

3. The method for identifying subsynchronous oscillation modes of a power system containing a grid-type converter based on matrix bundle and LM according to claim 2 is characterized by: The pole z i The solution process is as follows: 1) The matrix bundle is given by: [Y2]-λ[Y1]=[Z1][B]{[Z0]-λ[I]}[Z2] (6); Where [B] represents the residual matrix; [I] represents the n×n identity matrix; [Z0]=diag[z1,z2,…,z n ] (7); Where n represents the desired number of eigenvalues; L represents the beam parameter; 2) The steps of the matrix bundle method are as follows: ①Select L, satisfying n≤L≤Nn; ②Construct matrix [Y]: Wherein, the matrix [Y] is the matrix composed of the observation values ​​at each sampling time; ③Perform singular value decomposition on [Y]: [Y]=[U][S][V] T (11); Where [U] and [V] are unitary matrices, containing [Y][Y] respectively. T and [Y] T The eigenvector of [Y]; ④Construct matrices [V1] and [V2] and satisfy: V1=[v1 v2 v3 ... v n-1 ] (12); V2=[v2 v3 v4 … v n ] (13); In the formula, v n-1 is the n-1th right singular eigenvector; v n is the nth right singular eigenvector; ⑤Construct matrices [Y1] and [Y2]: [Y1]=[V1] T [V1] (14); <h2 style=";text-align:left;direction:ltr">[Y2]=[V2]<h2 style=";text-align:left;direction:ltr"> T <h2 style=";text-align:left;direction:ltr"> [V1] (15); ⑥ The pole zi to be solved is calculated as the eigenvalue of the matrix pair {[Y1]; [Y2]}.

4. The method for identifying subsynchronous oscillation modes of a power system containing a grid-type converter based on matrix bundle and LM according to claim 1 is characterized by: The Levenberg-Marquardt algorithm introduces a trust region, which changes the optimization problem from an unconstrained least squares problem to a constrained least squares problem, which is a nonlinear least squares problem. The specific steps are as follows: 1) Minimizing the sum of squares of errors is used as the expression of nonlinear least squares problem: In the formula, y i Is the system at t i The output measurement value at time, x is the vector composed of the identification parameters of each oscillation information obtained by the fitting function in the matrix bundle method; 2) Solve equation (19): The Taylor series expansion of the function f(x) at x0 is: min f(x0+Δx)=||f(x0)+J T Δx|| 2 (20); Where J is the Jacobian matrix, Δx is the small perturbation of the function f(x) at x0; In formula (21), x n ,y n are the identification parameter value and output measurement value at the nth sampling time respectively; Derivate the minimization function of equation (20), set the derivative to 0, and find the extreme value of the error function: JJ T Δx=-Jf(x0) (23); The iteration formula is: Δx=-(J T J) -1 J T f(x0) (24); Therefore, the matrix H is introduced to ensure the continued feasibility of iterative operations: H≈J T J+μI(25); In the formula, μ is the trust region, I is the unit matrix; The update iteration formula is: x k+1 =x k -(J T J+μI) -1 g k (26); In the formula, g k is the negative gradient direction, x k is the iteration value of the kth iteration; 3) Perform iterative operations: ① Given the initial value x0 and the initial optimization radius μ1: ②For the kth iteration, add the trust region to the Gauss-Newton method and solve: In formula (27), f(x k ) is the fitting function value of the kth iteration; Δx k is the step length; ③ Calculate and obtain the quantitative index ρ to measure the degree of approximation: ④. Then set the next optimization radius μ1′=2μ1, and let x k+1 =x k +Δx k Iterate Then set the next optimization radius μ1′=0.5μ1, let x k+1 =x k +Δx k Iterate ⑤. Determine whether the above algorithm converges. If not, return to step ②. If converged, end. The optimal solution obtained is the subsynchronous oscillation mode identification result of the power system.

5. A device for identifying subsynchronous oscillation modes of a power system containing a grid-type converter based on a matrix bundle and LM, using the method for identifying subsynchronous oscillation modes of a power system containing a grid-type converter based on a matrix bundle and LM according to claim 1, characterized in that: include: An eigenvalue extraction module is used to extract the eigenvalues ​​of the oscillation information identification parameters of a given signal by using a matrix bundle method, and use the vector composed of the extracted eigenvalues ​​as an initial value; The Levenberg-Marquardt algorithm module is used to iteratively optimize the vector composed of each eigenvalue starting from the initial value to obtain the optimal solution.

6. The subsynchronous oscillation mode identification device of a power system containing a grid-type converter based on a matrix bundle and LM according to claim 5 is characterized in that: The eigenvalue extraction module includes a fitting unit, a pole z i Solving for unit, eigenvalue λ i Solving unit and variable B with identification parameters i Solving unit; The fitting unit is used to implement the following functions: Construct a fitting function to fit the function of the given signal. The fitting function is expressed as: Where: A i represents the amplitude, σ i represents the attenuation factor, ω i Indicates frequency, Indicates phase; The function of the given signal is expressed as: y(t k )=y(k),k=0,1,…,N-1 (2); In the formula, y(k) represents the kth sample value among the samples with equal time intervals; N is the total number of samples; Rewrite equation (1) into the following complex exponential form: Where: λ i The characteristic value representing the identification parameter of the oscillation information; Simplify formula (3) and express it as: In formula (4), In formula (5), z i Indicates the system pole; The pole z i Solving for unit, eigenvalue λ i Solving unit and variable B with identification parameters i Substitute the solution result of the solving unit into the fitting function to obtain the fitted That is, the extracted eigenvalues.

7. The subsynchronous oscillation mode identification device of a power system containing a grid-type converter based on a matrix bundle and LM according to claim 6 is characterized in that: The pole z i Solver element for pole z i The solution is as follows: 1) The matrix bundle is given by: [Y2]-λ[Y1]=[Z1][B]{[Z0]-λ[I]}[Z2] (6); Where [B] represents the residual matrix; [I] represents the n×n identity matrix; [Z0]=diag[z1,z2,…,z n ] (7); n represents the desired number of eigenvalues; L represents the beam parameter; 2) The steps of the matrix bundle method are as follows: ①Select L, satisfying n≤L≤Nn; ②Construct matrix [Y]: The matrix [Y] is composed of the observation values ​​at each sampling time; ③Perform singular value decomposition on [Y]: [Y]=[U][S][V] T (11); where [U] and [V] are unitary matrices, containing [Y][Y] respectively. T and [Y] T The eigenvector of [Y]; ④Construct matrices [V1] and [V2] and satisfy: V1=[v1 v2 v3 …v n-1 ] (12); V2=[v2 v3 v4…v n ] (13); In the formula, v n-1 is the n-1th right singular eigenvector; v n is the nth right singular eigenvector; ⑤Construct matrices [Y1] and [Y2]: [Y1]=[V1] T [V1] (14); <h2 style=";text-align:left;direction:ltr">[Y2]=[V2]<h2 style=";text-align:left;direction:ltr"> T <h2 style=";text-align:left;direction:ltr"> [V1] (15); ⑥ The pole zi to be solved is calculated as the eigenvalue of the matrix pair {[Y1]; [Y2]}.

8. The subsynchronous oscillation mode identification device of a power system containing a grid-type converter based on a matrix bundle and LM according to claim 6 is characterized in that: The eigenvalue λ i The solving unit is used to obtain the characteristic value λ of the power system oscillation information identification parameter through the following formula: i : The variable B containing the identification parameters i Solving unit, used to solve B i , so that it satisfies that for all k The details are as follows: From formula (4), we can get: Formula (17) can be simplified as: ZB=Y(18); The matrix B is N×n, so equation (18) is solved by the least squares method.

9. The subsynchronous oscillation mode identification device of a power system containing a grid-type converter based on a matrix bundle and LM according to claim 6, characterized in that: The Levenberg-Marquardt algorithm module is used to implement the following functions: 1) Minimizing the sum of squares of errors is used as the expression of nonlinear least squares problem: In the formula, y i Is the system at t i The output measurement value at time, x is the vector composed of the identification parameters of each oscillation information obtained by the fitting function in the matrix bundle method; 2) Solve equation (19): The Taylor series expansion of the function f(x) at x0 is: min f(x0+Δx)=||f(x0)+J T Δx|| 2 (20); Where J is the Jacobian matrix; Derivate the minimization function of equation (20), set the derivative to 0, and find the extreme value of the error function: JJ T Δx=-Jf(x0) (23); The iteration formula is: Δx=-(J T J) -1 J T f(x0) (24); Therefore, the unit matrix is ​​introduced to ensure the continued feasibility of iterative operations: H≈J T J+μI(25); In the formula, μ is the trust region, I is the unit matrix; The update iteration formula is: x k+1 =x k -(J T J+μI) -1 g k (26); In the formula, g k is the negative gradient direction, x k is the iteration value of the kth iteration; 3) Perform iterative operations: ①Give the initial value x0 and the initial optimization radius μ1; ②For the kth iteration, add the trust region to the Gauss-Newton method and solve: ③ Calculate and obtain the quantitative index ρ to measure the degree of approximation: ④. Then set the next optimization radius μ′1=2μ1, and let x k+1 =x k +Δx k Iterate Then set the next optimization radius μ1′=0.5μ1, let x k+1 =x k +Δx k Iterate ⑤. Determine whether the above algorithm converges. If not, return to step ②. If converged, end. The optimal solution obtained is the subsynchronous oscillation mode identification result of the power system.

10. A storage medium, characterized in that: The storage medium is a computer-readable storage medium on which a computer program is stored. When the computer program is executed by a processor, the method for identifying subsynchronous oscillation modes of a power system containing a grid-type converter based on matrix bundles and LM as described in any one of claims 1 to 4 is implemented.

Images