Data-Driven Virtual Inertia Control Method of Doubly Fed Wind Turbine

Abstract

This paper presents a data-driven virtual inertia control method for doubly fed induction generator (DFIG)-based wind turbine to provide inertia support in the presence of frequency events. The Markov parameters of the system are first obtained by monitoring the grid frequency and system operation state. Then, a data-driven state observer is developed to evaluate the state vector of the optimal controller. Furthermore, the optimal controller of the inertia emulation system is developed through the closed solution of the differential Riccati equation. Moreover, a differential Riccati equation with self-correction capability is developed to enhance the anti-noise ability to reject noise interference in frequency measurement process. Finally, the simulation verification was performed in Matlab/Simulink to validate the effectiveness of the proposed control strategy. Simulation results showed that the proposed virtual inertia controller can adaptively tune control parameters online to provide transient inertia supports for the power grid by releasing the kinetic energy, so as to improve the robustness and anti-interference ability of the control system of the wind power system.

1. Introduction

The increasing penetration of wind power generation into the power system is becoming an important trend [1,2]. DFIG-based wind turbine is an important member of the commercial wind turbine family. The DFIG-based wind turbine is decoupled from power grid by back-to-back power converter, so that there exists no inherent response from wind turbine in the presence of grid frequency event [3,4]. With the increasing penetration of wind power generation, the impact of wind power system on the power system becomes more and more evident. Therefore, grid codes throughout the world pose the requirements for power regulation to support the frequency stability of the power grid [5,6].

Inertia emulation methods have been frequently proposed to support transient inertia responses in previous work, which can be divided into energy storage-based inertia emulation and rotor kinetic energy-based inertia emulation [7,8]. In [9], the authors provide the calculation method of the equivalent virtual inertia time constant reflecting the effective energy storage of the DFIG and put forward the variable parameter virtual inertia emulation method based on DFIG effective energy storage, so as to improve the ability of DFIG to participate in system frequency regulation. In [10], a DFIG-based wind turbine with hydrogen energy storage is proposed, where the hydrogen energy storage system is adopted to provide inertia support for power grid so as to improve the stability of power grid. In [11], a control method based on removable virtual resistors cooperated with reconfiguration of battery energy storage unit’s control structure is proposed. By installing removable virtual resistors on the energy storage device, the inertia support ability of the DFIG is further enhanced and the low-voltage ride through capability of wind turbine is improved. The different virtual inertia control methods based on energy storage devices adopted in [9,10,11] can provide effective inertia support for the power grid, but there still exists several drawbacks. For instance, energy storage (ES) devices are expensive, and the operation performance of ES devices is sensitive to the environment. Compared with energy storage-based methods, the rotor kinetic energy-based inertia emulation method can provide inertia support by controlling the rotational speed of wind turbine. An inertia control method based on dynamic equations and adaptive fuzzy algorithm is proposed in [12], which can adaptively adjust the inertia control gain within a wide wind speed range to mitigate the transient frequency fluctuation of power system. However, the dynamic response is relative slow. In [13], the relationship among virtual inertia, rotor speed, and grid frequency change are analyzed, and a new virtual inertia emulation method is developed, where the MPPT curve can be switched to virtual inertia control curve according to frequency deviation so as to release kinetic energy and provide inertia support for the power grid. However, this method tends to increase the computational burdens of the control system.

The data-driven virtual inertia emulation methods were developed to improve operation performance of inertia control strategies. In [14], the damping ratio of eigenvalues is selected as the objective function, and a dynamic linearized mathematical model is established to replace the control parameter optimization model. Moreover, the model parameters are solved by offline data, which improves the small signal stability of the system. In [15], an approximate dynamic programming method is proposed for online parameters tuning of traditional PD virtual inertia controller, where the knowledge of system model is not necessary. This method can realize the online optimal parameter tuning, which enhances the control performance of DFIG for frequency regulation, although these data-driven inertia control methods may address effectively transient frequency fluctuation. However, due to the wind farm operating in a very complex external environment, there exists measurement noise that will cause vibration of DFIG shell in the frequency measurement phase. If the impact of measurement noise is not fully considered, this will cause the algorithm to be susceptible to noise interference and decrease the anti-disturbance ability of system. In [16], an integral stability analysis method for stochastic parameter DFIG system considering virtual inertia control is proposed. On the basis of the traditional DFIG virtual inertia control, the functional analytical relationship between stochastic parameters and elements in the state matrix is derived, and the stochastic parameter interconnected system model considering Wiener noise is established, which improves the identification ability of the power system. However, the existing methods mainly focus on the influence of measurement noise on the inertia support capability of wind turbines, but it fails to reduce the computational burdens of the data-driven controller.

