**1. Introduction**

Recently, the increasing consumption of electrical energy, depletion of fossil fuels and the environmental problems related to using the non-renewable sources have promoted a growing interest in renewable energies [1,2]. Wind power as a renewable source represents an important and a promising solution for electrical demand and has recently gained more attention for its significant merits of cleanness and resource abundances. Various wind energy configurations are produced from the intensive studies and researches that are carried out in wind systems. One of the most common configurations is the grid-connected Doubly Fed Induction Generators (DFIGs) equipped with variable speed Wind Turbines (WTs). This configuration was widely installed in wind industry due to its significant merits such as independent control of active and reactive powers, low converters' costs and mechanical stress reduction [2].

A DFIG is a wound rotor generator in which the stator windings are directly connected to the grid. The rotor windings are connected to the grid through a back-to-back power converter that is composed from the Rotor Side Converter (RSC) and Grid Side Converter (GSC) components. This configuration allows the converter to handle the fraction of 20% to 30% of the total power. Therefore, the losses in the converter are lower compared to a system where the converter has to handle the full power leading to an improvement in the total e fficiency [2,3].

However, the switching of Insulated-Gate Bipolar Transistors (IGBTs) at both the RSC and GSC components is usually performed by Sinusoidal Pulse Width Modulation (SPWM) strategy, which generates the high frequency content in the utility current. A simple L-filter or an LCL-filter is adopted to reduce the Total Harmonic Distortion (THD) of the utility current [4]. To meet the grid code requirements, the LCL-filter is predominant in reducing the utility current harmonics. Indeed, it can lead to a better attenuation of harmonics using small values of inductances. However, the resonance phenomena of the LCL-filter must be damped properly in order to prevent the possible instability of these systems [4]. One way for solving this problem is employing a passive damping circuit. The passive damping is built by inserting a resistor branch in series or in parallel with the inductors or capacitors of the LCL-filter [4].

The conventional control scheme of the grid connected DFIG wind turbine system is built based on a vector control method [5–8]. Here, the Stator Flux Orientation (SFO) scheme is selected for DFIG power regulation at the RSC. The SFO-based vector control strategy enables a decoupled regulation of the active and reactive powers that is flowing between the DFIG and the grid. The active and reactive powers' control strategy is performed by regulating the converter currents using Proportional-Integral (PI) controllers [5–8]. On the other hand, the Voltage Oriented Control (VOC) is adopted at the GSC, which the DC-link voltage is maintained to constant value [9].

Zhou and Blaabjerg [5] proposed a frequency-domain based approach to attain the PI controllers' gains for the inner and outer loops at both RSC and GSC circuits. The gains of the PI controllers that verify the design condition are selected for control loop. However, Hamane et al. [6] performed the direct power control of the DFIG using both PI and sliding mode controllers. The synthesis of the PI controller is based on an algebraic pole compensation method. Moreover, Hamane et al. [7] developed a comparative study of PI, RST, sliding mode and fuzzy supervisory controllers, and the synthesis of PI and polynomial RST controllers are tuned based on the pole compensation and pole assignment methods, respectively. In addition, the indirect vector control strategy is adopted, in which the matrix converter is used at the wound rotor of the DFIG instead of the conventional back-to-back converter in [8]. The active/reactive powers and current components are regulated by using the PI controllers that are designed using the pole placement technique. Although these methods have presented good performances, the main drawback for this type of control is that the performance of the DFIG system highly depends on a proper tuning of the PI controller gains. However, since the parameters of PI controllers depend on the precise mathematical representation, the control schemes become prone to error. In addition, the execution of tuning process can be time-consuming, and the optimal gains may not be obtained. Hence, proposing a systematic approach to find the best parameters setting of PI controllers is an interesting task and the metaheuristics-based hard optimization theory may be considered a promising solution.

In recent years, an enormous variety of metaheuristics optimization algorithms has been applied to solve complex and hard problems in science and engineering fields. Bekakra and Attous [10] presented the Particle Swarm Optimization (PSO) algorithm to obtain the optimal gains of PI controllers for the indirect control of the active and reactive powers' loops of a DFIG. Integral Absolute Error (IAE), Integral Time-weighted Absolute Error (ITAE) and Integral Square Error (ISE) performance criteria were selected as objective functions. Vieira et al. [11] applied the Genetic Algorithm (GA) to obtain the gains of the PI controllers at the RSC, where the active and reactive powers loops were directly controlled. In addition, an indirect power control strategy was achieved using GA-tuned integral sliding mode controllers, where the GA was used to tune the parameters of the direct and quadrature rotor current controllers [12]. Moreover, Assareh et al. [13] proposed a hybrid GA along with a gravitational search algorithm to attain the optimal gains of the PI controller for the torque regulation of a DFIG-based WT system. In all of these studies, the optimization algorithms are applied to obtain only the PI controllers' gains of the active and reactive powers at the RSC without discussion about the DC-link voltage regulation loop at the GSC and their e ffect on the overall control system.

Therefore, this work tends to apply a unified Thermal Exchange Optimization (TEO) metaheuristic algorithm to optimize the gains of PI controllers for the outer-loops in the classical vector control scheme of a DFIG-based wind energy system. In particular, this paper deals with the PI controllers tuning of the active and reactive power control loops in the RSC. It also treats the PI controller optimization-based design for the DC-link voltage loop in the GSC component.While the inner current loops at both RSC and GSC are tuned according to classical methods. The global TEO metaheuristic that has been proposed by Kaveh and Dadras in 2017 [14] is mainly adopted to tune the controllers' parameters of the modeled DFIG-based Wind Energy Converter System (WECS). Moreover, other well-known methods such as PSO, GA, Harmony Search Algorithm (HSA), Water Cycle Algorithm (WCA) and Grasshopper Optimization Algorithm (GOA) are used for an empirical comparison study. Moreover, a statistical analysis is performed to check the significance of each algorithm by using nonparametric tests such as Friedman's rank and Bonferroni–Dunn's test. The main contribution of this work is to provide a systematic and less complex procedure based on an advanced TEO algorithm in order to design and tune all outer-loops PI controllers for the well-known vector control scheme of a DFIG-based wind energy system. The classical trials-errors based methods of PI controller tuning are no longer used and the design time is further reduced. The drawbacks of the classical tuning methods are significantly reduced.

