METHOD FOR ANALYZING STABILITY OF PERMANENT MAGNET SYNCHRONOUS GENERATOR-BASED WIND TURBINES CONNECTED TO WEAK POWER GRID CONSIDERING INFLUENCE OF POWER CONTROL

Abstract

Provided is a method for analyzing the stability of a PMSG-WT connected to a weak power grid considering the influence of power control. New energy power generation mostly uses a perturbation and observation (P&O) method for maximum power point tracking, and nonlinear discontinuous links therein make stability analysis difficult. The present application analyzes the stability of the PMSG-WT connected to the weak grid system based on a describing function method, and fully considers the nonlinear discontinuous links in the power loop, thus making the analysis result more accurate. At the same time, the describing function method is a method that can quantitatively calculate the frequency and amplitude of oscillation. The analysis method of the present application can provide a powerful and good reference for oscillation suppression and controller design.

Claims

1. A method for analyzing stability of a PMSG-WT connected to a weak power grid considering influence of power control, wherein the PMSG-WT comprises a wind turbine, a generator, an MSC, a DC capacitor, a GSC, a filter, an MSC controller, and a GSC controller; the MSC controller comprises a power loop, a rotating speed loop, and a machine-side current loop; the GSC controller comprises a voltage loop and a grid-side current loop; and the method comprises the following steps: step 1: acquiring main parameters of the PMSG-WT, establishing mathematical models of the wind turbine, the generator and the MSC, the MSC controller, the DC capacitor, the GSC and the filter, and the GSC controller, respectively, performing linearization in a dq coordinate system, and calculating steady-state operation parameters to obtain small-signal models of the wind turbine, the generator and the MSC, the MSC controller, the DC capacitor, the GSC and the filter, and the GSC controller; step 2: modeling the power loop in the MSC controller based on a describing function method, wherein its mathematical expression is: ${\omega }_g^{\mathrm{ref}}=\frac{\epsilon }{T_P}\int \mathrm{sgn}\left(P^{\mathrm{ref}}-P_n\right)\mathrm{sgn}\left(P_n-P_{n-1}\right)\mathrm{sgn}\left({\omega }_{g,n}-{\omega }_{g,n-1}\right)\mathrm{dt}$ where ε represents a perturbation step length of the power loop, T.sub.p represents a control period of the power loop, P.sup.ref represents a reference value of the output power of the PMSG-WT, P represents an output power of the PMSG-WT, P.sub.n represents an output power of the PMSG-WT at a current sampling moment n, P.sub.n−1 represents an output power of the PMSG-WT at a previous sampling moment, ω.sub.g represents a rotating speed of the generator, ω.sub.g.sup.ref represents a reference value of the rotating speed of the generator, ω.sub.g,n represents a rotating speed of the generator at the current sampling moment, and ω.sub.g,n−1 represents a rotating speed of the generator at the previous sampling moment; sgn(x) is a sign function, when x≥0, sgn(x)=1, and when x<0, sgn(x)=−1; considering an actual power-rotating speed curve of the PMSG-WT, then: $\mathrm{sgn}\left(P_n-P_{n-1}\right)\mathrm{sgn}\left({\omega }_{g,n}-{\omega }_{g,n-1}\right)=\{\begin{array}{cc}1& \left({\omega }_g\le {\omega }_{\mathrm{mpp}}\right)\\ -1& \left({\omega }_g>{\omega }_{\mathrm{mpp}}\right)\end{array}$ where ω.sub.mpp represents a rotating speed of the generator at a maximum power point; and then, the power loop model can be simplified as: ${\omega }_g^{\mathrm{ref}}=\{\begin{array}{cc}\frac{\epsilon }{T_P}\int \mathrm{sgn}\left(P^{\mathrm{ref}}-P_n\right)\mathrm{dt},& \left({\omega }_g\le {\omega }_{\mathrm{mpp}}\right)\\ -\frac{\epsilon }{T_P}\int \mathrm{sgn}\left(P^{\mathrm{ref}}-P_n\right)\mathrm{dt},& \left({\omega }_g>{\omega }_{\mathrm{mpp}}\right)\end{array}$ wherein the sign function in the formula can be modeled by using a describing function, and the describing function thereof is: $N\left(A\right)=\frac{4}{\pi A}$ in the formula, A represents an amplitude of an input signal; step 3: considering the influence of the weak AC power grid, combining the linear parts of the small-signal models of the weak power grid and the power loop with the small-signal models established in step 1, and deriving a transfer function G(s) of the linear part of a system; and step 4: drawing G(s) and −1/N(A) curves in a complex plane, and analyzing the stability of the system based on the describing function method, wherein the method is specifically that if G(s) contains a pole of the right half plane, the system must be unstable; if the G(s) does not contain the pole of the right half plane, then judging the stability of the system through the relationship between a G(s) curve and a −1/N(A) curve: a. if the G(s) curve does not surround the −1/N(A) curve, the system is stable and does not oscillate; b. if the G(s) curve and the −1/N(A) curve intersect, the system is critically stable, at this time, the system oscillates at a constant amplitude and a constant frequency, and the frequency and the amplitude of the oscillation can be calculated by the following formula: $\{\begin{array}{c}{G}_{\mathrm{Im}}\left({\omega }_0\right)=0\\ N\left(A_0\right)=-1/G_{\mathrm{Re}}\left({\omega }_0\right)\end{array}$ where G(jω)=G.sub.Re(ω)+jG.sub.Im(ω), G.sub.Re means to find the real part of a complex number, G.sub.Im means to find the imaginary part of the complex number, ω.sub.0 represents an angular frequency of the oscillation, and A.sub.0 represents the amplitude of the oscillation; and c. if the G(s) curve surrounds the −1/N(A) curve, the system is unstable.

2. The method according to claim 1, wherein in step 1, the small-signal models of the wind turbine, the generator and the MSC, and the MSC controller are established as follows: the mathematical model of the wind turbine is established as:
sJω.sub.g=T.sub.M−T.sub.e−Bω.sub.g in the formula, J represents rotational inertia of an equivalent concentrated mass block of the wind turbine and the generator, T.sub.m represents a mechanical torque of the generator, T.sub.e represents an electromagnetic torque of the generator, B represents a self-damping coefficient, it is deemed that B=0 herein, and s represents a parameter introduced by the Laplace transform; linearization is performed on the model to obtain:
sJΔω.sub.g=ΔT.sub.m−ΔT.sub.e the electromagnetic torque of the generator is:
T.sub.e= 3/2n.sub.pψ.sub.fi.sub.qr n.sub.p represents the number of the pole pairs of the generator, i.sub.qr represents stator current of the q-axis generator, and ψ.sub.f represents a permanent magnet flux linkage of the generator, and this formula is linearized to obtain:
ΔT.sub.e= 3/2n.sub.pψ.sub.fΔi.sub.qr the mechanical torque of the generator is:
T.sub.m=B.sub.tω.sub.g in the formula, $B_t=\frac{{\mathrm{dT}}_m}{d{\omega }_g}{.Math.}_{{\omega }_g={\omega }_g^{*}}$  represents a linearization constant of the mechanical torque of the wind turbine, and ω.sub.g* represents a steady-state value of the rotating speed of the generator, which is a specific constant at an equilibrium point; and this formula is linearized to obtain:
ΔT.sub.m=B.sub.tω.sub.g from this, the small-signal model of the wind turbine can be obtained, ${\mathrm{\Delta \omega }}_g=\frac{3n_p{\psi }_f}{2\left(B_t-\mathrm{sJ}\right)}\Delta i_{\mathrm{qr}}=\left[\begin{array}{cc}0& \frac{3n_p{\psi }_f}{2\left(B_t-\mathrm{sJ}\right)}\end{array}\right].Math.\Delta i_{\mathrm{dqr}}$ $\mathrm{setting}$ $G_{\mathrm{iq}\omega }=\left[\begin{array}{cc}0& \frac{3n_p{\psi }_f}{2\left(B_t-\mathrm{sJ}\right)}\end{array}\right],\mathrm{then}{\mathrm{\Delta \omega }}_g=G_{\mathrm{iq}\omega }.Math.\Delta i_{\mathrm{dqr}},{\mathrm{\Delta \omega }}_e=n_p{G}_{\mathrm{iq}\omega }.Math.\Delta i_{\mathrm{dqr}};$ the mathematical model of the generator and the MSC is established as: $\left[\begin{array}{cc}{\mathrm{sL}}_s+R_s& -{\omega }_e{L}_s\\ {\omega }_e{L}_s& {\mathrm{sL}}_s+R_s\end{array}\right].Math.\left[\begin{array}{c}{i}_{\mathrm{dr}}\\ i_{\mathrm{qr}}\end{array}\right]=\left[\begin{array}{c}0\\ {\omega }_e{\psi }_f\end{array}\right]-\left[\begin{array}{c}{d}_{\mathrm{dr}}\\ d_{\mathrm{qr}}\end{array}\right].Math.u_{\mathrm{dc}}$ in the formula, R.sub.s and L.sub.s represent rotor resistance and armature inductance of the generator, respectively, ω.sub.e represents an electrical angular speed of the rotor, and ω.sub.e=n.sub.pω.sub.g.; and the model is linearized to obtain the small-signal model of the generator and MSC:
Z.sub.dqr.Math.Δi.sub.dqr=−U*.sub.dc.Math.Δd.sub.dqr−D*.sub.dqr.Math.Δu.sub.dc+G.sub.ωe.Math.Δω.sub.e where $Z_{\mathrm{dqr}}=\left[\begin{array}{cc}{\mathrm{sL}}_s+R_s& -{\omega }_e^{*}{L}_s\\ {\omega }_e^{*}{L}_s& {\mathrm{sL}}_s+R_s\end{array}\right],\Delta i_{\mathrm{dqr}}=\left[\begin{array}{c}\Delta i_{\mathrm{dr}}\\ \Delta i_{\mathrm{qr}}\end{array}\right],\Delta d_{\mathrm{dqr}}=\left[\begin{array}{c}\Delta d_{\mathrm{dr}}\\ \Delta d_{\mathrm{qr}}\end{array}\right],D_{\mathrm{dqr}}^{*}=\left[\begin{array}{c}{D}_{\mathrm{dr}}^{*}\\ D_{\mathrm{qr}}^{*}\end{array}\right],G_{\omega e}=\left[\begin{array}{c}{L}_s{I}_{\mathrm{qr}}^{*}\\ {\psi }_f-L_s{I}_{\mathrm{dr}}^{*}\end{array}\right];$  capital letters and superscripts * represent steady-state components of corresponding lowercase variables, and Δ represents a small-signal component of the corresponding variable; the mathematical model of the MSC controller is established as: $\{\begin{array}{c}{d}_{\mathrm{dr}}^c{U}_{\mathrm{dc}}^{*}=-H_{\mathrm{cr}}\left(0-i_{\mathrm{dr}}^c\right)+{\omega }_e{L}_s{i}_{\mathrm{qr}}^c\\ d_{\mathrm{qr}}^c{U}_{\mathrm{dc}}^{*}=-H_{\mathrm{cr}}\left(-H_{\omega }\left({\omega }_g^{\mathrm{ref}}-{\omega }_g\right)-i_{\mathrm{qr}}^c\right)-{\omega }_e{L}_s{i}_{\mathrm{dr}}^c+{\omega }_e{\psi }_f\end{array}$ where $H_{\mathrm{cr}}=K_{\mathrm{cpr}}+\frac{K_{\mathrm{cir}}}{s},K_{\mathrm{cpr}}$  and K.sub.cir represent a proportional parameter and an integral parameter of machine-side current loop PI control, respectively, $H_{\omega }=K_{\omega p}+\frac{K_{\omega i}}{s},$  and K.sub.ωp and K.sub.ωi represent the proportional parameter and the integral parameter of rotating speed loop PI control, respectively; superscript c represents a dq coordinate system of the MSC controller; ω.sub.g.sup.ref represents the reference value of the rotating speed of the generator; and by linearizing the model, the small-signal model of the MSC controller can be obtained as: $\Delta d_{\mathrm{dqr}}^c=\frac{1}{U_{\mathrm{dc}}^{*}}\left(G_{\mathrm{cr}}.Math.\Delta i_{\mathrm{dqr}}^c+G_{\omega g}\left({\omega }_g^{\mathrm{ref}}-{\omega }_g\right)+G_{\omega e}.Math.{\mathrm{\Delta \omega }}_e\right)$ $\mathrm{where}$ $G_{\mathrm{cr}}=\left[\begin{array}{cc}{H}_{\mathrm{cr}}& {\omega }_e^{*}{L}_s\\ -{\omega }_e^{*}{L}_s& H_{\mathrm{cr}}\end{array}\right],G_{\omega g}=\left[\begin{array}{c}0\\ H_{\mathrm{cr}}{H}_{\omega }\end{array}\right],$ ω*.sub.e represents the steady-state value of the electrical angular speed of the rotor; affected by the disturbance of the rotating speed, the dq coordinate system of the MSC controller will have a phase angle difference with the dq coordinate system of the generator during the dynamic process; the electrical angle of the generator is:
θ.sub.e=n.sub.pω.sub.g/s the formula is linearized to obtain: ${\mathrm{\Delta \theta }}_e=\frac{n_p{\mathrm{\Delta \omega }}_g}{s}=\frac{3n_p^2{\psi }_f}{2s\left(B_t-\mathrm{sJ}\right)}\Delta i_{\mathrm{qr}}=H_{\theta e}\Delta i_{\mathrm{qr}}$ in the formula, $H_{\theta e}=\frac{3n_p^2{\psi }_f}{2s\left(B_t-\mathrm{sJ}\right)};$ therefore, the variable conversion relationship between the dq coordinate system of the MSC controller and the dq coordinate system of the generator is: $\{\begin{array}{c}\Delta i_{\mathrm{dqr}}^c=G_{\theta e}^i.Math.\Delta i_{\mathrm{dqr}}\\ \Delta d_{\mathrm{dqr}}^c=G_{\theta e}^d.Math.\Delta i_{\mathrm{dqr}}+\Delta d_{\mathrm{dqr}}\end{array}$ $\mathrm{where}$ $G_{\theta e}^i=\left[\begin{array}{cc}1& H_{\theta e}{I}_{\mathrm{qr}}^{*}\\ 0& 1-H_{\theta e}{I}_{\mathrm{dr}}^{*}\end{array}\right],G_{\theta e}^d=\left[\begin{array}{cc}0& H_{\theta e}{D}_{\mathrm{qr}}^{*}\\ 0& -H_{\theta e}{D}_{\mathrm{dr}}^{*}\end{array}\right];$ then, the small-signal model of the MSC controller is: $G_{\theta e}^d.Math.\Delta i_{\mathrm{dqr}}+\Delta d_{\mathrm{dqr}}=\frac{1}{U_{\mathrm{dc}}^{*}}\left(G_{\mathrm{cr}}.Math.G_{\theta e}^i.Math.\Delta i_{\mathrm{dqr}}+G_{\omega g}\left({\omega }_g^{\mathrm{ref}}-{\omega }_g\right)+G_{\omega e}.Math.{\mathrm{\Delta \omega }}_e\right).$

3. The method according to claim 2, wherein in step 1, the establishment process of the small-signal model of the DC capacitor is as follows: the mathematical model of the DC capacitor is established as:
sC.sub.dcu.sub.dc=i.sub.dc2−i.sub.dc1=1.5(d.sub.dri.sub.dr+d.sub.qri.sub.qr)−1.5(d.sub.dgi.sub.dg+d.sub.qgi.sub.qg) in the formula, C.sub.dc represents the DC capacitor, i.sub.dc1 represents grid-side DC current, i.sub.dc2 represents machine-side DC current, and by linearizing the model, the small-signal model of the DC capacitor can be obtained as:
sC.sub.dcΔu.sub.dc=1.5(D*.sub.dqr.sup.T.Math.Δi.sub.dqr+I*.sub.dqr.sup.T.Math.Δd.sub.dqr)−1.5(D*.sub.dqg.sup.T.Math.Δi.sub.dqg+I*.sub.dqg.sup.T.Math.Δd.sub.dqg) in the formula, $I_{\mathrm{dqr}}^{*}=\left[\begin{array}{c}{I}_{\mathrm{dr}}^{*}\\ I_{\mathrm{qr}}^{*}\end{array}\right],\Delta i_{\mathrm{dqg}}=\left[\begin{array}{c}\Delta i_{\mathrm{dg}}\\ \Delta i_{\mathrm{qg}}\end{array}\right],\Delta d_{\mathrm{dqg}}=\left[\begin{array}{c}\Delta d_{\mathrm{dg}}\\ \Delta d_{\mathrm{qg}}\end{array}\right],D_{\mathrm{dqg}}^{*}=\left[\begin{array}{c}{D}_{\mathrm{dg}}^{*}\\ D_{\mathrm{qg}}^{*}\end{array}\right],I_{\mathrm{dqg}}^{*}=\left[\begin{array}{c}{I}_{\mathrm{dg}}^{*}\\ I_{\mathrm{qg}}^{*}\end{array}\right];$  and the capital letters and the superscripts * represent the steady-state components of corresponding lowercase variables, and Δ represents the small-signal component of the corresponding variable.

4. The method according to claim 3, wherein in step 1, the establishment process of the small-signal models of the GSC and the filter, and the GSC controller is as follows: the mathematical model of the GSC and the filter is established as: $\left[\begin{array}{cc}{\mathrm{sL}}_f& -{\omega L}_f\\ {\omega L}_f& {\mathrm{sL}}_f\end{array}\right].Math.\left[\begin{array}{c}{i}_{\mathrm{dg}}\\ i_{\mathrm{qg}}\end{array}\right]=\left[\begin{array}{c}{d}_{\mathrm{dg}}\\ d_{\mathrm{qg}}\end{array}\right].Math.u_{\mathrm{dc}}-\left[\begin{array}{c}{u}_{\mathrm{dg}}\\ u_{\mathrm{qg}}\end{array}\right]$ in the formula, L.sub.f represents a filtering inductance, ω represents a power frequency angular frequency, ω=100π rad/s, and u.sub.dg and u.sub.qg represent a d-axis voltage and a q-axis voltage of a grid-connected point, respectively; by linearizing the model, the small-signal model of the GSC and the filter can be obtained as:
Z.sub.f.Math.Δi.sub.dqg=U*.sub.dc.Math.Δd.sub.dqg+D*.sub.dqg.Math.Δu.sub.dc−Δu.sub.dqg where $Z_f=\left[\begin{array}{cc}{\mathrm{sL}}_f& -\omega L_f\\ \omega L_f& {\mathrm{sL}}_f\end{array}\right],\Delta i_{\mathrm{dqg}}=\left[\begin{array}{c}\Delta i_{\mathrm{dg}}\\ \Delta i_{\mathrm{qg}}\end{array}\right],\Delta d_{\mathrm{dqg}}=\left[\begin{array}{c}\Delta d_{\mathrm{dg}}\\ \Delta d_{\mathrm{qg}}\end{array}\right],D_{\mathrm{dqg}}^{*}=\left[\begin{array}{c}{D}_{\mathrm{dg}}^{*}\\ D_{\mathrm{qg}}^{*}\end{array}\right],\Delta u_{\mathrm{dqg}}=\left[\begin{array}{c}\Delta u_{\mathrm{dg}}\\ \Delta u_{\mathrm{qg}}\end{array}\right];$  the capital letters and the superscripts * represent the steady-state components of corresponding lowercase variables, and Δ represents the small-signal component of the corresponding variable; the mathematical model of the GSC controller is established as: $\{\begin{array}{c}{d}_{\mathrm{dg}}^c{U}_{\mathrm{dc}}^{*}=u_{\mathrm{dg}}^c+H_{\mathrm{cg}}\left(-H_v\left(U_{\mathrm{dcref}}-u_{\mathrm{dc}}\right)-i_{\mathrm{dg}}^c\right)-\omega L_f{i}_{\mathrm{qg}}^c\\ d_{\mathrm{qg}}^c{U}_{\mathrm{dc}}^{*}=u_{\mathrm{qg}}^c+H_{\mathrm{cg}}\left(0-i_{\mathrm{qg}}^c\right)+\omega L_f{i}_{\mathrm{dg}}^c\end{array}$ where $H_{\mathrm{cg}}=K_{\mathrm{cpg}}+\frac{K_{\mathrm{cig}}}{s},$  K.sub.cpg and K.sub.cig represent a proportional parameter and an integral parameter of grid-side current loop PI control, respectively, $H_v=K_{\mathrm{vp}}+\frac{K_{\mathrm{vi}}}{s},$  K.sub.vp and K.sub.vi represent a proportional parameter and an integral parameter of grid-side voltage loop PI control, respectively, and U.sub.dcref represents a reference value of the DC voltage; in the GSC controller, a phase-locked loop is used to keep the fan synchronous with the power grid; the superscript c represents the dq coordinate system of the GSC controller; by linearizing the model, the small-signal model of the GSC controller can be obtained as:
Δd.sub.dqg.sup.c.Math.U*.sub.dc=Δu.sub.dqg.sup.c+G.sub.cgΔi.sub.dqg.sup.c+G.sub.udcΔu.sub.dc where $G_{\mathrm{cg}}=\left[\begin{array}{cc}-H_{\mathrm{cg}}& -\omega L_f\\ \omega L_f& -H_{\mathrm{cg}}\end{array}\right],G_{\mathrm{udc}}=\left[\begin{array}{c}{H}_{\mathrm{cg}}{H}_v\\ 0\end{array}\right];$ in addition, the dynamics of the phase-locked loop should also be considered in the GSC, and its mathematical model is: $\theta =H_{\mathrm{pll}}.Math.\frac{1}{s}.Math.u_{\mathrm{qg}}^c$ where $H_{\mathrm{pll}}=K_{\mathrm{ppll}}+\frac{K_{\mathrm{ipll}}}{s},$  K.sub.ppll and K.sub.ipll represent a proportional parameter and an integral parameter of a phase-locked loop PI controller, respectively, and u.sub.qg.sup.c represents the q-axis voltage of the grid-connected point in the dq coordinate system of the GSC controller; the model is linearized to obtain: $\Delta \theta =H_{\mathrm{pll}}.Math.\frac{1}{s}.Math.\Delta u_{\mathrm{qg}}^c$ wherein there is a certain deviation between the dq coordinate system of the system and the dq coordinate system of the controller, and the two coordinate systems can be converted to each other through the following equation: $\{\begin{array}{c}\Delta x_d^c=\Delta x_d+\mathrm{\Delta \theta }.Math.X_q^{*}\\ \Delta x_q^c=\Delta x_q+\mathrm{\Delta \theta }.Math.X_d^{*}\end{array}$ in the formula, the variables Δx.sub.d and Δx.sub.q may represent the output current Δi.sub.dg and Δi.sub.qg of the GSC, the voltages Δu.sub.dg and Δu.sub.qg of the grid-connected point, or the output duty ratios Δd.sub.dg and Δd.sub.qg of the grid-side controller, and X*.sub.q and X*.sub.d represent steady-state components I*.sub.qg, I*.sub.dg, U*.sub.qg, U*.sub.dg, D*.sub.qg, D*.sub.dg; from this, the small-signal model of the phase-locked loop can be derived, namely:
Δθ=G.sub.pll.Math.Δu.sub.qg where $G_{\mathrm{pll}}=\frac{H_{\mathrm{pll}}}{s+U_{\mathrm{dg}}^{*}{H}_{\mathrm{pll}}},$  and then the relationship between the dq coordinate system of the controller and the dq coordinate system of the system can be obtained as: $\{\begin{array}{c}\Delta u_{\mathrm{dqg}}^c=G_{\mathrm{pll}}^u.Math.\Delta u_{\mathrm{dqg}}\\ \Delta i_{\mathrm{dqg}}^c=G_{\mathrm{pll}}^i.Math.\Delta u_{\mathrm{dqg}}+\Delta i_{\mathrm{dqg}}\\ \Delta d_{\mathrm{dqg}}^c=G_{\mathrm{pll}}^d.Math.\Delta u_{\mathrm{dqg}}+\Delta d_{\mathrm{dqg}}\end{array}$ $\mathrm{where}{G}_{\mathrm{pll}}^u=\left[\begin{array}{cc}1& G_{\mathrm{pll}}{U}_{\mathrm{qg}}^{*}\\ 0& 1-G_{\mathrm{pll}}{U}_{\mathrm{dg}}^{*}\end{array}\right],G_{\mathrm{pll}}^i=\left[\begin{array}{cc}0& G_{\mathrm{pll}}{I}_{\mathrm{qg}}^{*}\\ 0& -G_{\mathrm{pll}}{I}_{\mathrm{dg}}^{*}\end{array}\right],G_{\mathrm{pll}}^d=\left[\begin{array}{cc}0& G_{\mathrm{pll}}{D}_{\mathrm{qg}}^{*}\\ 0& -G_{\mathrm{pll}}{D}_{\mathrm{dg}}^{*}\end{array}\right];$  and then the small-signal model of the GSC controller is converted into:
(G.sub.pll.sup.d.Math.Δu.sub.dqg+Δd.sub.dqg).Math.U*.sub.dc=G.sub.pll.sup.u.Math.Δu.sub.dqg+G.sub.cg(G.sub.pll.sup.i.Math.Δu.sub.dqg+Δi.sub.dqg)+G.sub.udcΔu.sub.dc.

5. The method according to claim 4, wherein the step 3 specifically comprises: the AC weak power grid is expressed by series equivalent inductance of an ideal voltage source, and its mathematical model is established as: $\left[\begin{array}{cc}{\mathrm{sL}}_g& -\omega L_g\\ \omega L_g& {\mathrm{sL}}_g\end{array}\right].Math.\left[\begin{array}{c}{i}_{\mathrm{dg}}\\ i_{\mathrm{qg}}\end{array}\right]=\left[\begin{array}{c}{u}_{\mathrm{dg}}\\ u_{\mathrm{qg}}\end{array}\right]-\left[\begin{array}{c}{u}_{\mathrm{ds}}\\ u_{\mathrm{qs}}\end{array}\right]$ in the formula, L.sub.g represents the equivalent inductance of the weak grid, u.sub.ds and u.sub.qs represent ideal voltage source voltages of the d axis and the q axis, respectively, and i.sub.dg and i.sub.qg represent the d axis current and the q axis current of the AC port of the GSC, respectively; and the formula is linearized to obtain:
Z.sub.g.Math.Δi.sub.dqg=Δu.sub.dqg in the formula, $Z_g=\left[\begin{array}{cc}{\mathrm{sL}}_g& -\omega L_g\\ \omega L_g& {\mathrm{sL}}_g\end{array}\right];$ the output power of the PMSG-WT is:
P=1.5(i.sub.dgu.sub.dg+i.sub.qgu.sub.qg) the formula is linearized to obtain:
ΔP=1.5(I*.sub.dqg.sup.T.Math.Z.sub.g+U*.sub.dqg.sup.T).Math.Δi.sub.dqg in the formula, $U_{\mathrm{dqg}}^{*}=\left[\begin{array}{c}{U}_{\mathrm{dg}}^{*}\\ U_{\mathrm{qg}}^{*}\end{array}\right];$  by combining the small-signal models of the wind turbine, the generator and the MSC, the MSC controller, the DC capacitor, the GSC and the filter, and the GSC controller with the linear part of the power loop, the transfer function G(s) of the linear part of the system can be obtained as: $G\left(s\right)=-\frac{1.5\epsilon }{T_{P}s}\left(I_{\mathrm{dqg}}^{*T}{Z}_g+U_{\mathrm{dqg}}^{*T}\right)M_5^{-1}\left(D_{\mathrm{dqg}}^{*}+G_{\mathrm{udc}}\right)M_7^{-1}{M}_8.Math.\frac{1}{1+T_{f^S}}.\frac{1}{1+1.5T_{P}s}$ $M_8=-D_{\mathrm{dqr}}^{*T}.Math.M_1^{-1}.Math.G_{\omega g}+I_{\mathrm{dqr}}^{*T}.Math.M_3$ $M_7=\left(\frac{2s{C}_{\mathrm{dc}}}{3}+I_{\mathrm{dqg}}^{*T}{M}_6-I_{\mathrm{dqr}}^{*T}{M}_4+D_{\mathrm{dqg}}^{*T}{M}_5^{-1}\left(D_{\mathrm{dqg}}^{*}+G_{\mathrm{udc}}\right)+D_{\mathrm{dqr}}^{*T}{M}_1^{-1}{D}_{\mathrm{dqr}}^{*}\right)$ $M_6=\left(\left(\left(G_{\mathrm{cg}}{G}_{\mathrm{pll}}^i-G_{\mathrm{pll}}^d{U}_{\mathrm{dc}}^{*}+G_{\mathrm{pll}}^u\right)Z_g+G_{\mathrm{cg}}\right)M_5^{-1}\left(D_{\mathrm{dqg}}^{*}+G_{\mathrm{udc}}\right)+G_{\mathrm{udc}}\right)/U_{\mathrm{dc}}^{*}$ $M_5=-G_{\mathrm{cg}}+Z_f-(\left(G_{\mathrm{cg}}{G}_{\mathrm{pll}}^i-G_{\mathrm{pll}}^d{U}_{\mathrm{dc}}^{*}+G_{\mathrm{pll}}^u\right)Z_g+Z_g$ $M_4=-M_2{M}_1^{-1}{D}_{\mathrm{dqr}}^{*}/U_{\mathrm{dc}}^{*}$ $M_3=\left(-M_2{M}_1^{-1}{G}_{\omega g}+G_{\omega g}\right)/U_{\mathrm{dc}}^{*}$ $M_2=G_{\omega e}{G}_{iq\omega }{n}_p+G_{\mathrm{cr}}{G}_{\theta e}^i-G_{\omega g}{G}_{iq\omega }-U_{\mathrm{dc}}^{*}{G}_{\theta e}^d$ $M_1=G_{\mathrm{cr}}{G}_{\theta e}^i-G_{\omega g}{G}_{iq\omega }-U_{\mathrm{dc}}^{*}{G}_{\theta e}^d+Z_{\mathrm{dqr}}$ in the formula, T.sub.f represents the period of a power sampling filter, 1/(1+T.sub.fs) represents the delay of the power sampling filter, and 1/(1+1.5T.sub.ps) represents the delay of the controller and the PWM.

Description

BRIEF DESCRIPTION OF DRAWINGS

[0066] FIG. 1 shows a topological structure (A) of a PMSG-WT and its controller structure (B and D represent coordinate conversion, C represents a phase-locked loop, E represents an MSC controller, and F represents a GSC controller);

[0067] FIG. 2 shows a distribution diagram of G(s) poles;

[0068] FIG. 3 shows G(s) and −1/N(A) curves when an L.sub.g value is changed;

[0069] FIG. 4a shows a simulation verification waveform (a) when Lg=0.1 mH;

[0070] FIG. 4b shows FFT spectrum analysis (b) when Lg=0.1 mH;

[0071] FIG. 5a shows the simulation verification waveform (a) when Lg=0.4 mH; and

[0072] FIG. 5b shows FFT spectrum analysis (b) when Lg=0.4 mH.

DESCRIPTION OF EMBODIMENTS

[0073] The topological structure of the PMSG-WT and its controller according to the present application are shown in FIG. 1. The PMSG-WT includes a wind turbine, a generator, an MSC, a DC capacitor, a GSC, and a filter. The wind turbine captures wind energy and converts the wind energy into mechanical energy, both the MSC and the GSC are two-level voltage source converters, wherein the MSC converts the alternating current output by the permanent magnet synchronous generator into direct current, and the GSC inverts the direct current into power frequency alternating current and integrates the power frequency alternating current into the power grid. Both the MSC and the GSC adopt a vector control method in a dq coordinate system, and convert three-phase voltage and current under the abc frame into the voltage and current under d-axis and q-axis through dq transformation. The MSC controller includes a power loop (P&O), a rotating speed loop (H.sub.ω), and a machine-side current loop (H.sub.cr); and the GSC controller includes a voltage loop (H.sub.v) and a grid-side current loop (H.sub.cg). In addition, the PMSG-WT also includes a phase-locked loop (PLL) and coordinate conversion links (abc/dq and dq/abc). Hereinafter, the present application will be further described in combination with specific embodiments:

[0074] In one embodiment of the present application, main parameters of the system are shown in Table 1.

TABLE-US-00001 TABLE 1 Main parameters of the system Parameter Value Equivalent inductance of a weak 0.1 mH power grid L.sub.g Voltage amplitude of an ideal 690{square root over (2)}/{square root over (3)} voltage source of the weak power grid U.sub.s Filtering inductance L.sub.f 1.2 mH DC capacitor C.sub.dc 2 mF Armature inductance of the 1.5028 mH generator L.sub.s Rotor resistance of generator R.sub.s 0.0224 Ω Permanent magnet flux 5.3445 Wb linkage ψ.sub.f Power frequency angular 100π rad/s frequency ω The number of pole pairs of the 48 generator n.sub.p Equivalent rotational inertia J 1417769.4 Perturbation step length of a 0.0001 rad/s power loop ε Control period of the power 0.1 ms loop T.sub.p Reference value of output power 0.688 MW of the PMSG-WT P.sup.ref Proportional parameter of a 10902 rotating speed loop K.sub.ωp Integral parameter of the rotating 17875 speed loop K.sub.ωi Proportional parameter of a 3.13 machine-side current K.sub.cpr Integral parameter of the 109.8 machine-side current loop K.sub.cir Proportional parameter of a 19.05 grid-side current loop K.sub.cpg Integral parameter of the 176 grid-side current loop K.sub.cig Proportional parameter of a 1.6 grid-side voltage loop K.sub.vp Integral parameter of the 28 grid-side voltage loop K.sub.vi Proportional parameter of the 2.1 phase-locked loop K.sub.ppll Integral parameter of the 36.6 phase-locked loop K.sub.ipll Reference value of DC voltage 1400 V U.sub.dcref Steady-state value of rotating 3 rad/s speed of the generator ω.sub.g* Linearization constant of a −2.742 × 10.sup.5 mechanical torque of the wind turbine B.sub.t Period of a power sampling 0.01 filter T.sub.f

[0075] In the embodiment of the present application, step 1, acquiring main parameters of the PMSG-WT, wherein the main parameters as shown in Table 1, establishing mathematical models of the wind turbine, the generator and the MSC, the MSC controller, the DC capacitor, the GSC and the filter, and the GSC controller, respectively, performing linearization in a dq coordinate system, and calculating steady-state operation parameters to obtain small-signal models of the wind turbine, the generator and the MSC, the MSC controller, the DC capacitor, the GSC and the filter, and the GSC controller:

[00036]$\begin{array}{c}\mathrm{sJ}{\mathrm{\Delta \omega }}_g=\Delta T_m-\Delta T_e\\ \Delta T_e=-\frac{3}{2}{n}_p{\psi }_f\Delta i_{\mathrm{qr}}\\ \Delta T_m=B_t{\mathrm{\Delta \omega }}_g\\ Z_{\mathrm{dqr}}.Math.\Delta i_{\mathrm{dqr}}=-U_{\mathrm{dc}}^{*}.Math.\Delta d_{\mathrm{dqr}}-D_{\mathrm{dqr}}^{*}.Math.\Delta u_{\mathrm{dc}}+G_{\omega e}.Math.{\mathrm{\Delta \omega }}_e\\ \Delta d_{\mathrm{dqr}}^c=\frac{1}{U_{\mathrm{dc}}^{*}}\left(G_{\mathrm{cr}}.Math.\Delta i_{\mathrm{dqr}}^c+G_{\omega g}\left({\omega }_g^{\mathrm{ref}}-{\omega }_g\right)+G_{\omega e}.Math.{\mathrm{\Delta \omega }}_e\right)\end{array}$$\{\begin{array}{c}\Delta i_{\mathrm{dqr}}^c=G_{\theta e}^i.Math.\Delta i_{\mathrm{dqr}}\\ \Delta d_{\mathrm{dqg}}^c=G_{\theta e}^d.Math.\Delta i_{\mathrm{dqr}}+\Delta d_{\mathrm{dqg}}\end{array}$${\mathrm{sC}}_{\mathrm{dc}}\Delta u_{\mathrm{dc}}=1.5\left({D_{\mathrm{dqr}}^{*}}^T.Math.\Delta i_{\mathrm{dqr}}+{I_{\mathrm{dqr}}^{*}}^T.Math.\Delta d_{\mathrm{dqr}}\right)-1.5\left({D_{\mathrm{dqg}}^{*}}^T.Math.\Delta i_{\mathrm{dqg}}+{I_{\mathrm{dqg}}^{*}}^T.Math.\Delta d_{\mathrm{dqg}}\right)$$\begin{array}{c}{Z}_f.Math.\Delta i_{\mathrm{dqg}}=U_{\mathrm{dc}}^{*}.Math.\Delta d_{\mathrm{dqg}}+D_{\mathrm{dqg}}^{*}.Math.\Delta u_{\mathrm{dc}}-\Delta u_{\mathrm{dqg}}\\ \Delta d_{\mathrm{dqg}}^c.Math.U_{\mathrm{dc}}^{*}=\Delta u_{\mathrm{dqg}}^c+G_{\mathrm{cg}}\Delta i_{\mathrm{dqg}}^c+G_{\mathrm{udc}}\Delta u_{\mathrm{dc}}\end{array}$$\{\begin{array}{c}\Delta u_{\mathrm{dqg}}^c=G_{\mathrm{pll}}^u.Math.\Delta u_{\mathrm{dqg}}\\ \Delta i_{\mathrm{dqg}}^c=G_{\mathrm{pll}}^i.Math.\Delta u_{\mathrm{dqg}}+\Delta i_{\mathrm{dqg}}\\ \Delta d_{\mathrm{dqg}}^c=G_{\mathrm{pll}}^d.Math.\Delta u_{\mathrm{dqg}}+\Delta d_{\mathrm{dqg}}\end{array}$

[0076] s represents a parameter introduced by the Laplace transform, J represents rotational inertia of an equivalent concentrated mass block of the wind turbine and the generator, ω.sub.g represents a rotating speed of the generator, T.sub.m represents a mechanical torque of the generator, and T.sub.e represents an electromagnetic torque of the generator. n.sub.p represents the number of pole pairs of the generator, ψ.sub.f represents a permanent magnet flux linkage of the generator, and B.sub.t represents a linearization constant of the mechanical torque of the wind turbine.

[00037]$Z_{\mathrm{dqr}}=\left[\begin{array}{cc}{\mathrm{sL}}_s+R_s& -{\omega }_e^{*}{L}_s\\ {\omega }_e^{*}{L}_s& {\mathrm{sL}}_s+R_s\end{array}\right],\Delta i_{\mathrm{dqr}}=\left[\begin{array}{c}\Delta i_{\mathrm{dr}}\\ \Delta i_{\mathrm{qr}}\end{array}\right],\Delta d_{\mathrm{dqr}}=\left[\begin{array}{c}\Delta d_{\mathrm{dr}}\\ \Delta d_{\mathrm{qr}}\end{array}\right],D_{\mathrm{dqr}}^{*}=\left[\begin{array}{c}{D}_{\mathrm{dr}}^{*}\\ D_{\mathrm{qr}}^{*}\end{array}\right],G_{\omega e}=\left[\begin{array}{c}{L}_s{I}_{\mathrm{qr}}^{*}\\ {\psi }_f-L_s{I}_{\mathrm{dr}}^{*}\end{array}\right],$

R.sub.s and L.sub.s represent rotor resistance and armature inductance of the generator, respectively, ω.sub.e represents an electrical angular speed of the rotor, and ω.sub.e=n.sub.pω.sub.g. i.sub.dr and i.sub.qr represent the stator current of the generator in the dq coordinate system, d.sub.dr and d.sub.qr represent output duty ratios of the MSC controller in the dq coordinate system, and u.sub.dc represents a DC voltage.

[00038]$G_{\mathrm{cr}}=\left[\begin{array}{cc}{H}_{\mathrm{cr}}& {\omega }_e^{*}{L}_s\\ -{\omega }_e^{*}{L}_s& H_{\mathrm{cr}}\end{array}\right],G_{\omega g}=\left[\begin{array}{c}0\\ H_{\mathrm{cr}}{H}_{\omega }\end{array}\right],H_{\mathrm{cr}}=K_{\mathrm{cpr}}+\frac{K_{\mathrm{cir}}}{s},$

K.sub.cpr and K.sub.cir represent a proportional parameter and an integral parameter of machine-side current loop PI control, respectively,

[00039]$H_{\omega }=K_{\omega p}+\frac{K_{\omega i}}{s},$

K.sub.ωp and K.sub.ωi represent the proportional parameter and the integral parameter of rotating speed loop PI control, respectively, ω.sub.g.sup.ref represents the reference value of the rotating speed of the generator, ω*.sub.e represents the steady-state value of the electrical angular speed of the rotor, and the superscript c represents a dq coordinate system of the MSC controller.

[00040]$G_{\theta e}^i=\left[\begin{array}{cc}1& H_{\theta e}{I}_{\mathrm{qr}}^{*}\\ 0& 1-H_{\theta e}{I}_{\mathrm{dr}}^{*}\end{array}\right],G_{\theta e}^d=\left[\begin{array}{cc}0& H_{\theta e}{D}_{\mathrm{qr}}^{*}\\ 0& -H_{\theta e}{D}_{\mathrm{dr}}^{*}\end{array}\right],H_{\theta e}=\frac{3n_p^2{\psi }_f}{2s\left(B_t-\mathrm{sJ}\right)}.I_{\mathrm{dqr}}^{*}\left[\begin{array}{c}{I}_{\mathrm{dr}}^{*}\\ I_{\mathrm{qr}}^{*}\end{array}\right],\Delta i_{\mathrm{dqg}}=\left[\begin{array}{c}\Delta i_{\mathrm{dg}}\\ \Delta i_{\mathrm{qg}}\end{array}\right],\Delta d_{\mathrm{dqg}}=\left[\begin{array}{c}\Delta d_{\mathrm{dg}}\\ \Delta d_{\mathrm{qg}}\end{array}\right],D_{\mathrm{dqg}}^{*}=\left[\begin{array}{c}{D}_{\mathrm{dg}}^{*}\\ D_{\mathrm{qg}}^{*}\end{array}\right],I_{\mathrm{dqg}}^{*}=\left[\begin{array}{c}{I}_{\mathrm{dg}}^{*}\\ I_{\mathrm{qg}}^{*}\end{array}\right],$

C.sub.dc represents the DC capacitor, i.sub.dg and i.sub.qg represent d-axis current and q-axis current of an AC port of the GSC, respectively, and d.sub.dg and d.sub.qg represent the output duty ratios of the GSC controller in the dq coordinate system, respectively.

[00041]$Z_f=\left[\begin{array}{cc}{\mathrm{sL}}_f& -\omega L_f\\ \omega L_f& {\mathrm{sL}}_f\end{array}\right],\Delta i_{\mathrm{dqg}}=\left[\begin{array}{c}\Delta i_{\mathrm{dg}}\\ \Delta i_{\mathrm{qg}}\end{array}\right],\Delta d_{\mathrm{dqg}}=\left[\begin{array}{c}\Delta d_{\mathrm{dg}}\\ \Delta d_{\mathrm{qg}}\end{array}\right],D_{\mathrm{dqg}}^{*}=\left[\begin{array}{c}{D}_{\mathrm{dg}}^{*}\\ D_{\mathrm{qg}}^{*}\end{array}\right],\Delta u_{\mathrm{dqg}}=\left[\begin{array}{c}\Delta u_{\mathrm{dg}}\\ \Delta u_{\mathrm{qg}}\end{array}\right],$

L.sub.f represents filtering inductance, ω represents power frequency angular frequency, ω=100π rad/s, i.sub.dg and i.sub.qg represent the d-axis current and the q-axis current of the AC port of the GSC, respectively, d.sub.dg and d.sub.qg represent the output duty ratios of the GSC controller in the dq coordinate system, respectively, and u.sub.dg and u.sub.qg represent a d-axis voltage and a q-axis voltage of a grid-connected point, respectively.

[00042]$G_{\mathrm{cg}}=\left[\begin{array}{cc}-H_{\mathrm{cg}}& -\omega L_f\\ \omega L_f& -H_{\mathrm{cg}}\end{array}\right],G_{\mathrm{udc}}=\left[\begin{array}{c}{H}_{\mathrm{cg}}{H}_v\\ 0\end{array}\right],H_{\mathrm{cg}}=K_{\mathrm{cpg}}+\frac{K_{\mathrm{cig}}}{s},$

K.sub.cpg and K.sub.cig represent the proportional parameter and the integral parameter of grid-side current loop PI control, respectively,

[00043]$H_v=K_{\mathrm{vp}}+\frac{K_{\mathrm{vi}}}{s},$

and K.sub.vp and K.sub.vi represent the proportional parameter and the integral parameter of grid-side voltage loop PI control, respectively.

[00044]$G_{\mathrm{pll}}^u=\left[\begin{array}{cc}1& G_{\mathrm{pll}}{U}_{\mathrm{qg}}^{*}\\ 0& 1-G_{\mathrm{pll}}{U}_{\mathrm{dg}}^{*}\end{array}\right],G_{\mathrm{pll}}^i=\left[\begin{array}{cc}0& G_{\mathrm{pll}}{I}_{\mathrm{qg}}^{*}\\ 0& -G_{\mathrm{pll}}{I}_{\mathrm{dg}}^{*}\end{array}\right],G_{\mathrm{pll}}^d=\left[\begin{array}{cc}0& G_{\mathrm{pll}}{D}_{\mathrm{qg}}^{*}\\ 0& -G_{\mathrm{pll}}{D}_{\mathrm{dg}}^{*}\end{array}\right],G_{\mathrm{pll}}=\frac{H_{\mathrm{pll}}}{s+U_{\mathrm{dg}}^{*}{H}_{\mathrm{pll}}},H_{\mathrm{pll}}=K_{\mathrm{ppll}}+\frac{K_{\mathrm{ipll}}}{s},$

K.sub.ppll and K.sub.ipll represent the proportional parameter and the integral parameter of a phase-locked loop PI controller, respectively.

[0077] Capital letters and superscripts * represent steady-state components of corresponding lowercase variables, and the specific calculation method is as follows:

[00045]$\{\begin{array}{c}{I}_{\mathrm{qr}}^{*}=3P^{\mathrm{ref}}/n_p{\Psi }_f{\omega }_g^{*}\\ I_{\mathrm{dr}}^{*}=0\\ U_{\mathrm{dc}}^{*}=U_{\mathrm{dcref}}\\ D_{\mathrm{dr}}^{*}=\left({\omega }_e{I}_{\mathrm{qr}}^{*}{L}_s-R_s{I}_{\mathrm{dr}}^{*}\right)/U_{\mathrm{dc}}^{*}\\ D_{\mathrm{qr}}^{*}=\left(-{\omega }_e{L}_s{I}_{\mathrm{dr}}^{*}+{\omega }_e{\Psi }_f-R_s{I}_{\mathrm{qr}}^{*}\right)/U_{\mathrm{dc}}^{*}\\ U_{\mathrm{dg}}^{*}=\sqrt{\left({\left(1.5U_s\right)}^2+\sqrt{{\left(1.5U_s\right)}^4-{\left(3P^{\mathrm{ref}}\omega L_g\right)}^2}\right)/4.5}\\ U_{\mathrm{qg}}^{*}=0\\ I_{\mathrm{dg}}^{*}=P^{\mathrm{ref}}/\left(1.5U_{\mathrm{dg}}^{*}\right)\\ I_{\mathrm{qg}}^{*}=0\\ D_{\mathrm{dg}}^{*}=\left(-I_{\mathrm{qg}}^{*}\omega L_f+U_{\mathrm{dg}}^{*}\right)/U_{\mathrm{dc}}^{*}\\ D_{\mathrm{dg}}^{*}=\left(I_{\mathrm{dg}}^{*}\omega L_f+U_{\mathrm{dg}}^{*}\right)/U_{\mathrm{dc}}^{*}\end{array}$

[0078] Step 2: modeling the power loop in the MSC controller based on a describing function method, wherein its mathematical expression is:

[00046]${\omega }_g^{\mathrm{ref}}=\{\begin{array}{cc}\frac{\epsilon }{s{T}_P}\mathrm{sgn}\left(P^{\mathrm{ref}}-P_n\right),& \left({\omega }_g\le {\omega }_{\mathrm{mpp}}\right)\\ -\frac{\epsilon }{s{T}_P}\mathrm{sgn}\left(P^{\mathrm{ref}}-P_n\right),& \left({\omega }_g>{\omega }_{\mathrm{mpp}}\right)\end{array}$

[0079] the sign function in the formula can be modeled by using a describing function, and the describing function is:

[00047]$N\left(A\right)=\frac{4}{\pi A}$

[0080] Step 3: considering the influence of the weak AC power grid, the weak AC power grid model is:

Z.sub.g.Math.Δi.sub.dqg=Δu.sub.dqg

[0081] in the formula

[00048]$,Z_g=\left[\begin{array}{cc}{\mathrm{sL}}_g& -\omega L_g\\ \omega L_g& {\mathrm{sL}}_g\end{array}\right].$

The small-signal model of the output power of the PMSG-WT is:

ΔP=1.5(I*.sub.dqg.sup.T.Math.Z.sub.g+U*.sub.dqg.sup.T).Math.Δi.sub.dqg

[0082] combing the linear parts of the weak power grid model and the power loop model with the small-signal model of the PMSG-WT in step 1, and deriving a transfer function G(s) of the linear part of the system:

[00049]$G\left(s\right)=-\frac{1.5\epsilon }{T_{P}s}\left(I_{\mathrm{dqg}}^{*T}{Z}_g+U_{\mathrm{dqg}}^{*T}\right)M_5^{-1}\left(D_{\mathrm{dqg}}^{*}+G_{\mathrm{udc}}\right)M_7^{-1}{M}_8.Math.\frac{1}{1+T_{f}s}.Math.\frac{1}{1+T_{P}s}$$\mathrm{where}$$M_8=-D_{\mathrm{dqr}}^{*T}.Math.M_1^{-1}.Math.G_{\omega g}+I_{\mathrm{dqr}}^{*T}.Math.M_3$$M_7=\left(\frac{2s{C}_{\mathrm{dc}}}{3}+I_{\mathrm{dqg}}^{*T}{M}_6-I_{\mathrm{dqr}}^{*T}{M}_4+D_{\mathrm{dqg}}^{*T}{M}_5^{-1}\left(D_{\mathrm{dqg}}^{*}+G_{\mathrm{udc}}\right)+D_{\mathrm{dqr}}^{*T}{M}_1^{-1}{D}_{\mathrm{dqr}}^{*}\right)$$M_6=\left(\left(\left(G_{\mathrm{cg}}{G}_{\mathrm{plli}}-G_{\mathrm{plld}}{U}_{\mathrm{dc}}^{*}+G_{\mathrm{pllu}}\right)Z_g+G_{\mathrm{cg}}\right)M_5^{-1}\left(D_{\mathrm{dqg}}^{*}+G_{\mathrm{udc}}\right)+G_{\mathrm{udc}}\right)/U_{\mathrm{dc}}^{*}$$M_5=-G_{\mathrm{cg}}+Z_f-\left(G_{\mathrm{cg}}{G}_{\mathrm{plli}}-G_{\mathrm{plld}}{U}_{\mathrm{dc}}^{*}+G_{\mathrm{pllu}}\right)Z_g+Z_g$$M_4=-M_2{M}_1^{-1}{D}_{\mathrm{dqr}}^{*}/U_{\mathrm{dc}}^{*}$$M_3=\left(-M_2{M}_1^{-1}{G}_{\omega g}+G_{\omega g}\right)/U_{\mathrm{dc}}^{*}$$M_2=G_{\omega e}{G}_{iq\omega }{n}_p+G_{\mathrm{cr}}{G}_{\theta e}^i-G_{\omega g}{G}_{iq\omega }-U_{\mathrm{dc}}^{*}{G}_{\theta e}^d$$M_1=G_{\mathrm{cr}}{G}_{\theta e}^i-G_{\omega g}{G}_{iq\omega }-U_{\mathrm{dc}}^{*}{G}_{\theta e}^d+z_{\mathrm{dqr}}$

[0083] Step 4: analyzing the stability of the system. Firstly, a G(s) pole diagram is drawn, as shown in FIG. 2, it can be seen that G(s) does not contain the pole of the right half plane (the real part is greater than 0), so the first condition for the stability of the system is satisfied. G(s) and −1/N(A) images are drawn in a complex plane, as shown in FIG. 3, G(s) and −1/N(A) intersect, indicating that the system is in a critically stable state at this time. It can be seen from calculation that, when Lg=0.1 mH, the oscillation frequency of the system is about 129 rad/s (20.5 Hz), and the oscillation amplitude is about 61 kW; and when Lg=0.4 m, the oscillation frequency of the system is about 131 rad/s (20.9 Hz), the oscillation amplitude is about 56 kW. When the power grid strength decreases (Lg increases), the oscillation amplitude of the system decreases, indicating that under certain conditions, the increase of the equivalent impedance of the power grid is conducive to maintaining the stability of the system.

[0084] FIGS. 4a, 4b, 5a and 5b show system simulation and FFT spectrum analysis results when Lg=0.1 mH and Lg=0.4 mH, respectively. In FIGS. 4a and 4b, the DC component is 0.6889 MW, the corresponding oscillation frequency is 19.5 Hz, and the amplitude is 58 kW, which are basically consistent with the theoretical analysis results; and in FIGS. 5a and 5b, the DC component is 0.6889 MW, the corresponding oscillation frequency is 19.5 Hz, and the amplitude is 53 kW, which are basically consistent with the theoretical analysis results. The simulation results verify the effectiveness and accuracy of the analysis method.