To improve the ability of DFIG to provide inertia support for grid, as well as to reduce the impact of measurement noises on inertia control, this paper presents a data-driven optimal control strategy for inertia emulation control. The proposed controller consists of two complementary subsystems, the data-driven observer and the optimal controller based on the system Markov parameters. The data-driven state observer is constructed by monitoring the grid frequency and the operating status of DFIG to obtain the Markov parameters of the system. In addition, the optimal controller based on the Markov parameters of the system is developed by the closed solution of the differential Riccati equation with self-tuning capability. The main contributions of this paper are explained as follows.

• A novel data-driven virtual inertia emulation method is proposed, which is able to achieve the better control performance compared with traditional inertia control strategies.
• A data-driven state observer and an optimal controller based on system Markov parameters are developed. In this proposed controller, only real-time input and output data of the system is required to complete the relevant computations without using a detailed mathematical model.
• In this work, the differential Riccati equation with self-correction ability is introduced on the basis of the aforementioned optimal controller, which improves the anti-disturbance ability of system.

The rest of this paper is organized as follows. The basic principle of DFIG and the traditional virtual inertia control strategies are introduced in Section 2. In Section 3, the data-driven inertia control strategy is developed. In Section 4, simulation verification is performed to validate the proposed control method in Matlab/Simulink. The conclusions are drawn in Section 5.

2. System Description and Problem Formulation

2.1. System Description

Figure 1 shows the diagram of the traditional DFIG-based wind turbine, which is composed of doubly fed induction machine, rotor side converter, grid side converter, DFIG controller, and power grid. The rotor side converter and grid side converter can decouple DFIG from power grid and perform variable speed operation.

The power converter is enabled by control system to perform the following functions. (1) To implement the maximum power point tracking when the wind speed is lower than the rated wind speed. (2) To maintain the DC-link voltage of power converter. (3) To provide the reactive power support in the presence of grid faults.

The mechanical power extracted from wind turbine can be represented as (1) [17,18].

$$P_m=0.5\pi \rho R^2{C}_p(\lambda ,\beta )v^3$$

(1)

where P_m is the mechanical power extracted from the wind energy by the wind turbines, ρ = 1.22 kg/m³ is the air density in kilograms per cubic meter, R is the turbine radius, ω is the turbine angular speed, v is the wind speed, λ is the blade top speed ratio, and the λ = ωR/v. β is the blades pitch angle, and C_p(λ,β) is the power performance coefficient of the turbine, wherein its value is related to the blades pitch angle and the blade top speed ratio, which can be represented as (2):