The remainder of this paper is organized as follows. The mathematical model of DFIG based WT is presented in Section 2. The description of the vector control scheme for the DFIG system at both RSC and GSC is analyzed in Section 3. In addition, Section 4 is devoted to the formulation of the outer-loops PI controllers' tuning problem given as a constrained optimization problem to be solved by the proposed TEO algorithm. In Section 5, such a TEO algorithm is described and its pseudo-code for the software implementation is given. Section 6 presents the implementation and the validation of the proposed TEO-tuned PI controller approach. Concluding remarks are given in Section 7.

#### **2. Modelling of the DFIG Based Wind Energy Converter**

The configuration of the DFIG-based WECS is depicted in Figure 1. A WT is joined to the DFIG by means of a gearbox. The DFIG is an induction machine, which the stator windings are directly connected to the grid, while the rotor windings are connected to the grid thanks to a back-to-back converter.

**Figure 1.** Configuration of the Doubly Fed Induction Generator (DFIG)-based Wind Energy Converter System (WECS).

#### *2.1. Modelling of the Wind Turbine*

The rotor blade of the WT is responsible for catching the wind power and converting it into kinetic energy. The captured mechanical power *Pm* is given as follows [2]:

$$P\_m = \frac{1}{2} \rho \mathbb{C}\_p(\lambda, \beta) \pi R^2 V\_w^3 \tag{1}$$

where ρ is the air density, *Cp* is the power coefficient and it depends on both the tip-speed ratio λ and the blade pitch angle β, *R* is the turbine radius and *Vw* is the wind speed.

The tip speed ratio is defined as the ratio of the blade tip speed to the wind speed. It is given by:

$$
\lambda = \frac{\Omega\_t R}{V\_w} \tag{2}
$$

where Ω*t* is the angular speed of the WT.

#### *2.2. Modelling of the DFIG*

The dynamical model of the studied DFIG can be expressed by the following electrical equations [6–8]:

$$\begin{cases} V\_{ds} = R\_s i\_{ds} + \frac{d\varphi\_{ds}}{dt} - \alpha\_s \wp\_{qs} \\\ V\_{qs} = R\_s i\_{qs} + \frac{d\varphi\_{qs}}{dt} + \alpha\_s \wp\_{ds} \\\ V\_{dr} = R\_r i\_{dr} + \frac{d\varphi\_{dr}}{dt} - (\alpha\_s - \omega\_m)\wp\_{qr} \\\ V\_{qr} = R\_r i\_{qr} + \frac{d\varphi\_{qr}}{dt} + (\omega\_s - \omega\_m)\wp\_{dr} \end{cases} \tag{3}$$

where *Vs* and *is* are the stator voltage and current, *Vr* and *ir* are the rotor voltage and current, ϕ*s* and ϕ*r* are the stator and rotor flux linkages, *Rs* and *Rr* are the stator and rotor resistances, ω*s* and ω*m* are the stator and rotor angular frequencies, respectively. The subscripts *d* and *q* denote the direct and quadrature axis components, respectively. The stator and rotor flux linkages are defined as follows:

$$\begin{cases} \begin{aligned} \varphi\_{ds} &= L\_{s}i\_{ds} + L\_{m}i\_{dr} \text{ and } \varphi\_{qs} = L\_{s}i\_{qs} + L\_{m}i\_{qr} \\ \varphi\_{dr} &= L\_{r}i\_{dr} + L\_{m}i\_{ds} \text{ and } \varphi\_{qr} = L\_{r}i\_{qr} + L\_{m}i\_{qs} \end{aligned} \end{cases} \tag{4}$$

where *Ls*, *Lr* and *Lm* are the stator, rotor and magnetizing inductances, respectively.

#### *2.3. Modelling of the GSC and the DC-Link*

The GSC is connected to the grid through an LCL-filter. However, for better understanding the control of GSC, it is necessary to describe the model of the grid side system and DC-link parts. The mathematical formulation of the GSC in the *dq* synchronous frame is defined by Equation (5). In this equation, *LT* represents the sum of the converter *Li* and grid side *Lg* inductances. In fact, the LCL-filter can be approximated to an L-filter equal to the sum of the LCL-filter inductors [4]:

$$\begin{cases} L\_T \frac{di\_{\mathcal{g}d}}{dt} = -R\_T i\_{\mathcal{g}d} + \omega\_{\mathcal{g}} L\_T i\_{\mathcal{g}q} + \varepsilon\_{\mathcal{g}d} - V\_{df} \\\ L\_T \frac{di\_{\mathcal{g}q}}{dt} = -R\_T i\_{\mathcal{g}q} - \omega\_{\mathcal{g}} L\_T i\_{\mathcal{g}d} + \varepsilon\_{\mathcal{g}q} - V\_{qf} \end{cases} \tag{5}$$

where *Vd f* and *Vq f* are the components of the converter side voltage, *igd* and *igq* are the components of the grid currents, *egd* and *egq* are the components of the grid voltage, and *RT* is the sum of the converter and grid side resistors.

The grid active and reactive powers are expressed as [9]:

$$\begin{cases} P\_{\mathcal{S}} = \frac{3}{2} (e\_{d\mathcal{S}} i\_{d\mathcal{S}} + e\_{\eta\mathcal{S}} i\_{\eta\mathcal{S}}) \\\ Q\_{\mathcal{S}} = \frac{3}{2} (e\_{\eta\mathcal{S}} i\_{d\mathcal{S}} - e\_{d\mathcal{S}} i\_{\eta\mathcal{S}}) \end{cases} \tag{6}$$

#### **3. Vector Control of the DFIG-Based Wind Energy Converter**

The classical direct power control of the DFIG system is divided into the RSC and GSC control loops, where the control structure for the RSC and GSC components consist of two cascaded control layers. The inner PI layers are adopted to regulate the components of rotor and grid currents. On the other hand, the outer PI ones in the RSC are implemented to control the active and reactive powers while the outer PI loop in the GSC is adopted to maintain the DC-link voltage to its reference value. However, in this work, the inner PI loops in the RSC and GSC circuits are designed based on a pole assignment technique [8,15], while the outer PI ones are tuned based on the proposed metaheuristics algorithms-based approach.

#### *3.1. Control of the RSC*

The RSC circuit is responsible to independently regulate the active and reactive powers in order to extract the maximum available power. The vector control by the SFO scheme is employed for the power regulation of the DFIG system.

Considering that the electrical network is stable, the stator flux ϕ*s* is constant. In addition, for the medium and high power DFIG, the value of the stator resistance is very small and therefore it can be neglected [6,7]. Hence, based on these assumptions, the stator voltages and fluxes in Equation (3) and Equation (4) can be rewritten respectively as follows:

$$\begin{cases} \ V\_{ds} = 0\\ V\_{qs} = V\_s = \omega\_s \wp\_{ds} \end{cases} \tag{7}$$

$$\begin{cases} \; \boldsymbol{\varrho}\_{\rm ds} = \boldsymbol{\varrho}\_{\rm ds} = \boldsymbol{L}\_{\rm s} \dot{\boldsymbol{i}}\_{\rm ds} + \boldsymbol{L}\_{m} \dot{\boldsymbol{i}}\_{\rm dr} \\ \; \boldsymbol{\varrho}\_{\rm qs} = \boldsymbol{0} = \boldsymbol{L}\_{r} \dot{\boldsymbol{i}}\_{\rm qs} + \boldsymbol{L}\_{m} \dot{\boldsymbol{i}}\_{\rm qr} \end{cases} \tag{8}$$

The active and reactive stator powers and rotor voltages are given, respectively, by:

⎪⎪⎪⎨

$$\begin{cases} P\_s = \frac{3}{2} \left( V\_{ds} i\_{ds} + V\_{qs} i\_{qs} \right) = -\frac{3}{2} \frac{L\_m}{L\_s} V\_s i\_{qr} \\\ Q\_s = \frac{3}{2} \left( V\_{qs} i\_{ds} + V\_{ds} i\_{qs} \right) = \frac{3}{2} \left( \frac{V\_s^2}{\alpha\_3 L\_a} - \frac{L\_m}{L\_s} V\_s i\_{dr} \right) \end{cases} \tag{9}$$

$$\begin{cases} V\_{dr} = R\_r \dot{i}\_{dr} + \sigma L\_r \underbrace{\frac{\Delta\_{dr}}{\frac{di\_{qr}}{dt}} - \overbrace{\sigma \omega\_r L\_r \dot{i}\_{qr}}^{\Delta\_{dr}}}\_{\Delta\_{qr}} \\\\ V\_{qr} = R\_r \dot{i}\_{qr} + \sigma L\_r \frac{d \dot{i}\_{qr}}{dt} + \sigma \omega\_r L\_r \dot{i}\_{dr} + \frac{L\_m}{L\_s} \omega\_r \dot{q}\_{ds} \end{cases} \tag{10}$$