$$\left\{\begin{array}{l}{C}_p(\lambda ,\beta )=0.5176(\frac{116}{{\lambda }_i}-0.4\beta -5)e^{-21/{\lambda }_i}+0.0068\lambda \\ {\lambda }_i=\frac{1}{\frac{1}{\lambda +0.08\beta }-\frac{0.035}{{\beta }^3+1}}\end{array}\right.$$

(2)

In order to make the C_p(λ,β) reach the optimal value C_pmax(λ,β), the rotor speed of the DFIG is controlled by wind power converter to capture the optimal tip speed ratio λ_opt. The wind turbine parameter values provided by the wind turbine supplier and related documents [19] can be used to obtain the relationship between the wind energy utilization coefficient, the tip speed ratio, and the pitch angle at the optimal operating point of the wind turbine. In order to analyze the influence of virtual inertia on the wind turbine, this paper mainly analyzed the influence of virtual inertia control on C_p and λ when the pitch angle is near zero. Figure 2 shows the relationship curve between tip speed ratio and wind energy utilization coefficient when β = 0°. It can be seen that the maximum value of wind energy utilization coefficient is C_pmax ≈ 0.48 when the blade tip speed ratio reaches the maximum value, and the optimal value of blade tip speed ratio is λ_opt ≈ 8.1. At this time, the wind energy utilization factor is the highest. When the gird load suddenly increases, the system frequency will decrease. Then, the rotational speed of wind turbine is reduced to release the kinetic energy of the rotor, where the wind turbine can be operated in region 1. At this time, the value of tip speed ratio λ will decrease. Therefore, C_p will decrease from the highest value C_pmax to the left along the curve of the relationship tip speed ratio and the wind energy utilization coefficient. On the other hand, when the gird load decreases, the system frequency will increase. Then, the rotational speed can be increased to support frequency regulation, where the wind turbine can be operated in region 2. At this time, the value of tip speed ratio λ will increase. Therefore, C_p will decrease from the highest value C_pmax to the right along the curve of the relationship tip speed ratio and the wind energy utilization coefficient.

The model of the wind energy conversion system with DFIG is described as follows:

$$\left\{\begin{array}{l}{u}_s=-R{i}_s+p{\psi }_s+j{\omega }_s{\psi }_s\\ {\psi }_s=-L_s{i}_s+L_m{i}_r\\ T_e=n_p{L}_m{(i_{sq}{i}_{rd}-i_{sd}{i}_{rq})}^{}{{}^{}}_{}\\ P_s=u_{sd}{i}_{sd}+u_{sq}{i}_{sq}\end{array}\right.$$

(3)

where ψ_s is the stator flux, P_s is the stator active power, u_s is the stator voltage, i_s is the stator current, and i_r is the rotor current. p is the differential operator, L_s is the stator inductance, L_m is the stator mutual inductance, ω_s is the synchronous speed, and R_s is the stator resistance. T_e is the electromagnetic torque of the motor and n_p is the pole-pairs number of the motor.

2.2. PD Virtual Inertia Controller

The wind turbine has the ability to adaptively address frequency response if a frequency power outer loop is added to the DFIG control system. Figure 3 shows the control structure. The power outer loop includes two parts: the proportional part acts on the dynamic frequency deviation (Δf), and the differential part acts on the frequency deviation change rate (dΔf/dt). In the initial stage when the power grid is out of balance, the system frequency responds quickly, and the differential part plays a major role at this time. When the system frequency drops or rises seriously, the proportional part will play a major role. The active power compensated to DFIG can be given as (4).

$$P_c=K_d\frac{d\Delta f}{dt}+K_p\Delta f$$

(4)

where P_c is the inertia response reference of active power, K_d is the differential parameter, K_p is the proportional parameter, Δf is the frequency offset, and dΔf/dt is the frequency rate.

In order to improve the transient performance of the system frequency, the equivalent inertia of the power system is increased by releasing the rotor energy storage and simulating the inertia response of the synchronous generator. The inertia response of the power system can be given as (5) [20].

$$\Delta P=2H\frac{d\Delta f}{dt}+D\Delta f$$

(5)

where H is the inertia time constant, and D is the damping coefficient. Considering the inertia response of the wind turbine, the traditional PD virtual inertia control is given as (6).

$$\Delta P_c=\left(2H+K_d\right)\frac{d\Delta f}{dt}+\left(D+K_p\right)\Delta f$$

(6)

It is easy to perform the PD-based virtual inertia control strategy in practical operation. However, the PD-based virtual inertia control strategy fails to provide more effective inertia support according to specific conditions in actual wind power systems. In (5), the frequency change rate measured from the power grid is susceptible to the effect of measurement noises. These measurement noises are, in the majority, caused by aging and parameter drift of measuring equipment and/or the difference of electromechanical characteristics between different equipment, leading to the difference of frequency measurements [21]. Due to the existence of measurement noise in the process of grid frequency measurement, the DFIG inertia emulation system has stronger uncertainty, and therefore it is urgent to develop an accurate DFIG inertia emulation system and to address the frequency change rate measurement noise so as to improve control effect of virtual inertia. To address the above problems, this paper presents a data-driven optimal control strategy for inertia emulation control.

3. Data-Driven Virtual Inertia Control Method

In this section, a data-driven virtual inertia control method is developed to perform inertia control. Figure 4 shows the diagram of data-driven inertia emulation controller, which is composed of state observer and the optimal controller based on the system Markov parameters. In this proposed inertia controller, a data-driven state observer based on the Markov parameters of the system is established, which can effective solution to estimate the state vector x_c(k) of the optimal controller [22,23]. Then, the optimal controller of the inertia emulation system is designed through the closed solution of the differential Riccati equation with self-correction capability.

3.1. Description of the Optimal Problem

The power reference value P_ref in the wind turbine control system are changed by traditional virtual inertia control, which makes the DFIGs show controllable inertia response characteristics when the system frequency fluctuates. In order to improve the inertia control performance, a data-driven inertia control method is proposed in this paper. The DFIG inertia emulation system at time k can be given in a state-space form as (7).

$$\left\{\begin{array}{c}x\left(k+1\right)=Ax\left(k\right)+Bu\left(k\right)\\ y\left(k\right)=Cx\left(k\right)\end{array}\right.$$

(7)

where x(k) is the n-dimensional state vector, u(k) is the m-dimensional control input vector, and y(k) is the l-dimensional measurement output vector. A, B, C are the parameter matrix. Then, the LQ optimal control is used to solve the control function.

$$u\left(k\right)=f\left[u\left(0\right),u\left(1\right),\cdots ,u\left(k-1\right),y\left(0\right),y\left(1\right),\cdots ,y(k-1)\right]$$

(8)

The performance index function is given as (9).

$$\min J=\epsilon \left\{y{\left(N\right)}^{T}Qy\left(N\right)+{\displaystyle \sum _{k=0}^{N-1}y{\left(k\right)}^{T}Qy\left(k\right)+u{\left(k\right)}^{T}Ru\left(k\right)}\right\}$$

(9)

where Q is a positive semi-definite symmetric weight matrix and R is a positive definite symmetric weight matrix.

When the state space model is unknown, the Markov parameter M_i can be obtained by using the online input and output data. In addition, the appropriate weights Q_W and R_W are selected to solve the optimal control function as shown in (8), which can be given as (10).

$$u\left(k\right)=f\left[Q_W,R_W,y\left(k-1\right),u\left(k-1\right),\left\{M_i\right\}\right]$$

(10)

where u is input vector of the DFIG inertia emulation system that is the grid frequency deviation Δf, and it is also the output of the wind power control system. y is output vector of the DFIG inertia emulation system that is the compensated active power ΔP_c, and it is also the input of the wind power control system. Q_W is positive semi-definite symmetric weight matrix, R_W is positive definite symmetric weight matrix, M_i = CA⁽ⁱ⁻¹⁾B; M_i is the Markov parameter of the wind energy conversion system, where i = 1, 2, …, n.

The optimal control function is obtained by solving (11)

$$\min J=\epsilon \left\{\Delta f{\left(N\right)}^T{Q}_W\Delta f\left(N\right)+{\displaystyle \sum _{k=0}^{N-1}\left(\Delta f{\left(k\right)}^T{Q}_W\Delta f\left(k\right)+\Delta P_c{\left(k\right)}^T{R}_W\Delta P_c\left(k\right)\right)}\right\}$$

(11)

to reach the minimum. Therefore, the DFIG inertia emulation system can make the frequency f accurately track the rated value f_N when the power supply load is unbalance so as to provide the inertia support.

3.2. Data-Driven Algorithm

The state vector of the DFIG inertia emulation system for p(p > 0) times can be obtained as (12).

$$\left\{\begin{array}{c}x\left(k+p\right)=A^{p}x\left(k\right)+B_p{u}_p\left(k\right)\\ y_p\left(k\right)=O_{p}x\left(k\right)+T_p{u}_p\left(k\right)\end{array}\right.$$

(12)

where u_p(k) and y_p(k) are respectively derived from u(k) and y(k), starting at first p-step input and output data form a column vector as (13) and (14).

$$u_p\left(k\right)={\left[\begin{array}{ccc}{u}_1\left(k\right)& \cdots & u_p\left(k+p-1\right)\end{array}\right]}^T$$

(13)

$$y_p\left(k\right)={\left[\begin{array}{ccc}{y}_1\left(k\right)& \cdots & y_p\left(k+p-1\right)\end{array}\right]}^T$$

(14)

The matrix B_p in (12) is an n × pm controllability matrix, O_p is a pl × n observability matrix, and T_p is a pl × pm Toeplitz matrix of the system Markov parameters, where

$$\begin{gathered} B_{\mathrm{p}}=\left[\begin{array}{cccc}{A}^{\mathrm{p}-1}B& \cdots & AB& B\end{array}\right],O_{\mathrm{p}}=\left[\begin{array}{c}C\\ CA\\ ⋮\\ C{A}^{\mathrm{p}-1}\end{array}\right] \\ {T}_{\mathrm{p}}=\left[\begin{array}{ccccc}0& \cdots & \cdots & \cdots & 0\\ CB& 0& \ddots & & ⋮\\ CAB& CB& \ddots & \ddots & ⋮\\ ⋮& ⋮& \ddots & \ddots & ⋮\\ C{A}^{\mathrm{p}-2}B& C{A}^{\mathrm{p}-3}B& \cdots & CB& 0\end{array}\right] \end{gathered}$$

Lemma 1.

[24] Assuming that the observable system model, the controllable matrix B_p and the diagonal constant matrix T_p are known; as long as pm ≥ n, there is a matrix M makes A^p + MO_p = 0 holds.

Applying Lemma 1, (12) can be estimated by input and output as (15).

$$\left\{\begin{array}{l}x\left(k+p\right)=\left(B_p+M{T}_p\right)u_p\left(k\right)-M{y}_p\left(k\right)\hfill \\ y_p\left(k+p\right)=P_1{u}_p\left(k\right)+T_p{u}_p\left(k+p\right)+P_2{y}_p\left(k\right)\hfill \end{array}\right.$$

(15)

where

$$P_1=O_{\mathrm{p}}\left(B_{\mathrm{p}}+M_i{T}_{\mathrm{p}}\right)$$

$$P_2=-O_{\mathrm{p}}{M}_i$$

$$M_i=-{\mathrm{A}}_i{}^{\mathrm{p}}{O}_{\mathrm{p}}{}^{+}\frac{{\partial }^2\mathsf{\Omega }}{\partial u^2}$$

In order to obtain Markov parameters of the data-driven system, the input and output data column vectors are formed into matrices Y and V in the form of (16) and (17), and the relationship equation between the input and output data can be given as (18) [25].

$$Y=\left[\begin{array}{ccc}{y}_p\left(k+p\right)& y_p\left(k+p+1\right)\hspace{1em}\cdots & y_p\left(k+p+L\right)\end{array}\right]$$

(16)

$$V=\left[\begin{array}{cccc}{u}_p\left(k\right)& u_p\left(k+1\right)& \cdots & u_p\left(k+L\right)\\ u_p\left(k+p\right)& u_p\left(k+p+1\right)& \cdots & u_p\left(k+p+L\right)\\ y_p\left(k\right)& y_p\left(k+1\right)& \cdots & y_p\left(k+L\right)\end{array}\right]$$

(17)

$$\left[\begin{array}{ccc}{P}_1& T_p& P_2\end{array}\right]=Y{V}^T{\left(V{V}^T\right)}^{+}$$

(18)

P₁, P₂, and matrix T_p can be obtained by solving (18). When p = N + 1, the Markov parameter M_i = CA⁽ⁱ⁻¹⁾B of the DFIG inertia emulation system can be extracted from T_p.

Lemma 2.

[26] The closed solution formula of the difference Riccati equation is given as (19).

$$X\left(k\right)=C{\left(k\right)}^{T}Q\left(k\right)C\left(k\right)-C{\left(k\right)}^{T}Q\left(k\right)S\left(k\right)·{\left(R\left(k\right)+S{\left(k\right)}^{T}Q\left(k\right)S\left(k\right)\right)}^{-1}S{\left(k\right)}^{T}Q\left(k\right)C\left(k\right)$$

(19)

where

$$C\left(k\right)=\left[\begin{array}{c}CA\\ C{A}^2\\ ⋮\\ C{A}^{N-k}\end{array}\right],S\left(k\right)=\left[\begin{array}{ccccc}0& \cdots & \cdots & \cdots & 0\\ CB& 0& \ddots & & ⋮\\ CAB& CB& \ddots & \ddots & ⋮\\ ⋮& ⋮& \ddots & \ddots & ⋮\\ C{A}^{N-K-1}B& C{A}^{N-K-2}B& \cdots & CB& 0\end{array}\right],S\left(N\right)=0$$

R(k) = diag (R, R,⋯, R) and Q(k) = diag (Q, Q,⋯, Q) are both N − k + 1 dimensional matrices.

For the inertia emulation system with given Markov parameters, the closed solution of the differential Riccati equation is introduced. Applying Lemma 2, the data-driven controller is given as (20).

$$u\left(k\right)=G_w\left(k\right)x_c\left(k\right)$$

(20)

where G_w(k) is the gain of data-driven controller, and x_c(k) is the state vector of data-driven controller, which can be given as (21) and (22).

$$G_w\left(k\right)=-\left(R_W+\theta {\left(k+1\right)}^T\mathsf{\Omega }\left(k+1\right)\right.{\left.·\theta \left(k+1\right)\right)}^{-1}\theta {\left(k+1\right)}^T\Omega \left(k+1\right)$$

(21)

$$x_c\left(k\right)=C\left(k+1\right)Ax\left(k\right)=O_{p}x\left(k\right)=\left[\begin{array}{c}CA\\ C{A}^2\\ ⋮\\ C{A}^{N-k}\end{array}\right]x\left(k\right)$$

(22)

where

$$\mathsf{\theta }\left(k+1\right)=C\left(k+1\right)B={\left[\begin{array}{cccc}{M}_1& M_2& \cdots & M_{N-k}\end{array}\right]}^T$$

$$\mathsf{\Omega }\left(k+1\right)=Q\left(k+1\right)-Q\left(k+1\right)S\left(k+1\right)·{\left(R\left(k+1\right)+S{\left(k+1\right)}^{T}Q\left(k+1\right)S\left(k+1\right)\right)}^{-1}·S{\left(k+1\right)}^{T}Q\left(k+1\right)$$

$$S\left(k+1\right)=\left[\begin{array}{ccccc}0& \cdots & \cdots & \cdots & 0\\ M_1& 0& \ddots & & ⋮\\ M_2& M_1& \ddots & \ddots & ⋮\\ ⋮& ⋮& \ddots & \ddots & ⋮\\ M_{N-k-1}& M_{N-k-2}& \cdots & M_1& 0\end{array}\right]$$

S(N) = 0, and R(k + 1) = diag (R_W, R_W,⋯, R_W) and Q(k + 1) = diag (Q_W, Q_W,⋯, Q_W) are both N − k dimensional matrices.

3.3. State Vector Estimation

In (22), the state vector x_c(k) of the optimal controller based on data is related with A, C, and the state variable x(k). In the case of unknown system model, this paper presents a new method to construct a state observer to estimate x_c(k), which can be estimated by taking the state vector x(k) of (15) into (22) when p = N + 1.

$$\begin{array}{l}{\widehat{x}}_c\left(k\right)=\left[\begin{array}{cc}{0}_l& I_{\left(N-k\right)l}\end{array}\right]O_p\widehat{x}\left(k\right)\\ =\left[\begin{array}{cc}{0}_l& I_{\left(N-k\right)l}\end{array}\right]\left(P_1{u}_p\left(k-p\right)+P_2{y}_p\left(k-p\right)\right)\end{array}$$

(23)

x_c(k) is the column vector of ${\widehat{x}}_c\left(k\right)$ from row l + 1 to row l(n − k + 1), so that the state vector of the data-driven optimal controller can be estimated by the input and output vectors and P₁ and P₂.

3.4. Self-Correction of the State Observer under Measurement Noise

To improve the robustness of control system with immunity to measurement noises in practical operation, the control system with capability of anti-interference and noise reduction is developed in this section. The measurement noise is emulated by adding random noise signals, which is introduced into Riccati equation and the data-based state observer. In addition, system frequency is measured by the phase-locked loop. Then, an optimal controller with self-correction ability under correlated noise is developed. The DFIG inertia emulation system with measurement noise can be given as (24).

$$\left\{\begin{array}{c}x\left(k+p\right)=A^{p}x\left(k\right)+B_p{u}_p\left(k\right)+\mathsf{\Gamma }{\omega }_0\left(\mathrm{k}\right)\\ y_p\left(k\right)=O_{p}x\left(k\right)+T_p{u}_p\left(k\right)+v\left(k\right)\end{array}\right.$$

(24)

where ω₀(k) and v(k) are the input noise and measurement noise, respectively, and A^p, O_p, and Γ are the state transition matrix, measurement matrix, and input noise matrix, respectively, which contain the unknown parameter vector θ.

Assuming that ω₀(k) and v(k) are the correlated white noises with zero mean, variances Q₀ and R, and correlated matrix S, then (25) can be obtained.

$$E\left\{\left[\begin{array}{c}{\omega }_0\left(k\right)\\ v\left(k\right)\end{array}\right]\left[\begin{array}{cc}{\omega }_0{}^T\left(k\right)& v^T\left(k\right)\end{array}\right]\right\}=\left[\begin{array}{cc}{Q}_0& S\\ S^T& R\end{array}\right]$$

(25)

where E represents mathematical expectation, and T represents transpose.

When the θ is known, (24) can be represented as (26).

$$x\left(k+p\right)=A^{p}x\left(k\right)+B_p{u}_p\left(k\right)+\mathsf{\Gamma }{\omega }_0\left(k\right)+J\left[y_p\left(k\right)-O_{p}x\left(k\right)-T_p{u}_p\left(k\right)-v\left(k\right)\right]$$

(26)

where, take J = ΓSR⁻¹, and define:

$$\mathsf{\Phi }=A^p-O_p,\overline{\omega }\left(k\right)=\mathsf{\Gamma }{\omega }_0\left(k\right)-Jv\left(k\right)$$

and (26) can be represented as (27).

$$x\left(k+p\right)=\mathsf{\Phi }x\left(k\right)+J{y}_p\left(k\right)+\left(B_p-J{T}_p\right)u_p\left(k\right)+\overline{\omega }\left(k\right)$$

(27)

Then, white noise ω(k) is defined as (28).

$$\omega \left(k\right)={\omega }_0\left(k\right)-S{R}^{-1}v\left(k\right)$$

(28)

(27) is transformed into (29).

$$x\left(k+p\right)=\mathsf{\Phi }x\left(k\right)+J{y}_p\left(k\right)+\left(B_p-J{T}_p\right)u_p\left(k\right)+\mathsf{\Gamma }\omega \left(\mathrm{k}\right)$$

(29)

and the variance of ω(k) is Q = Q₀-SR⁻¹S^T.

In this way, (22) is rewritten as (30).

$$x\left(k+p\right)=\left(B_p+M{T}_p-J{T}_p\right)u_p\left(k\right)-\left(J-M\right)y_p\left(k\right)+\mathsf{\Gamma }\omega \left(k\right)+Mv\left(k\right)$$

(30)

When p = N + 1, x_c(k) can be obtained by taking (30) into (21).

$$\begin{array}{l}{\widehat{x}}_c\left(k\right)=\left[\begin{array}{cc}{0}_l& I_{\left(N-k\right)l}\end{array}\right]O_p\widehat{x}\left(k\right)\\ =\left[\begin{array}{cc}{0}_l& I_{\left(N-k\right)l}\end{array}\right]\left[\left(P_1-J{T}_p\right)u_p\left(k-p\right)+\left(P_2+J\right)y_p\left(k-p\right)+O_p\mathsf{\Gamma }\omega \left(k\right)-P_{2}v\left(k\right)\right]\end{array}$$

(31)

3.5. Algorithm Flow of Markov Riccati Controller

Figure 5 shows the flowchart of the proposed data-driven inertia control algorithm. When the data-driven controller is integrated into the DFIG inertia emulation system, the input and output data are collected to construct matrix Y and V. In order to obtain the accurate Markov parameters, the rows related to input vector in V must be linearly independent. Therefore, when selecting the length L of the data, it is necessary to ensure that the obtained data set is sufficient. The specific algorithm flow is given as follows.

(1)

The matrix Y and V is formulated by system output data and control input data collected in DFIG system according to (16) and (17).

(2)

The P₁, P₂, and matrix T_P are obtained by solving (18). Then, the Markov parameters M_i (i = 1, 2,⋯, n) of DFIG inertia emulation system is obtained from T_P by taking p = N + 1.

(3)

A set of data-based control gains G_w(k) is derived by substituting M_i into (21).

(4)

The estimated value of state vector is obtained by substituting P₁, P₂, and related input and output data into (31), and the state vector x_c(k) of optimal controller is obtained by taking the estimated value of state vector’s second row to (n − k + 1)l row column vector.

(5)

The control input u(k) is calculated according to (20).

(6)

The new output data y(k) is obtained by taking the control input u(k) into wind power system, and the next moment state vector x_c(k + 1) could be estimated by substituting the new control input u(k) and system output y(k) into Equation (31); then, the next moment control input u(k + 1) is obtained, which is cycled in turn.

4. Simulation Verification

In order to validate the effectiveness of the proposed data-driven inertia control method, we established a simulation model of a 1.5 MW DFIG-based wind turbine in Matlab/Simulink. The parameters of the system model are given in

Table 1. Simulation parameters of DFIG-based wind turbine.

Parameter	Value
Air density, ρ	1.22 kg/m3
DFIG rated voltage, Veref	575 V
DFIG rated power, Peref	1500 kW
Inertia of generator, Hg	0.685 s
Wind turbine inertia constant, HWT	4.32 s
Magnetizing inductance, Lm	2.9 (p.u.)
Stator self-inductance, Ls	0.18 (p.u.)
Rotor self-inductance, Lr	0.16 (p.u.)
Rated wind speed, vr	11 m/s
Sampling time of data driven controller, T	0.00001 s
QW	1
RW	0.0001
N	8
p	9
L	37

. Considering simulation time and accuracy of the model [27], we chose the length of the control horizon as N = 8, p = N + 1 = 9, and L = 3p + 10 = 37. Moreover, considering the switching loss of power converter, we chose the sampling period of PWM in power converter as T = 100 us [28,29]. Since the control bandwidth of the proposed data-driven controller is within the control bandwidth of the DFIG power controller, it thus can be implemented in the existing DFIG power controller with same sampling time. In addition, the dead band and filters were installed in the data acquisition link of the data-driven controller, where they aimed to filter out invalid signals.

The inertia emulation ability of wind turbine is important to provide inertia response similar with synchronous generators to reduce the frequency change rate of the system in the presence of frequency event. Therefore, under the condition of 11 m/s wind speed, the simulation study investigated the transient inertia support capability of the proposed control method under three different conditions, namely, system frequency decrease, system frequency increase, and system frequency decrease with measurement noises. In order to validate the advantages of the proposed method, we provided the comparative analysis to investigate the frequency regulation performance among three control methods under the same conditions.

NO Inertia Control: DFIG was operated without inertia emulation capability.

PD Inertia Control: DFIG was operated with PD-based inertia control to provide transient inertia support.

The proposed data-driven inertia control: DFIG was operated with the proposed inertia controller to provide transient inertia support.

4.1. Case I: System Frequency Decrease

When the system frequency deviation is higher than 0.03 Hz, the proposed virtual inertia control method is activated. The specific calculation process of the controller is given as follows. According to the sampled data, we can obtain the Markov parameters of the inertia emulation system.

$${\left\{M_i\right\}}_1^8=\left\{0.0142\hspace{1em}0.0145\hspace{1em}2.5696\hspace{1em}0.0155\hspace{1em}0.0211\hspace{1em}-0.0139\hspace{1em}-0.0091\hspace{1em}7.9556\right\}$$

According to step (3), a set of control gain G_w(k) (k = 1, 2,⋯, 8) is obtained as follows:

$$\begin{array}{l}{G}_w\left(1\right)=\left[-88.8433\hspace{1em}-54.3848\hspace{1em}44.8816\hspace{1em}-38.6886\hspace{1em}-60.5618\hspace{1em}50.2962\hspace{1em}-2.9208\hspace{1em}33.1088\right]\\ G_w\left(2\right)=\left[-82.7040\hspace{1em}-54.1597\hspace{1em}21.0978\hspace{1em}-47.1021\hspace{1em}-25.7684\hspace{1em}39.0785\hspace{1em}12.7832\right]\\ G_w\left(3\right)=\left[-58.4073\hspace{1em}-53.7086\hspace{1em}17.9921\hspace{1em}-28.2684\hspace{1em}-29.6122\hspace{1em}43.0438\right]\\ G_w\left(4\right)=\left[-16.8306\hspace{1em}-33.7387\hspace{1em}51.3757\hspace{1em}-22.0870\hspace{1em}-10.9116\right]\\ G_w\left(5\right)=\left[-37.0653\hspace{1em}-48.9030\hspace{1em}52.2306\hspace{1em}-30.8350\right]\\ G_w\left(6\right)=\left[-54.9232\hspace{1em}-42.3938\hspace{1em}33.4975\right]\\ G_w\left(7\right)=\left[-20.5507\hspace{1em}-8.0355\right]\\ G_w\left(8\right)=\left[-47.1136\right]\end{array}$$

According to step (4), the first state vector is obtained as follows:

$$x_c\left(1\right)={\left[0.0673\hspace{1em}0.0577\hspace{1em}0.0453\hspace{1em}0.0425\hspace{1em}0.0329\hspace{1em}0.0188\hspace{1em}0.0223\hspace{1em}0.0320\right]}^T$$

According to step (5), G_w(k) and x_c(k) are multiplied to obtain the control input ΔP_c:

$$\Delta P_c=\left\{-0.8779\hspace{1em}-0.9564\hspace{1em}-0.7580\hspace{1em}-0.2048\hspace{1em}-0.4259\hspace{1em}-0.4624\hspace{1em}-0.1846\hspace{1em}-0.3171\right\}$$

The next moment input is obtained by taking the calculated ΔP_c into the DFIG power compensation system, and cycle in turn. The specific simulation results are shown in Figure 6.

Figure 6a shows the comparative analysis of DFIG output active power under different control methods when system frequency decreases. It can be seen that the output active power of PD inertia control and the proposed data-driven inertia control was higher than that of NO inertia control. The active power was increased by DFIG by reducing the rotational speed.

Figure 6b shows the comparative analysis of DFIG rotor speed under different control methods when system frequency decreased. It can be seen that the rotational speed of wind turbine was reduced, so as to release the rotational kinetic energy and provide inertia support for power system.

Figure 6c shows the comparative analysis of DFIG reference active power under different control strategies when system frequency decreased. It can be seen that the active power reference of inertia control under the proposed data-driven inertia control and PD inertia control was increased to provide inertia support for power system.

Simulation results showed that compared with PD inertia control, the proposed data-driven inertia control method was able to activate more kinetic energy, which thus provided effective inertia supports for the power grid.

4.2. Case II: System Frequency Increase

When the system frequency deviation is higher than 0.03 Hz, the data-driven virtual inertia controller is activated. Figure 7 shows the simulation results regarding this case.

Figure 7a shows the comparison of DFIG output active power under different control methods when system frequency increased. It can be seen that the output active power of PD inertia control and the proposed data-driven inertia control was lower than that of NO inertia control.

Figure 7b shows the comparison of DFIG rotor speed under different control methods when system frequency increased. It can be seen that compared with PD inertia control, the use of proposed data-driven inertia control can effectively increase the rotor speed of DFIG, thereby suppressing the output power of DFIG and providing inertia support for the power system.

Figure 7c shows the comparison of DFIG reference active power under different control strategies when system frequency increased. It can be seen that the reference active power of the inertia control under the proposed data-driven inertia control and PD inertia control were reduced after system frequency increased, which in turn affected the DFIG rotor speed and provided inertia support for the power system.

Simulation results showed that when system frequency increased, the proposed data-driven inertia control method can contribute to transient frequency stability by increasing rotational speed of wind turbine.

4.3. Case III: System Frequency Decrease with Measurement Noise

The aim of this case was to investigate the anti-interference ability of inertia emulation controller with measurement noises, thus NO inertia control strategies were not related with this issue. To emulate the measurement noise, we added a noise signal with a random fluctuation range of ±0.02 to the measured frequency feedback signal to simulate the noise amplification effect of frequency sampling in the DFIG system. When the system frequency deviation is higher than 0.03 Hz, the data-driven virtual inertia control method is activated. Figure 8 shows the simulation results.

Figure 8a shows that the comparison of system frequency under different control methods in the presence of noise and system frequency decreases. It can be seen that the performance of the PD inertia control was significantly affected by measurement noise, which caused the system frequency fluctuation from 50 to 49.82 Hz with obvious glitches. When the control method was the proposed data-driven inertia control, the system frequency drop was relatively improved compared with PD inertia control, and it was decreased from 50 to 49.85 Hz with slightly smaller glitches. After the data-driven virtual inertia control method was adopted, the lowest point of frequency drop was only 0.02 Hz lower than the case without noise signal.

Figure 8b shows the comparison of DFIG output active power under different control methods in the presence of noise and system frequency decreases. When the system had measurement noise, the inertial response speed and active power output of PD inertia control appeared to be weakened.

Figure 8c shows the comparison of DFIG rotor speed under different control methods in the presence of noise and system frequency decreases.

Simulation results showed that when system frequency decreased with measurement noise, compared with the traditional PD inertia control, the proposed data-driven inertia control method can provide effective inertia support to the power grid in the presence of measurement noise. Moreover, the rotational speed of wind turbine and active power output were merely affected by measurement noise.

5. Conclusions

This paper presents a novel data-driven virtual inertia emulation method to provide the transient inertia supports for power grid in the presence of frequency events. In this proposed inertia controller, a data-driven state observer based on the Markov parameters of the system was established. Then, the optimal controller of the inertia emulation system through the closed solution of the differential Riccati equation with self-correction capability was developed. Compared with traditional inertia emulation methods, the proposed method is able to adaptively adjust the active power output of DFIG through the real-time input and output data of the system. Moreover, the developed controller has good robustness, allowing it to reject the effects of measurement noise. The simulation results show that the proposed data-driven inertia emulation method can effectively provide the transient inertia support for power grid under different operation conditions, which further contributes to the stability and security of the power grid in the presence of frequency events.