where σ = *Lr* − *L*<sup>2</sup> *m*/*Ls* is the machine leakage coe fficient and ω*r* = ω*s* − ω*m* is the slip angular frequency.

As shown in Equation (9), the stator active and reactive powers are independent from each other. In addition, the components of stator power are linearly varying with the direct and quadrature rotor currents. This powers' regulation is thus performed using PI controllers for the *d* and *q* axis components of the rotor currents.

The outer control loops for the stator active and reactive powers are tuned based on the proposed metaheuristics-based procedure. Hence, the design of PI controllers for the inner rotor current loops is detailed according to Figure 2.

**Figure 2.** Proportional-Integral (PI) control scheme for the inner current loop of the RSC.

Based on a pole placement method as follows:

$$R\_{Pl}(s) = K\_{pir} \left( 1 + \frac{1}{T\_{\bar{\imath}ir\,\text{S}}} \right) \tag{11}$$

where *Kpir* and *Tiir* are the proportional gain and time constant of the rotor PI controllers, respectively.

The transfer function that describes the relationship between the voltages and currents dynamics is given as follows [8,15]:

$$R\_{\rm plant}(s) = \frac{i\_{rdq}}{V\_{rdq}} = \frac{1}{R\_r + \sigma L\_r s} \tag{12}$$

Since the same transfer functions for the direct and quadrature rotor currents are considered, the closed-loop transfer function *Rcli*(*s*) between the reference rotor current *i*∗*rdq* and the actual one *irdq* is written as:

$$R\_{cli}(s) = \frac{\dot{I}\_{rdq}}{\dot{I}\_{rdq}^\*} = \frac{K\_{prir}(T\_{iir}s + 1)}{T\_{iir}\sigma L\_r s^2 + T\_{iir}(R\_r + K\_{pir})s + K\_{prir}} \tag{13}$$

By equating the denominator of the closed-loop transfer function of Equation (13) under the general form of a second order system, the parameters of PI rotor current controllers can be found as follows [8,15]:

$$\begin{cases} \begin{array}{c} \mathcal{K}\_{\text{pir}} = 2\mathcal{\xi}\omega\iota\_{\text{n}}\sigma L\_{r} - R\_{r} \\ T\_{\text{iir}} = \frac{2\mathcal{\xi}}{\omega\iota\_{\text{n}}} - \frac{R\_{r}}{\alpha\omega\iota\_{\text{n}}^{2}L\_{r}} \end{array} \tag{14}$$

where ξ and ω*n* are the damping coefficient and natural frequency of the desired closed-loop reference model, respectively.

For the design purpose, the discussion is made with respect to the choice of the damping coefficient ξ and natural frequency ω*n*. Since the inner current loop in the cascade control scheme will have a much larger bandwidth than the one used in the outer-loop, a regulation at 1/20 of the switching frequency *fsw*,*RSC* is retained. For the outer loop, a regulation between 1/20 and 1/10 of the inner loop bandwidth is made [5]. Then, feed forward compensations Δ*dr* and <sup>Δ</sup>*qr* should be added back to generate the desired rotor voltages *<sup>V</sup>*<sup>∗</sup>*dr*and *<sup>V</sup>*<sup>∗</sup>*qr*.

#### *3.2. Control of the GSC*

The VOC strategy is applied to control the GSC, which usually contains one outer PI control loop that regulates the DC-link voltage regardless of the magnitude and the direction of the rotor power. Two inner PI current loops that regulate the direct and quadrature gird currents are also included [16,17]. In addition, the passive damping strategy is employed to mitigate the resonance problem. The mathematical representation of GSC is built on the approximated model of the LCL-filter, which is used to tune the PI grid current controllers. Whereas, the whole transfer function of the LCL-filter with passive damping is taken into account for the stability analysis and investigation [4]. The block diagram of the grid current loops is shown in Figure 3.

**Figure 3.** Design of the PI controller current loop of grid side converter.

To implement the VOC scheme, the *d*-axis is aligned with the grid vector voltage. Therefore, this leads to a *d*-axis grid voltage equal to its magnitude, and the *q*-axis voltage which is then equal to zero. Hence, the grid power expressions can be expressed as [16,17]:

$$\begin{cases} P\_{\mathcal{S}} = \frac{3}{2} e\_{d\mathcal{S}} i\_{d\mathcal{S}} \\ Q\_{\mathcal{S}} = -\frac{3}{2} e\_{d\mathcal{S}} i\_{d\mathcal{S}} \end{cases} \tag{15}$$

Equation (15) decides the active power and consequently the DC-link voltage is controlled via the direct current whereas the quadrature current is used to regulate the reactive power that exchanges with the grid [17]. With the VOC approach, the dynamic equations of the grid currents are rewritten as follows:

$$\begin{split} c\_{d\xi} &= R\_T i\_{d\xi} - \overbrace{\alpha\_{\mathcal{S}} L\_T i\_{q\mathcal{S}}}^{\Delta\_{d\xi}} + L\_T \frac{d i\_{d\xi}}{dt} + V\_{df} \\ 0 &= R\_T i\_{q\mathcal{S}} + \overbrace{\alpha\_{\mathcal{S}} L\_T i\_{d\xi}}^{\Delta\_{d\xi}} + L\_T \frac{d i\_{q\mathcal{S}}}{dt} + V\_{qf} \end{split} \tag{16}$$

By applying the same methodology of PI rotor current controllers design in the RSC, the gains of PI controllers for the grid currents dynamics are found as [15,18]:

$$\begin{cases} \begin{cases} \mathcal{K}\_{\text{pig}} = 2\xi \omega\_n L\_T - R\_T\\ \mathcal{K}\_{\text{iig}} = \omega\_n^2 L\_T \end{cases} \tag{17} $$

where ω*n* = 847.80 rad/ sec and ξ = 0.707.

#### **4. PI Controllers Tuning Problem Formulation**

In the PI control framework, appropriate values of *Kp* and *Ki* gains are generally obtained by empirical methods and trials-errors based procedures [19]. These nonsystematic and challenging tasks become more di fficult and time-consuming, especially for the complex and large-scale systems like the studied DFIG-based WECS. So, the idea of formulating the *Kp* and *Ki* gains' selection as an optimization problem is a promising solution. Such a control problem can be nonlinear, non-smooth or even non-convex and can be e ffectively solved thanks to advanced metaheuristics [10,11]. In this work, three PI controllers for the outer-loops at both RSC and GSC components are considered for the optimization process. These PI controllers for the DC-link voltage, active and reactive powers' dynamics are systematically tuned thanks to the proposed TEO metaheuristic which is described in the following. Figure 4 gives the proposed optimization-based tuning scheme of the PI controllers for the DFIG-based WECS.

**Figure 4.** The proposed control algorithm applied to the DFIG system.

The decision variables of the formulated optimization problem are obviously the proportional and integral gains of PI controllers for the active/reactive powers and DC-link loops. They are specified for the optimization process as follows:

$$\mathbf{x} = \begin{bmatrix} \mathbf{K}\_{\mathbf{p}\mathbf{P}\_{s'}} \mathbf{K}\_{\mathbf{i}\mathbf{P}\_{s'}} \mathbf{K}\_{\mathbf{p}\mathbf{Q}\_{s'}} \mathbf{K}\_{\mathbf{i}\mathbf{Q}\_{s'}} \mathbf{K}\_{\mathbf{p}\mathbf{d}c'} \mathbf{K}\_{\mathbf{i}\mathbf{d}c} \end{bmatrix}^{\mathrm{T}} \in \mathbf{S} \subseteq \mathbb{R}\_{+}^{6} \tag{18}$$

where S = *x* ∈ <sup>R</sup>6+,*xlow* ≤ *x* ≤ *<sup>x</sup>up*denotes the bounded search space for all PI parameters.

⎧⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎨⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎩

However, the defined objective functions are minimized taking into account a scope of time-domain restrictions. These are concerning to the maximum overshoot δmax, steady-state error *Ess*, rise time *tr* and/or settling time *ts* of the closed-loop system step-response [20]. So, the tuning issue associating the PI controllers of the DFGI-based WECS can be expressed as follows:

$$\begin{aligned} & \text{Minimize } f\_m(\mathbf{x}), \; m \in \{IAE, ISE, ITSE, ITE\} \\ & \mathbf{x} = \begin{bmatrix} K\_{pP\_s}, K\_{iP\_s}, K\_{pQ\_s}, K\_{iQ\_s}, K\_{pdc}, K\_{idc} \end{bmatrix}^T \in \mathcal{S} \subseteq \mathbb{R}\_+^6 \\ & \text{subject to }: \\ & g\_1(\mathbf{x}) = \delta\_{P\_s} - \delta\_{P\_s}^{\max} \le 0 \\ & g\_2(\mathbf{x}) = \delta\_{Q\_s} - \delta\_{Q\_s}^{\max} \le 0 \\ & g\_3(\mathbf{x}) = \delta\_{dc} - \delta\_{dc}^{\max} \le 0 \\ & K\_{pj, \min} \le K\_p \le K\_{pj, \max} \\ & K\_{ij, \min} \le K\_i \le K\_{jf, \max}, j \in \{P\_s, Q\_s, dc\} \end{aligned} \tag{19}$$

where *fm* : <sup>R</sup>6+ → R are the cost functions, *gq* : <sup>R</sup>6+ → R are the problem's inequality constraints, δ*dc*, δ*Ps* and <sup>δ</sup>*Qs* are the overshoots of the controlled DC-voltage, active and reactive powers dynamics, respectively, δmax *dc*, δmax *Ps*and δmax *Qs*denote their maximum given value.

*Math. Comput. Appl.* **2019**, *24*, 77

Since the optimization problem (19) is a multi-objective type, weighted sum-based aggregation functions are used as follows to handle with the multiple costs [21]:

$$f\_{\rm IAE}(\mathbf{x}) = w\_{P\_s} \int\_0^T \left| e\_{P\_s}(\mathbf{x}, t) \right| dt + w\_{Q\_s} \int\_0^T \left| e\_{Q\_s}(\mathbf{x}, t) \right| dt + w\_{dc} \int\_0^T \left| e\_{d\mathbf{c}}(\mathbf{x}, t) \right| dt \tag{20}$$

$$f\_{\rm ISE}(\mathbf{x}) = w\_{\rm P\_s} \int\_0^T e\_{\rm P\_s}^2(\mathbf{x}, t)dt + w\_{\rm Q\_s} \int\_0^T e\_{\rm Q\_s}^2(\mathbf{x}, t)dt + w\_{\rm dc} \int\_0^T e\_{\rm dc}^2(\mathbf{x}, t)dt\tag{21}$$

$$f\_{ITSE}(\mathbf{x}) = w p\_s \int\_0^T t e\_{P\_s}^2(\mathbf{x}, t) dt + w\_{Q\_s} \int\_0^T t e\_{Q\_s}^2(\mathbf{x}, t) dt + w\_{dc} \int\_0^T t e\_{dc}^2(\mathbf{x}, t) dt \tag{22}$$

$$f\_{ITAE}(\mathbf{x}) = w\_{P\_s} \int\_0^T t \Big| e\_{p\_s}(\mathbf{x}, t) \Big| dt + w\_{Q\_s} \int\_0^T t \Big| e\_{Q\_s}(\mathbf{x}, t) \Big| dt + w\_{dc} \int\_0^T t \Big| e\_{dc}(\mathbf{x}, t) \Big| dt \tag{23}$$

where *T* denotes the total simulation time, *wPs* > 0, *wQs* > 0 and *wdc* > 0 are the weighting coefficients of the aggregation functions satisfying *wPs* + *wQs* + *wdc* = 1, and *ej*(.), *j* ∈ {*Ps*, *Qs*, *dc*} denote the tracking errors between the plant output and the relative set-point values, i.e., *ePs*(*<sup>x</sup>*, *t*) = *P*∗*s* − *Ps*(*<sup>x</sup>*, *t*), *eQs*(*<sup>x</sup>*, *t*) = *Q*<sup>∗</sup>*s* − *Qs*(*<sup>x</sup>*, *t*) and *edc*(*<sup>x</sup>*, *t*) = *<sup>V</sup>*<sup>∗</sup>*dc* − *Vdc*(*<sup>x</sup>*, *t*).

The judgment matrix method is used to determine the weighting coefficients of the objective functions [22,23]. Such a method grades all objective functions based on the importance of each one. After calculating the eigenvalues of such a matrix, the weighting coefficients of functions (20)–(23) can be chosen as *wPs* = 0.6370, *wQs* = 0.2583 and *wdc* = 0.1047.

#### **5. Thermal Exchange Optimization Algorithm**

The Thermal Exchange Optimization (TEO) algorithm is a novel metaheuristic inspired by Newton's law of the cooling [3]. In this population-based metaheuristic, each agen<sup>t</sup> is modeled as a cooling object. By associating another agen<sup>t</sup> as a surrounding fluid, a heat transfer and thermal exchange occurs between them. The new temperature of each agen<sup>t</sup> is considered as a new position of the potential solution in the search space.

In a *d*-dimensional search space and at the *kth* iteration time, each object of the population is described by its temperature *Tki* = *Tki*,1, *<sup>T</sup>ki*,2, ... , *<sup>T</sup>ki*,*<sup>d</sup>*, (*i*, *k*) ∈ 1, *Npop* × [0, *Niter*]. In order to enhance the optimization performances, the TEO algorithm uses a Thermal Memory (TM) to store a number of historically best vectors as well as their related fitness. Therefore, these stored solution vectors are added to the population and the same numbers of current worst objects are deleted. After that, a growing order of solutions is retained according to their cost function values. The *Npop* sorted objects are equally divided into two groups of environment and cooling objects as shown in Figure 5. The environment objects are denoted as *Tk*1 , *Tk*2 , ... , *<sup>T</sup>kNpop*/2 while the cooling ones are *<sup>T</sup>kNpop*/2+<sup>1</sup>, *<sup>T</sup>kNpop*/2+<sup>2</sup>, ... , *TkNpop* .

**Figure 5.** Pairs of environment and cooling objects.

Since the object has lower η can exchange the temperature slightly, the TEO metaheuristic proposed a similar formula to evaluate the value of η for each object as described in Equation (24):

$$\eta = \frac{\text{Cost (object)}}{\text{Cost(worst object)}} \tag{24}$$

Over the course of iterations, the value *t* = *k*/*Niter* increases as follows leading to improving the exploration mechanism [3]. In order to enhance the TEO exploration capacity, Equation (25) has been proposed to avoid the trapping in local optima and update the environmental temperature as follows:

$$T\_{env,i}^{k+1} = (1 - \mathcal{U}|0,1\rangle(c\_1 + c\_2(1-t)))T\_{env,i}^k \tag{25}$$

where *c*1 and *c*2 are the controlling variables chosen as 0 or 1, U{0, 1} is a uniformly random number, *T<sup>k</sup> env*,*i* is the previous environmental temperature, (1 − *t*) is proposed to reduce the randomness by closing to the last iterations. When the decreasing of the randomness is not required, the value of *c*2 is set to zero.

The new temperatures for each object, i.e., either cooling objects or environmental ones, are updated as follows:

$$T\_i^{k+1} = T\_{env,i}^k + \left(T\_i^k - T\_{env,i}^k\right)e^{-\eta t} \tag{26}$$

In the TEO formalism, another formulation is suggested to improve the exploration ability. A control parameter 0 ≤ *pro* ≤ 1 is introduced and determines whether a component of each cooling object must be changed or not. For each object, such a parameter is compared to *Rand*(*i*),*<sup>i</sup>* = 1, 2, ... , *Npop*, which is a random number uniformly distributed between 0 and 1. If *Rand*(*i*) < *pro*, one dimension of the *ith* object is chosen randomly and its value is regenerated as follows [3,24]:

$$T\_{\left(i,j\right)} = T\_{j,\text{min}} + \mathcal{U}(0,1) \times \left(T\_{j,\text{max}} - T\_{j,\text{min}}\right) \tag{27}$$

where *<sup>T</sup>*(*<sup>i</sup>*,*j*) is the *jth* component of the *ith* object, *Tj*,max and *Tj*,min are the upper and lower bounds of the *jth* component, respectively.

Finally, the steps of the proposed TEO algorithm are summarized as follows:

**Step 1.** Randomly initialize the temperature for all objects *T*0*i*, *i* = 1, 2, ... , *Npop*.

**Step 2.** Calculate the fitness of each search object.

**Step 3.** Save some *T* best vectors and their related cost values in the TM.

**Step 4.** Add the saved solutions and remove the same numbers of the worst objects.

**Step 5.** Arrange the objects according to their related fitness in an ascending order.

**Step 6.** Divide the objects into two equal groups: environment and cooling objects.

**Step 7.** Calculate the parameters η and *t*.

**Step 8.** Change the environment temperatures by Equation (25).

**Step 9.** Update the temperatures according to Equations (26) and (27).

**Step 10.** Check the termination criterion and repeat the iterations.

To describe the proposed metaheuristics-based tuning strategy of the PI controllers in the studied DFIG-based energy conversion system, a detailed flowchart is shown in Figure 6.

**Figure 6.** Flowchart of the proposed Thermal Exchange Optimization (TEO)-tuned PI controllers' parameters.

#### **6. Simulation Results and Discussions**

#### *6.1. Execution of the Metaheuristic Algorithms*

The proposed metaheuristics-based direct power control of the studied DFIG-based WECS is built using the MATLAB/Simulink environment. The respective nominal parameters of the DFIG, power converters at both sides, L-filter and LCL-filter are presented in Table 1. Specifically, this work deals with the PI controllers tuning for the active and reactive powers loops in the RSC circuit, as well as it treats with the PI controller optimization-based selection of the DC-link voltage loop in the GSC component. It is worth indicating that the inner current loops in the RSC and GSC are tuned according to the pole assignment method. Since the grid current closed-loop of the LCL-filter is unstable, the passive damping method is applied as a solution to ensure the stability of such current dynamics. To investigate the e ffectiveness and superiority of the proposed TEO method, the well-known global algorithms GA, PSO, HSA, WCA and GOA are considered for comparison purposes.



All reported algorithms are independently run 10 times. The termination criterion is set as a maximum number of iteration reached *Niter* = 100 for a population size of the thermal agents equal to *Npop* = 50. All the values of the used common parameters are kept equal. All algorithms have been executed on a PC computer with Core TM i5-7200U CPU and 2.5 GHz/8.00 GB RAM. The specific control parameters of the reported algorithms are listed in Table 2.



The optimization problem (19) is minimized for various performance criteria, i.e., IAE, ISE, ITAE and Integral Time-weighted Square Error (ITSE), and under time-domain operational constraints. These performance indexes are calculated within two conditions, which are considered as follows:



To find the accurate values of the decision variables *x* = *KpPs*,*KiPs*,*KpQs*,*KiQs*,*Kpdc*,*Kidc* , the bounded search domain is initially set to 0 ≤ *x* ∈ R<sup>6</sup> + ≤ 105. Several runs of the algorithm are performed within this initial limitation. The results indicated that the search space could be reduced to 0 ≤ *x* ∈ R<sup>6</sup> + ≤ 400. Therefore, the TEO algorithm explores within a smaller search domain for the next runs and more precise solution output vector is achievable.

Table 3 gives the statistical results attained by the introduced algorithms under minimizing the cost functions described in Equations (20)–(23). The worst, average, best values and the standard deviation (STD) experimental results are summarized in the Table 3, where the optimal mean value of each index is highlighted in bold and underlined. In addition, the Elapsed Time (ET) is also listed, which is defined as the time that any algorithm requires finding the best solution. It can be clearly observed that the proposed TEO produces very competitive solutions with the reported algorithms. Appendix A summarizes the obtained gains of the PI controllers for each proposed optimization method. Indeed, these tuned PI controllers' gains lead to the best transient and steady state responses of the entire reported algorithms.


**Table 3.** Statistical results of optimization problem (19) over 10 independent runs.

In addition, Figure 7 shows the convergence histories of the mean objective function values of IAE, ISE, ITAE and ITSE criteria, respectively. Indeed, it is shown that the proposed TEO metaheuristic for both IAE and ISE indices outperforms the other reported methods in terms of the fastness and non-premature convergence as well as the solutions' quality. The TEO-based method for ITAE and ITSE criteria gives the best solution as a second order after the GA one. From these results, the superiority of the TEO metaheuristic is still shown in terms of exploitation and exploration capabilities for local and global searches. This further justifies the use of such a global metaheuristic to systematic and easy design of the direct power control of the DFIG-based WECS. The Box-and-Whisker plots of the mean objective function values are shown in Figure 8. From these results, one can observe that the boxplots of the TEO algorithm for all given performance criteria are generally lower and narrower than other algorithms. This confirms the high reproductivity of the TEO algorithm toward finding the optimal values of the solution.

In the remainder of the results, the obtained gains of the optimized PI controllers are used to assess the frequency-and time-domains performances of the outer-loops at the RSC and GSC power components. Moreover, comparisons with the classical pole placement [8,15], frequency response [5] symmetrical optimum [27], Ziegler–Nichols [28] and Tyreus–Luyben [29] tuning methods are made for the DC-link voltage dynamics as shown in Table 4 and Figure 9. The closed-loop response of the DC-link voltage dynamics is investigated under variable voltage profile. The DC-link voltage reference is varied up and down at various step levels. Figure 9 describes such a response around the final set-point value for di fferent tuned PI controllers. Although, in real systems, the DC-link voltage is required to be constant. This scenario is proposed to check the capability of the proposed TEO algorithm under di fferent circumstances. The aim is to show the di fference between the classical tuning methods and the proposed optimization-based one. Referring to this result, the ITAE-based TEO algorithm can regulate the DC-link voltage dynamics with higher performance compared to the other algorithms.

**Figure 7.** Convergence rates comparison: (**a**) IAE criterion; (**b**) ISE criterion; (**c**) ITAE criterion; (**d**) ITSE criterion.

**Figure 8.** Box-and-Whisker plot of optimization problem (19): (**a**) IAE criterion; (**b**) ISE criterion; (**c**) ITAE criterion; (**d**) ITSE criterion.


**Table 4.** Time-domain performances for controlled DC-link voltage under step changes scenario.

**Figure 9.** DC-link voltage responses against step changes for different-tuned PI controllers.

The voltage tracking is achieved with best performance in terms of response precision and fastness, i.e., the steady-state error and rise/settling time metrics are minimal as depicted in Table 4. The performance in terms of response damping is acceptable for the proposed TEO-tuned PI control as shown in Figure 9. Referring to the numerical results in Table 4, one can clearly observe that the proposed TEO algorithm showed the best rise/settling times and steady-state error except the overshoot index, which the Tyreus–Luyben method gained the minimum overshoot. Since several runs were executed to obtain the PI controllers' gains with the classical Ziegler–Nichols and Tyreus–Luyben methods, the tuning process becomes a tedious and time-consuming task. These minor degradations of the overshoot performance do not influence the effectiveness of the proposed TEO-based tuning method with respect to the systematization of the synthesis procedure, the simplicity of implementation and the superiority in other performance indices. The proposed TEO algorithm found the optimal gains of the PI controllers within a reduced computation time and under operational constraints compared to the reported trials-errors based tuning procedures which generally give local solutions for the formulated control problem.

Regarding to the grid current dynamics, the same tuning methods with adding the technical optimum instead of the symmetrical ones are used to select and adjust the PI controllers' gains. All obtained results are reported in Table 5.


**Table 5.** Time-domain performances for controlled quadrature grid current under a step change scenario.

The control performance of the introduced tuning methods is investigated by supposing that there is a step change in the q-axis grid current reference at time *t* = 0.2 sec. The reference q-axis grid current is changed from zero to −0.2 pu. However, the *q*-axis reference grid current is usually set at zero to achieve a unity power factor. Figure 10 shows the quadrature grid current response under step changes. It can be observed that the Tyreus–Luyben method gives superior performance in comparison with the other reported methods. Here, it is important to mention that the gains of the PI controllers for the inner current loops at both RSC and GSC are selected thanks to the classical methods and not with the proposed optimization algorithms. Since the plant models of the inner loops are available, the use of the pole-placement, Tyreus–Luyben and Ziegler–Nichols methods is well adapted instead of metaheuristics-based methods which do not require models of systems to be controlled. The Bode plot of the grid current loop is presented in Figure 11.

**Figure 10.** Quadrature grid current responses against a step change for di fferent-tuned PI controllers.

**Figure 11.** Bode plot of the open loop grid current loop: with and without passive damping based technical optimum tuning.

The quality of the three phase grid currents is investigated through measuring the THD of the AC grid currents. As depicted in Figure 12. It is obviously shown that the use of the LCL-filter achieves a better attenuation with a THD value about 4.05% for which the IEEE519-1992 standard limits are respected. Further analysis can be made to compare the THD of the AC grid currents in the case of using an L-filter based structure. According to the results of Figure 12, it is clear that the THD for adopting the LCL-filer is smaller than that using the L-filer type at the same value of inductance.

**Figure 12.** Harmonics spectrum of the AC grid currents.

Referring to the RSC control results, the stator active and reactive powers can be regulated by controlling the q-axis rotor current. The reference value of stator active power is generated from a Maximum Power Point Tracking (MPPT) strategy to extract the maximum power from the wind. In addition, the value of reference reactive power is set at di fferent step levels to check the performance of the proposed PI controllers. As shown in Figure 13, both active and reactive powers track e ffectively their reference values with good performance in terms of speed and damping dynamics. In addition, the decoupling between the active and reactive powers is perfectly assured and the DFIG extracts the maximum available power, which is approximately about 1.5 MW. Figure 14 demonstrates the control performance of the PI controller for the reactive power loop around its final set-point value based on di fferent tuned PI methods. The aim is to show the di fference between the optimization tuning-based methods. On the other hand, both left and right sides of Figure 15 show the tracking performance of the direct and quadrature rotor currents, respectively. Such tracking dynamics are perfectly performed and improved thanks to the proposed TEO-based tuning and control method.

**Figure 13.** Time-domain variations of the TEO-tuned PI controlled active and reactive powers.

**Figure 14.** Performance comparison of the controlled reactive power under different PI controllers' tuning methods.

**Figure 15.** Time-domain variations of the TEO-based control of the direct and quadrature rotor currents.

#### *6.2. Statistical Analysis and Comparison*

In this part, the mean execution values related to the different optimization criteria are sorted to assess the best operating one according to its average objective function performance. Moreover, a statistical comparison based on the nonparametric Friedman and Bonferroni–Dunn is carried out by using these mean performances [30,31]. The average ranks for all the proposed methods based on the four performances indices are provided in Table 6. One can note that the proposed TEO metaheuristic has worthily attained the lowest average ranks compared to the remaining methods. A statistical analysis has been performed to highlight the importance of the TEO-based tuning method over other algorithms [32].


**Table 6.** Average Rank based statistical analysis of mean performances.

The Friedman test for six competitor algorithms (*m* = 6) and four indices (*l* = 4) provides the computed value χ2*F* = 12.14 of the χ-distribution. The critical value of such a distribution with degrees of freedom *m* − 1 = 5 and at confidence level α = 0.05 is equal to χ25,0.05 = 11.07. Since the above computed score is greater than this statistical value, the null hypothesis is declined. Moreover, for the Iman–Davenport test [30,31], the statistic is distributed with *m* − 1 = 5 and (*m* − 1) × (*l* − 1) = 15 degrees of freedom. For this test, the null hypothesis is also rejected as the calculated value *FF* = 4.64 is greater than the critical value *F*5,15,0.05 = 2.9 at the same significant level of confidence. All these statistical results indicate that there are significant differences among the performances of the reported algorithms for the optimization problem (19). Thus, the proposed TEO is found to be the most effective one having the best average Friedman ranking as given in Table 6. However, a Bonferroni–Dunn post-hoc test is called to investigate whether or not the proposed TEO algorithm is significantly better than another algorithm at the above considered level of confidence [32]. The corresponding Critical Differences (CD) of the reported algorithms at the confidence levels α = 0.05 (95% significance level) and α = 0.1 (90% significance level) are computed as CD0.05 = 3.17 and CD0.1 = 2.88, respectively. Figure 16 illustrates the graphical representation of the Bonferroni–Dunn test considering the TEO as the control algorithm [32,33].

**Figure 16.** Graphical representation of Bonferroni–Dunn's test for problem (19).

The bar of the TEO algorithm is the lowest high among the reported algorithms and the heights of the bars corresponding to the PSO algorithm. The GOA method violates the horizontal lines of significant levels. This reveals that the TEO algorithm performs at least significantly better than these two algorithms over the solutions equality.

#### *6.3. Computational Time E*ffi*ciency*

In this subsection, the average elapsed times of Table 3 are used to assess the computational efficiency of the reported algorithms. Such an algorithmic property of resource usage is quantified by the Computational Time Efficiency (CTE) metric defined as follows [34]:

$$\text{CTE} = \frac{\text{ETI}\_{\text{algorithm}}}{\text{ETI}\_{\text{total}}} (\text{\textdegree) }\tag{28}$$

where *ETalgorithm* and *ETtotal* denote the average and total elapsed times for the same index.

According to the ET measures of problem (19), the CTE for each algorithm over the reported optimization criterion is summarized in Table 7.



In terms of computational e fficiency, it can be observed from these results that the proposed TEO algorithm attained the second rank for IAE and ITAE criteria and the third and fourth ranks for ISE and ITSE criteria, respectively. Roughly, the HSA algorithm presents the best computational e fficiency in terms of time resource usage. The GOA metaheuristic achieved the best average elapsed time as shown in Figure 17.

**Figure 17.** Average elapsed time over all criteria for the introduced algorithms.

#### *6.4. Sensitivity Analysis*

In the TEO formalism, the parameters *c*1 and *c*2 are introduced to reduce the randomness of the algorithm over the course of iterations and to improve the exploration capabilities. The size of thermal memory TM is considered as a predefined parameter, which can enhance the performances of the algorithm without increasing the computational cost [3,24]. However, *pro* is a user parameter that improves the global search capacity. For illustration purposes, the resolution of problem (19) is tested by changing these above main parameters. As shown in Table 8, the case with the parameters set (*c*1 = 1, *c*2 = 1) shows the superiority of the optimization results compared to other scenarios. In addition, it can be seen from Table 9 that the user parameter *TM* = 10 leads to obtain the best solutions. Finally, Table 10 investigates the e ffect of the parameter *pro* on the algorithm performance. The convergence of the TEO algorithm is always guaranteed under all these control parameters variations which shows the insensitivity of the proposed physics-inspired metaheuristic.


**Table 8.** Comparison under different values of *c*1 and *c*2: IAE criterion.

**Table 9.** Comparison under different values of *TM*: IAE criterion.


**Table 10.** Comparison under different values of *pro*: IAE criterion.


## **7. Conclusions**

This paper proposes an intelligent metaheuristics-based design procedure to tune the outer-loop PI controllers for the direct power control scheme of DFIG-based energy conversion systems. In this study, only outer PI loops in the RSC and GSC circuits are optimized for the active and reactive powers and DC-link voltage dynamics. Since the reported classical techniques of PI controller tuning are tedious, time-consuming and not systematic, a TEO-based approach has been proposed and successfully implemented. The PI controllers tuning problem is firstly formulated as a constrained optimization program under nonlinear and non-smooth operational constraints. The introduced TEO algorithm is then employed separately to minimize several time-domain performance criteria such as IAE, ISE, ITAE and ITSE indices as objective functions. The proposed TEO-tuned PI controllers methodology improves the performance and robustness of the controlled DFIG-based energy converter in terms of rising time, settling time, and steady state indices. The classical trials-errors based methods of PI controllers tuning are no longer used and the design time is further reduced. To evaluate the performance superiority of the proposed TEO-based approach in finding the global minimum value of the objective function for various performance indices, a comparison study with the PSO, GA, HSA, GSO and WCA is performed. The demonstrative results exhibit that the proposed TEO gives very complete results in terms of global search capabilities, robustness and non-premature convergence. Finally, the statistical analysis is achieved by using the Friedman's rank and Bonferroni-Dunn's test. The corresponding results show that the proposed TEO-based method is a promising alternative approach for controlling the DFIG system by systematically tuning the unknown PI controllers' parameters efficiently.

**Author Contributions:** Conceptualization, M.M.A. and S.B.; Methodology, S.B.; Software, M.M.A.; Validation, M.M.A.; Formal analysis, S.B.; Investigation, S.B.; Resources, M.M.A. and S.B.; Writing—original draft preparation, M.M.A.; Writing—review and editing, S.B.; Visualization, S.B.; Supervision, S.B.

**Funding:** This research received no external funding.

**Acknowledgments:** The authors would like to express their sincere appreciation to the anonymous reviewers for their valuable and constructive comments, which greatly improve the quality of this paper.

**Conflicts of Interest:** The authors declare no conflict of interest.
