A Reproducible Benchmarking Framework for Differential Evolution Variants in Wind Turbine Controller Tuning

1. Introduction

The regulation of nonlinear dynamical systems operating under varying conditions remains a fundamental challenge in control engineering. Strong nonlinearities, coupled subsystems, and operating-region-dependent dynamics often limit the effectiveness of analytical tuning rules and linearization-based design methodologies. In such environments, controller performance may degrade significantly when disturbance rejection, operating-point variations, and multi-objective requirements must be simultaneously addressed.

These challenges are particularly relevant in modern energy systems, aerospace applications, and complex electromechanical plants, where high-performance regulation must be maintained despite nonlinear interactions between subsystems and varying environmental conditions [1,2]. In such scenarios, achieving robust disturbance rejection, stability across multiple operating regimes, and satisfactory transient behavior simultaneously becomes a demanding control-design task.

Large-scale wind turbines constitute a representative example of nonlinear multi-regime systems. Their dynamics arise from the interaction of aerodynamic, mechanical, and electrical subsystems, leading to strong nonlinear behavior and operating-region-dependent control objectives [3,4,5,6]. In particular, during high-load operation, control objectives are primarily associated with disturbance rejection around regulated operating points rather than classical reference tracking. These characteristics increase tuning complexity and reduce the effectiveness of conventional linear design tools. Consequently, the systematic regulation of wind turbines provides a demanding and realistic context for evaluating advanced controller-tuning methodologies [7,8].

Despite the nonlinear nature of wind turbine dynamics, industrial control systems frequently rely on relatively simple feedback controllers due to their ease of implementation, robustness, and interpretability. Among these, the Proportional–Integral–Derivative (PID) controller remains the dominant industrial standard. Its widespread adoption is largely attributed to its intuitive structure, low implementation cost, and the existence of well-established tuning rules such as the Ziegler–Nichols and Cohen–Coon methods [9].

However, the structural simplicity of PID controllers can limit achievable performance in systems exhibiting high-order dynamics, significant nonlinearities, or strong coupling between subsystems. In such situations, the limited number of adjustable parameters may restrict the ability to shape closed-loop dynamics effectively. To overcome these limitations, several extensions of the PID structure have been proposed in the literature to increase design flexibility and improve dynamic compensation capabilities.

Among these extended structures, the Proportional Integral Derivative Accelerative (PIDA) controller, originally introduced in [10], incorporates an additional zero in the controller transfer function through the inclusion of an accelerative term. Compared with the conventional PID structure, this additional term increases the degrees of freedom available for shaping the closed-loop dynamics. In this sense, PIDA and related higher-order derivative controllers have been reported to improve transient response and disturbance-rejection capability in processes with demanding dynamic characteristics, such as high-order, delayed, or oscillatory plants [11,12].

From a physical perspective, this feature is relevant to wind turbine control because rapid wind-induced variations can lead to transient speed deviations and oscillatory responses in the regulated system. The accelerative action may therefore provide additional dynamic compensation, potentially contributing to reduced overshoot, faster regulation around the operating point, and improved disturbance attenuation. This is especially important in above-rated operation, where pitch control is the main mechanism for regulating generator speed and mitigating aerodynamic disturbances [13].

However, the increased dynamic compensation capability of the PIDA controller also introduces additional tuning challenges. If the accelerative term is not properly filtered or tuned, it may increase sensitivity to high-frequency components, measurement noise, or aggressive control actions, which can affect closed-loop robustness and actuator effort [14,15,16]. Therefore, the practical use of PIDA controllers requires a systematic tuning strategy capable of balancing performance improvement, stability preservation, and robustness under different operating conditions. As a result, PIDA controllers have attracted attention in applications where enhanced dynamic compensation is desirable [17,18,19]. Nevertheless, unlike PID controllers, PIDA structures lack consolidated and systematically validated tuning methodologies, particularly in nonlinear multi-regime environments. This absence of structured tuning approaches motivates the exploration of optimization-based methods.

In recent decades, metaheuristic optimization algorithms have emerged as effective alternatives for controller tuning in nonlinear and high dimensional systems [20,21,22]. Among these techniques, Differential Evolution (DE) has received considerable attention due to its algorithmic simplicity, robustness, and strong performance in continuous search spaces [23,24,25]. Since its original formulation, DE has evolved into a broad family of variants derived from modifications to mutation and crossover strategies. These variants exhibit distinct exploration–exploitation characteristics and convergence behaviors, leading to performance differences across optimization problems [26,27,28,29]. Additionally, sensitivity to manually tuned control parameters in classical variants has motivated the development of adaptive formulations such as CODE, JADE, jDE, and SaDE [30,31], as well as empirical investigations into parameter-setting strategies [32].

Despite the widespread application of DE in controller tuning, existing comparative studies often evaluate only limited subsets of algorithm variants or rely on mathematical benchmark functions that do not fully represent the dynamic complexity of nonlinear control systems. Furthermore, differences in experimental protocols—including random initialization, simulation environments, cost-function definitions, search-space settings, and implementation procedures—can hinder reproducibility and complicate direct comparison among reported results [33,34,35]. Consequently, a methodological gap remains in the availability of controlled and reproducible benchmarking procedures for systematically evaluating a broad set of classical DE variants in realistic nonlinear controller-tuning problems. Under heterogeneous experimental conditions, it becomes difficult to distinguish performance differences associated with the algorithmic configuration from those introduced by the experimental design itself.

To address this methodological gap, this work proposes a controlled and reproducible DE-based benchmarking framework for nonlinear multi-region controller-tuning problems. The framework establishes a structured experimental protocol based on fixed-seed initialization, predefined search bounds, consistent performance-index definitions, and unified simulation procedures. Within each operating region, all DE variants are evaluated using the same optimization configuration, initial population, objective function, controller structure, and simulation scenario, while region-specific settings are consistently maintained across variants. A total of twenty classical DE variants are evaluated under this protocol, providing a common experimental basis for reproducible and traceable performance comparison.

The original contribution of this work lies in the unified application of this controlled and reproducible benchmarking framework to the evaluation of classical DE variants for PID and PIDA controller tuning in a nonlinear wind turbine system. By applying a common experimental protocol within each operating region, the framework reduces the influence of initialization and implementation-related inconsistencies, providing a controlled basis for examining the convergence behavior and resulting controller performance associated with different mutation and crossover strategies. In addition, the study evaluates the suitability of PID and PIDA controller structures across partial- and above-rated operating regions and subsequently examines their generalization performance in simulation using a measured SCADA wind-speed profile. This additional evaluation makes it possible to identify cases in which favorable tuning-stage performance does not translate into satisfactory behavior under different wind conditions, revealing signs of overfitting. Therefore, the reported results provide insight into the combined influence of the DE configuration, controller structure, and wind turbine operating region under the considered experimental scenarios.

The remainder of this article is organized as follows. Section 2 describes the wind turbine model and its operating regions. Section 3 presents the controller structure and the optimization problem formulation. Section 4 introduces the proposed benchmarking framework. Section 5 details the simulation setup and experimental protocol. Section 6 discusses the comparative results. Finally, Section 7 summarizes the main conclusions and outlines directions for future research.

2. Wind Turbine Modeling and Operating Regions

The wind turbine model considered in this work is the National Renewable Energy Laboratory (NREL) 5-MW reference turbine [36], which is widely adopted as a benchmark in large scale wind turbine control studies. The aerodynamic, mechanical, pitch actuator, and generator parameters follow the original NREL specification, with a slight deviation in the rated wind speed due to differences in the aerodynamic representation, yielding a comparable maximum power coefficient and a rated wind speed of $11.6\mathrm{m}/\mathrm{s}$. Figure 1 provides a schematic representation of the overall wind turbine system, depicting the interconnection of its subsystems and the associated input–output signals. The physical interpretation of these variables is reported in Table 1.

The controlled variable is the generator angular speed ${\omega }_g$, while the control inputs correspond to the pitch-angle reference ${\beta }^{*}$ and the electromagnetic torque reference ${\tau }_{em}^{*}$. The electrical power output $P_e$ depends on both the generator speed and electromagnetic torque, and its regulation is therefore determined by the operating-region-specific control strategy. Each subsystem of the wind turbine model is characterized by a specific set of parameters that govern its dynamic behavior and interactions within the overall system. In the following, all subsystems are systematically described. The numerical values of the model parameters, adopted from the NREL reference turbine, are summarized in Table 2, which provides a consolidated overview of the model parameterization.

2.1. Aerodynamic Model

The aerodynamic subsystem models the conversion of wind kinetic energy into mechanical power and torque. The aerodynamic power is computed using Equation (1).

$$P_t={\displaystyle \frac{1}{2}}\rho \pi r^2{v}^3{C}_p(\lambda ,\beta ),$$

(1)

where $\lambda ={\omega }_{t}r/v$ denotes the tip-speed ratio. The power coefficient $C_p(\lambda ,\beta )$ is computed using Equations (2) and (3) following the formulation proposed in [37]. The corresponding aerodynamic torque is then obtained as ${\tau }_t=P_t/{\omega }_t$.

$$C_p(\lambda ,\beta )=0.5176\left({\displaystyle \frac{116}{{\lambda }_i}}-0.4\beta -5\right)e^{-21/{\lambda }_i}+0.0068\lambda$$

(2)

$${\displaystyle \frac{1}{{\lambda }_i}}={\displaystyle \frac{1}{\lambda +0.08\beta }}-{\displaystyle \frac{0.035}{{\beta }^3+1}}$$

(3)

The auxiliary variable ${\lambda }_i$ represents an effective tip-speed ratio introduced for the empirical approximation of $C_p(\lambda ,\beta )$.

2.2. Pitch Actuator Model

The pitch actuator adjusts the blade pitch angle $\beta$ and is modeled as a simplified second-order system, given by:

$$\ddot{\beta }=-{\omega }_n^2\beta -2\zeta {\omega }_n\dot{\beta }+{\omega }_n^2{\beta }^{*},$$

(4)

where ${\beta }^{*}$ is the controller-generated pitch reference, ${\omega }_n$ the natural frequency, and $\zeta$ the damping ratio. Physical constraints are enforced according to Table 2.

2.3. Drivetrain Model

The mechanical drivetrain is modeled using a two-mass representation that captures the torsional interaction between the turbine and the generator through a flexible shaft and an ideal gearbox. The turbine and generator inertias, $J_t$ and $J_g$, rotate at angular speeds ${\omega }_t$ and ${\omega }_g$, respectively, and are coupled by a shaft with stiffness $K_s$ and damping $D_s$. The applied torques include the aerodynamic torque ${\tau }_t$ acting on the turbine, the shaft torque ${\tau }_s$ transmitted through the drivetrain, and the electromagnetic torque ${\tau }_{em}$ generated by the electrical machine. The drivetrain dynamics are described by Equations (5)–(7).

$$J_t{\displaystyle \frac{d}{dt}}\left({\omega }_t\right)={\tau }_t-{\tau }_s-D_t{\omega }_t$$

(5)

$$J_g{\displaystyle \frac{d}{dt}}\left({\omega }_g\right)={\displaystyle \frac{{\tau }_s}{N_g}}-{\tau }_{em}-D_g{\omega }_g$$

(6)

$${\tau }_s=K_s\left({\theta }_t-{\displaystyle \frac{{\theta }_g}{N_g}}\right)+D_s\left({\omega }_t-{\displaystyle \frac{{\omega }_g}{N_g}}\right)$$

(7)

2.4. Electrical Generator Model

The electrical subsystem is modeled using a simplified first-order representation of the electromagnetic torque dynamics, as described by Equation (8).

$${\dot{\tau }}_{em}=-{\displaystyle \frac{1}{{\tau }_g}}{\tau }_{em}+{\displaystyle \frac{1}{{\tau }_g}}{\tau }_{em}^{*},$$

(8)

The parameter ${\tau }_{em}^{*}$ denotes the torque reference, whereas ${\tau }_g$ represents the electrical time constant. The electrical power output is computed as follows:

$$P_e=\eta {\tau }_{em}{\omega }_g,$$

(9)

where $\eta$ denotes the overall electromechanical efficiency. Torque and electrical-power saturation effects are neglected in the present model, allowing transient values beyond nominal rated conditions during simulation.

2.5. Wind Turbine Operating Regions

Wind turbine operation is typically divided into four regions, as illustrated in Figure 2.

In Region 2, the control objective is to support Maximum Power Point Tracking (MPPT) operation by regulating the generator speed through the electromagnetic torque while maintaining the pitch angle at its optimal value. This strategy aims to keep the turbine operating near the aerodynamic condition associated with maximum power extraction under partial-load conditions. In Region 3, corresponding to operation above the rated wind speed, the objective is to maintain the generator speed and electrical power near their rated values by adjusting the blade pitch while holding the electromagnetic torque constant.

3. Controller Structure and Problem Formulation

3.1. PIDA Controller Structure and Formulation

The transfer function of the PIDA controller is given by Equation (10), which incorporates low-pass filters in the derivative and acceleration actions to limit high-frequency amplification and improve the practical implementation of these terms.

$$C\left(s\right)=K_P+K_I{\displaystyle \frac{1}{s}}+K_D{\displaystyle \frac{s}{T_{D}s+1}}+K_A{\displaystyle \frac{s^2}{{(T_{A}s+1)}^2}}$$

(10)

The parameters $K_P$, $K_I$, $K_D$, and $K_A$ denote the proportional, integral, derivative, and acceleration gains, respectively, while $T_D$ and $T_A$ are the time constants of the low-pass filters associated with the derivative and acceleration actions, ensuring a proper transfer function. To mitigate integrator windup in the presence of actuator saturation, the PIDA controller incorporates an anti-windup compensation scheme based on the back-calculation method. The resulting controller structure is shown in Figure 3.

The gain $K_B$ is defined in Equation (11) and is obtained as the geometric mean of the time constants associated with the integral, derivative, and acceleration actions [38]. This gain weights the difference between the computed control signal $U_C\left(s\right)$ and the saturated control signal $U_S\left(s\right)$, facilitating integrator discharge and mitigating windup effects when actuator saturation occurs.

$$K_B={\displaystyle \frac{1}{\sqrt[3]{\frac{K_D{K}_A}{K_I{K}_P}}}}$$

(11)

3.2. Optimization Problem Formulation for Controller Tuning

The DE algorithm is employed to search for the controller parameter vector that minimizes the selected objective function within the prescribed search bounds:

$$\underset{\mathsf{\theta }}{min}J\left(\mathsf{\theta }\right)$$

(12)

where $\mathsf{\theta }=[K_P,K_I,K_D,K_A,T_D,T_A]$ denotes the decision-variable vector for the PIDA controller. For the PID structure, the corresponding parameter vector is obtained by excluding the acceleration-related parameters $K_A$ and $T_A$. In both cases, the optimization problem is formulated as a bounded continuous search, with each controller parameter constrained within predefined limits.

The objective function is formulated in accordance with the specific control objectives and performance requirements of the regulation problem. The proposed framework supports integral performance indices, allowing the cost formulation to be aligned with the targeted transient and steady-state characteristics.

In Region 2, the Integral of Absolute Error (IAE), defined in Equation (13), is adopted. This criterion penalizes the accumulated tracking error and provides an appropriate performance measure for the relatively smooth and quasi-linear dynamics typically observed in this operating regime.

$$IAE\left(\mathsf{\theta }\right)={\int }_0^t\left|e\left(t\right)\right|dt$$

(13)

In Region 3, where the primary control objective shifts toward disturbance rejection under highly nonlinear aerodynamic conditions, the Integral of Time-weighted Absolute Error (ITAE), defined in Equation (14), is employed. By introducing a time-weighting factor, this index penalizes persistent errors more aggressively, thereby promoting faster attenuation of disturbances and improved steady-state regulation

$$ITAE\left(\mathsf{\theta }\right)={\int }_0^{t}t\left|e\left(t\right)\right|dt$$

(14)

4. Differential Evolution–Based Tuning Framework

4.1. Differential Evolution Fundamentals

Differential Evolution is a stochastic, population-based optimization algorithm based on mutation and crossover operators. Algorithm 1 presents the canonical workflow of the DE algorithm. It begins by defining the control parameters, including the population size $N_p$, the mutation scaling factor F, and the crossover probability $C_r$. An initial population of candidate solutions is then randomly generated within the prescribed search bounds. At each generation G, for each target vector ${\mathit{x}}_{i,G}$, a mutant vector ${\mathit{v}}_{i,G}$ is generated according to the selected mutation strategy. Subsequently, each mutant vector is crossed with its corresponding target vector to produce a trial vector ${\mathit{u}}_{i,G}$. A one-to-one selection scheme is then applied, whereby the trial vector replaces the target only if it yields a lower cost value. This iterative process continues until the maximum number of generations $N_g$ is reached or another predefined stopping criterion is satisfied.

Algorithm 1 Standard DE algorithm

1: Define control parameters: $N_p$, F, $C_r$, and stopping criterion 2: Initialize population ${\left\{x_{i,0}\right\}}_{i=1}^{N_p}$ randomly within the search bounds 3: Set generation counter $G\leftarrow 0$ 4: while stopping criterion not satisfied do 5:     for $i=1$ to $N_p$ do 6:         Generate ${\mathit{v}}_{i,G}$ according to the selected mutation strategy 7:         Generate trial vector ${\mathit{u}}_{i,G}$ by crossing ${\mathit{v}}_{i,G}$ with ${\mathit{x}}_{i,G}$ 8:         if $f\left({\mathit{u}}_{i,G}\right)\le f\left({\mathit{x}}_{i,G}\right)$ then 9:            ${\mathit{x}}_{i,G+1}\leftarrow {\mathit{u}}_{i,G}$ 10:         else 11:            ${\mathit{x}}_{i,G+1}\leftarrow {\mathit{x}}_{i,G}$ 12:         end if 13:     end for 14:     $G\leftarrow G+1$ 15: end while 16: return best solution found in the final population

Different variants of the DE algorithm arise from the specific design of its mutation and crossover operators, which directly shape the search dynamics and convergence behavior. Although the overall optimization workflow remains unchanged, alternative formulations of these operators lead to distinct algorithmic behaviors. Accordingly, DE variants are commonly denoted using the standard DE/x/y/z notation, where x identifies the mutation strategy, y denotes the number of difference vectors used in the mutation process, and z specifies the crossover scheme. The mutation operator typically employs one or two difference vectors ($y=1$ or $y=2$), while the crossover mechanism is commonly implemented using either a binomial (bin) or an exponential (exp) scheme.

In this study, a total of 20 classical DE variants are considered, obtained by combining five mutation strategies (best, rand, current-to-best, current-to-rand, and rand-to-best), two differentiation orders ($y=1,2$), and two crossover operators (bin and exp). The complete set of analyzed variants is summarized in Table 3, and the same variants are consistently applied in both Region 2 and Region 3 to ensure a fair and homogeneous benchmarking framework.

4.2. Controlled Initialization and Fair Comparison Strategy

Due to the stochastic nature of Differential Evolution (DE), optimization outcomes depend not only on the algorithmic configuration of each variant but also on the random processes involved during the search. Randomness affects several stages of DE, including initial-population generation, individual selection for mutation, and the probabilistic decisions associated with crossover. Consequently, when different DE variants are evaluated under independently generated stochastic conditions, part of the observed performance differences may be influenced by random sampling effects in addition to the mutation strategy, differentiation order, and crossover scheme.

To promote a controlled and reproducible comparison, this study adopts a fixed-seed random-number generation strategy. For each operating region and controller architecture, the same random seed and random-number generator configuration are used for all evaluated DE variants. This procedure ensures that each optimization experiment can be reproduced from a predefined stochastic initialization and that all variants are evaluated from a common initial population under the corresponding experimental scenario. By controlling the initialization procedure, the framework reduces one source of experimental variability and provides a more consistent basis for comparing the behavior of the different DE configurations.

Reproducibility is implemented in MATLAB through the rng function, using a fixed positive integer seed and a predefined generator type. Before executing each set of comparable optimization runs, the random number generator is initialized explicitly, ensuring that the initial population and subsequent stochastic operations can be reproduced exactly. In contrast, if the seed is not specified, MATLAB initializes the generator using time-dependent information from the system clock, which leads to non-reproducible results.

This controlled initialization strategy is particularly relevant because the study compares multiple DE variants applied to the same controller-tuning problem. Within each operating region, the wind profile, simulation horizon, controller structure, parameter bounds, objective function, population size, number of generations, stopping criterion, and initial population are kept unchanged across variants. Region-specific settings are therefore consistently applied to all DE configurations, providing a common experimental basis for the comparison. Under these conditions, differences in convergence behavior and resulting controller performance can be examined while reducing the influence of initialization and implementation-related inconsistencies.

The fixed-seed protocol is intended to ensure reproducibility and controlled comparison under the defined experimental conditions. Accordingly, the reported results characterize the behavior of the evaluated DE variants within the considered stochastic realization and operating scenarios rather than providing a statistical characterization across independent random realizations. Nevertheless, the protocol enables the reported controller parameters, convergence curves, and performance indices to be regenerated and provides a transparent basis for extending the benchmark to additional DE variants, controller structures, or repeated stochastic evaluations in future studies.

5. Simulation Setup and Experimental Design

This section outlines the simulation environment and experimental design adopted to ensure a consistent and reproducible application of the DE–based tuning framework in Regions 2 and 3 of the wind turbine, addressing the region-specific control problems defined in Section 2.5. All simulations were implemented in the Simulink environment of MATLAB R2024b, where the system model, control architecture, and optimization procedures were fully integrated. The simulation setup, operating conditions, performance metrics, and optimization parameters are described in the following subsections. To promote transparency and reproducibility, the complete implementation of the proposed framework, including the Simulink models and optimization scripts, is publicly available through an external repository [39].

5.1. Simulation Environment

The wind turbine model was implemented according to the mathematical formulation presented in Section 2. All simulations were performed using a variable-step solver with MATLAB/Simulink automatic solver selection. The relative tolerance was set to $1\times {10}^{-3}$, whereas the absolute tolerance, maximum step size, minimum step size, and initial step size were automatically determined by Simulink. Zero-crossing detection was configured using the local block settings with the nonadaptive algorithm. The same numerical solver configuration was maintained across all DE variants within each operating region to ensure consistent numerical integration conditions throughout the benchmarking procedure. The complete Simulink implementation of the wind turbine system, including all subsystems described in Section 2, is shown in Figure 4, which illustrates the block-level functional architecture of the nonlinear model and the corresponding signal flow among subsystems. The diagram represents the closed-loop operation of the system under normal operating conditions, where wind speed and reference signals are processed by the control strategy to generate pitch-angle and electromagnetic-torque references according to the current operating region.

To evaluate the proposed tuning framework, two independent simulation experiments were designed, each corresponding to a distinct operating region of the wind turbine.

In Region 2, the process variable is the generator angular speed, whose reference is a step sequence with values $[94.11,122.9,37.64]\mathrm{rad}/\mathrm{s}$, enforced through electromagnetic torque regulation. These values are selected to excite system nonlinearities by combining sub-nominal, nominal, and low-speed operating points. The corresponding wind-speed profile is $[7.5,11.4,3]\mathrm{m}/\mathrm{s}$, with transitions every $25\mathrm{s}$ and a total simulation time of $75\mathrm{s}$.

In Region 3, the process variable is held constant at the rated generator speed by regulating the blade pitch angle, with initial conditions ${\omega }_{t_i}=1.267\mathrm{rad}/\mathrm{s}$ and ${\omega }_{g_i}=122.9\mathrm{rad}/\mathrm{s}$. The wind-speed profile is defined as $[12,16,24.9,18]\mathrm{m}/\mathrm{s}$, introducing disturbances of varying amplitudes across the full operating range, with transitions every $125\mathrm{s}$ and a total simulation time of $500\mathrm{s}$.

All simulations were performed using the control structure shown in Figure 3, implemented in discrete time via the Backward–Euler formulation with a controller sampling time of $0.05\mathrm{s}$. This approach enables a unified control implementation across operating regions using a single adaptive algorithm, avoiding switches between continuous-time controllers. Moreover, it prevents the numerical slowdown typically associated with continuous-time tuning when the filter time constants $T_D$ and $T_A$ approach zero. In such cases, derivatives are approximated numerically in the discrete domain, thereby ensuring stable and computationally efficient simulations.

5.2. DE Parameter Settings and Reproducibility

This subsection summarizes the configuration of the DE algorithm used for controller tuning in both operating regions. To ensure reproducibility and fair benchmarking, a deterministic setup was adopted, with fixed random-number generation and identical initial populations across all DE variants. The main DE parameters, described in Section 4.1, are listed in Table 4.

The observed differences between operating regions stem from their distinct dynamic characteristics. In Region 2, the system exhibits relatively simple and less sensitive dynamics, resulting in a smoother objective-function landscape. Under these conditions, exploration–primarily governed by the DE mutation factor F and crossover rate $C_r$–is more effective and allows convergence to be achieved within a smaller number of generations.

In contrast, Region 3 is characterized by stronger nonlinearities and higher sensitivity to control actions, mainly due to the interaction between pitch dynamics and aerodynamic effects. Such characteristics typically lead to a more rugged and sensitive search space, in which excessive exploration may degrade convergence or induce unstable candidate solutions. Consequently, the admissible level of exploration must be reduced, resulting in a more constrained and exploitation-oriented search process.

Despite these differences, identical random-number generation settings are employed in all cases, using a fixed seed of 1 and the twister generator, in order to ensure reproducibility. The only stopping criterion considered is the completion of the predefined number of generations, and the cost function corresponds to that defined in SubSection 3.2. The search space for Region 2 is defined as

$$\{K_P,K_I,K_D,K_A\in [-4\times {10}^4,0];T_D,T_A\in [0,1]\},$$

whereas for Region 3 it is given by

$$\{K_P,K_I,K_D,K_A\in [-5,0];T_D,T_A\in [0,3]\}.$$

The search bounds were selected based on the order of magnitude of the control signals associated with each operating region.

In terms of computational cost, the optimization process is mainly determined by the number of iterations required to tune each DE variant. In this context, a simulation refers to the complete optimization process of one DE variant, whereas an iteration refers to the evaluation of a single candidate individual through a complete closed-loop Simulink execution of the nonlinear wind turbine model.

The number of iterations per simulation was determined from the population size and the maximum number of generations defined in Table 4. Based on these parameters, each DE variant requires approximately 5000 iterations in Region 2 and 12000 iterations in Region 3. Each iteration in Region 2 evaluates 75 s of turbine operation, while each iteration in Region 3 evaluates 500 s. Therefore, the computational demand is considerably higher in Region 3 because each simulation involves both a larger number of iterations and a longer simulation horizon for each evaluated individual.

Using a single processor core, the execution time of one complete simulation ranged from approximately 40 to 60 min in Region 2, depending on the evaluated DE variant. In Region 3, the execution time increased to approximately 3.5 to 4 h per variant. The developed framework also includes the option of enabling MATLAB’s Parallel Computing Toolbox. When this option is activated, the execution time can be reduced by up to four times, since several iteration loops can be evaluated simultaneously using a larger number of processor cores.

All simulations were executed in MATLAB/Simulink R2024b on a workstation equipped with an Intel(R) Core(TM) i7-10750H CPU @ 2.60 GHz, 6 cores, 32 GB RAM, an NVIDIA GeForce RTX 2060 GPU with 6 GB of dedicated memory, a 64-bit x64-based operating system, and Windows 11 version 24H2.

6. Results and Discussion

6.1. Region 2 Results

6.1.1. Optimization Performance of DE Variants in Region 2

The optimization performance of the DE variants in Region 2 is analyzed through their convergence speed. Figure 5a–e illustrate the evolution of the best objective-function value ($J_{\mathrm{Best}}$) over successive generations for all variants, grouped by mutation strategy.

As shown in Figure 5a,b,e, exploitation-oriented variants based on best, current-to-best, and rand-to-best mutation strategies with a single differentiation exhibit faster convergence, with $J_{\mathrm{Best}}$ stabilizing within a smaller number of generations. This behavior is mainly due to the stronger exploitation tendency characteristic of these variants. Here, new candidate solutions are produced around high-performing individuals, which confines the search to a more limited region of the solution space.

In contrast, exploration-oriented variants, particularly those based on rand-type strategies and/or two differentiations, tend to converge more slowly, as shown in Figure 5c,d. This slower convergence is associated with increased population diversity and a wider search region. Nevertheless, after a sufficient number of generations, the performance differences between exploration- and exploitation-based strategies become less pronounced, and further improvements in $J_{\mathrm{Best}}$ appear increasingly influenced by stochastic effects rather than by the mutation strategy alone.

The optimization results are summarized in Table 5, which presents the optimized controller parameters and the corresponding objective-function values obtained for each DE variant under the Region 2 configuration described in Section 5.

According to Table 5, the DE/best/1/bin variant attains the lowest objective value among the evaluated strategies, reaching 856.8178 at generation 94. However, this result does not imply slow convergence. By generation 25, the objective value has already decreased to 856.9840, indicating that the variant approaches a near-optimal solution within relatively few generations, followed by only marginal refinements. The final improvement amounts to only 0.1662 units (approximately 0.02%), which is consistent with the convergence behavior observed in Figure 5a.

Although DE/best/1/bin achieves the lowest objective value, a comparable value of 857.7840 at generation 94 is also reached by the DE/rand/1/bin variant. The main difference lies in the convergence dynamics, as this variant exhibits a slower rate of improvement, as shown in Figure 5d. For instance, at generation 93, the objective value is 860.267, whereas at generation 25 it remains at 900.2380, representing a reduction of 42.4540 units (4.44%). The remaining variants exhibit convergence patterns consistent with their respective mutation strategies, as previously discussed.

6.1.2. Closed-Loop Controller Performance in Region 2

The closed-loop dynamic response of the system and the corresponding control actions are shown in Figure 6a,b. The controllers considered in this analysis are the PIDA controller tuned using the DE variant that achieved the lowest objective-function value under the Region 2 fixed-seed configuration reported in Table 5, a PID controller tuned using the same variant, and a PID controller tuned using the MATLAB PID Tuner, which serves as the baseline for performance comparison. The corresponding controller parameters used in this comparison are summarized in Table 6.

Controller performance is quantitatively assessed using transient-response characteristics together with the integral performance index (IAE) associated with the cost function in this operative region. Owing to the nonlinear dynamics of the system and the presence of multiple reference changes, transient metrics vary across transitions. Therefore, overshoot and settling time are reported for each setpoint change defined in Section 5.1 and illustrated in Figure 6a. These metrics, along with the IAE index, are summarized in Table 7.

As indicated in Table 7, several PIDA controllers achieve slightly lower IAE values than the PID controller tuned with DE/best/1/bin, while all evaluated DE-tuned controllers yield lower IAE values than the baseline controller under the considered Region 2 scenario. In particular, the PIDA controller tuned with DE/best/1/bin achieves a 0.2 % improvement over the PID tuned with the same variant according to the IAE criterion, which is used as the cost function in Region 2, and a 2.93 % improvement relative to the baseline controller.

Nevertheless, these controllers exhibit noticeable differences in their transient response characteristics. In general, the baseline controller presents no overshoot, indicating a more damped response; however, this also results in slower dynamics. This behavior is reflected in the settling times, where the baseline controller consistently shows larger values for all setpoint changes. In contrast, PIDA-DE/best/1/bin exhibits a slight advantage over PID-DE/best/1/bin in terms of overshoot, while maintaining very similar settling times.

The limited performance gap between the PID and PIDA controllers in Region 2 can be attributed to the low-order, quasi-linear behavior dominating this operating regime. As discussed in Section 2.4, the electrical generator is modeled as a first-order linear system and acts as the dominant actuator in Region 2. Under these conditions, the additional acceleration action of the PIDA controller does not yield systematic performance improvements, and a well-tuned PID controller can achieve comparable closed-loop behavior.

This interpretation is consistent with the control signals shown in Figure 6b, where the DE-tuned PID exhibits a more aggressive control action, accelerating convergence toward the reference while maintaining a transient response comparable to that of the PIDA controller.

Finally, under the considered Region 2 scenario, the baseline controller yields a higher IAE value than all evaluated DE-tuned controllers. This difference is consistent with the fact that the DE-based approach directly optimizes controller parameters on the nonlinear model under the Region 2 operating conditions, whereas the baseline controller is obtained using MATLAB PID Tuner. Overall, these results support the effectiveness of the proposed DE-based tuning framework for the evaluated Region 2 configuration.

6.2. Region 3 Results

6.2.1. Optimization Performance of DE Variants in Region 3

Region 3 is analyzed following the same methodology adopted for Region 2, with particular emphasis on convergence speed under a more demanding and nonlinear operating regime. Figure 7a–e illustrate the evolution of the best objective value ($J_{\mathrm{Best}}$) over successive generations for all DE variants, grouped according to their mutation strategy.

Compared to Region 2, the convergence behavior observed in Region 3 exhibits more pronounced differences among the evaluated DE variants. As shown in Figure 7a,b,e, variants based on best-individual mutation strategies combined with a single differentiation order tend to converge more rapidly, reaching stable objective values within a smaller number of generations.

This behavior reflects a stronger emphasis on exploitation, whereby the search process remains concentrated around the best-performing individuals, limiting large parameter variations and promoting faster convergence under the considered operating conditions.

In contrast, the convergence curves shown in Figure 7c,d indicate that exploration-oriented variants generally exhibit slower convergence. In Region 3, this effect is more pronounced than in Region 2, which can be attributed to the increased nonlinearity and sensitivity of the optimization landscape. The broader exploration induced by these variants leads to larger parameter perturbations, which in some cases result in poorer objective values and less consistent convergence behavior, thereby delaying convergence.

The optimization results are summarized in Table 8, which reports the optimized parameter vectors and corresponding objective values obtained for each DE variant under the parameters defined in Section 5. All variants were tuned using the PIDA controller structure.

Under the considered Region 3 fixed-seed configuration, the DE/best/1/bin variant achieves the lowest recorded objective-function value, reaching 37840.68 at generation 184. Similar to Region 2, this variant exhibits only minor improvements in the cost function during later generations. In fact, at generation 91, the objective value is 37849.00, indicating that before reaching half of the total generations, the solution is already close to its final value, with a difference of 8.32 units (approximately 0.02%).

In contrast to Region 2, two additional variants achieve objective-function values very close to the lowest recorded value. The first is DE/current-to-best/1/bin, reaching 37840.70 at generation 188. Similar to DE/best/1/bin, this variant exhibits only minor improvements in later generations, reaching 37847.10 at generation 103, which corresponds to an improvement of approximately 0.02%. The second variant is DE/rand-to-best/1/bin, reaching an objective value of 37840.96 at generation 191. As in the previous cases, the improvement after generation 119 is approximately 0.02%.

In this region, the three variants with the lowest recorded objective-function values also exhibit comparable convergence speeds. These variants are based on mutation strategies that use the best individual. For the considered Region 3 scenario, the results suggest that an appropriate balance between exploration and exploitation is important for convergence behavior. The remaining variants exhibit convergence patterns consistent with their respective mutation strategies, as previously discussed.

6.2.2. Closed-Loop Controller Performance in Region 3

Figure 8a,b show the system responses and corresponding control signals under Region 3 operating conditions. The results are presented for the PID and PIDA controllers tuned using DE/best/1/bin, which achieved the lowest objective-function value under the Region 3 fixed-seed configuration reported in Table 8, together with a PID controller tuned using the MATLAB PID Tuner, which serves as the baseline. The corresponding controller parameters used in this comparison are summarized in Table 9.

The system response in Region 3 is quantified using integral performance indices. Unlike Region 2, classical transient-response metrics are not considered, as the primary control objective is disturbance rejection rather than reference tracking. In this operating region, metrics such as overshoot and settling time are largely influenced by the disturbance magnitude rather than by the intrinsic system dynamics. Accordingly, only integral performance indices are reported in Table 9.

As shown in Table 9, the PIDA controller outperforms the PID controller in Region 3, in contrast to the behavior observed in Region 2. Using DE/best/1/bin, which yielded the lowest objective-function value under the considered Region 3 configuration, the PIDA controller achieves an improvement of approximately 64 % in the ITAE criterion–used as the cost function in Region 3–relative to the PID controller tuned with the same DE variant, and about 92 % compared with the baseline controller.

This marked performance gain is attributed to the increased dynamic complexity of Region 3, where the blade pitch angle becomes the dominant control input. According to Equations (1)–(3), the pitch angle directly influences the dominant nonlinear aerodynamic component of the system. Under these conditions, the system exhibits higher-order and strongly nonlinear dynamics, which justify the use of a more complex control structure. The additional zeros introduced by the acceleration action of the PIDA controller enable more effective compensation of nonlinearities, as reflected in Figure 8b, where the PIDA controller exhibits faster and less oscillatory control actions than the PID.

As in Region 2, the evaluated DE-tuned controllers yield lower objective-function values than the baseline under the considered Region 3 scenario; however, the performance gap is more pronounced in this region. This can be attributed to the fact that the DE-based framework directly optimizes the controller gains on the nonlinear model under the defined Region 3 operating conditions, which is reflected in lower ITAE values under the evaluated disturbance-rejection scenario.

These results support the suitability of the proposed DE-based tuning framework for the considered Region 3 configuration.

6.3. Validation under Real Wind Profile

To validate the controllers obtained in Section 6.1.2 and Section 6.2.2, a real wind-speed profile measured by a Supervisory Control and Data Acquisition (SCADA) system from a wind turbine was used [40]. Specifically, the wind profile recorded on January 7, 2018, shown in Figure 9, was selected.

Both the baseline controller and the selected DE-tuned controllers were evaluated using the wind profile shown in Figure 9, which spans wind-speed ranges suitable for analyzing transitions between operating regions.

The generator angular speed and its reference, obtained with the baseline controllers, are shown in Figure 10a.

Although the controlled variable is the generator speed, the most relevant output for power production is the active electrical power, defined by Equation (9) and shown in Figure 10b.

The control signals corresponding to Figure 10a,b are shown in Figure 10c,d. As observed in Figure 10c, the electromagnetic torque presents a smoother behavior than the pitch angle shown in Figure 10d. In contrast, the pitch signal exhibits more pronounced oscillations, which are transferred to the generator-speed response and contribute to the oscillatory behavior observed in Figure 10a.

As discussed in Section 6.2, the baseline PID controller in Region 3 exhibits a relatively slow dynamic response. Consequently, the control action tends to oscillate in an attempt to compensate for wind-induced speed deviations and recover the rated operating condition. This behavior not only increases the oscillatory content of the system output, but also has a physical implication for the turbine. Repeated pitch variations may increase the activity and wear of the pitch actuator, while oscillations in generator speed may introduce additional mechanical stress in the drivetrain and rotor shaft. Therefore, the analysis of the control signals confirms that the baseline PID controller does not only present poorer numerical performance, but also a less desirable dynamic behavior from an operational and mechanical perspective.

The same evaluation scenario was subsequently applied using the PIDA-DE/best/1/bin controller in Region 2 and the PIDA-DE/best/1/bin controller in Region 3. The resulting generator angular speed is shown in Figure 11a. It is worth noting that the reference signals computed by the MPPT strategy exhibit slight variations depending on the selected controller, as they are influenced by internal model variables.

Figure 11a reveals a noteworthy behavior. Under high-load operating conditions, the controller tuned for Region 3 outperforms the baseline controller. This improvement is also reflected in the pitch-angle signal shown in Figure 11d, which is clearly smoother than the corresponding signal obtained with the baseline PID controller.

Although the DE-tuned PIDA controller shows a clear improvement in Region 3, its performance is comparatively less favorable in Region 2. In this operating region, the generator-speed response exhibits small oscillations, particularly around $0.5·{10}^4\mathrm{s}$ and during the interval between $3·{10}^4\mathrm{s}$ and $4·{10}^4\mathrm{s}$.

This response suggests a highly oscillatory control behavior in Region 2, with a consequent degradation in the quality of the control action. The corresponding electromagnetic torque control effort is shown in Figure 11c. As observed, the control signal exhibits repeated oscillatory corrections under this operating condition, which are then reflected in the generator-speed response shown in Figure 11a. From a physical perspective, this behavior is undesirable because persistent oscillatory control actions may increase dynamic loading on the drivetrain and rotor shaft.

In accordance with Equation (9), the active electrical power also exhibits oscillatory behavior in Region 2, as shown in Figure 11b. This response is not favorable for power generation, since fluctuations in generator speed and electromagnetic torque are directly reflected in the generated electrical power. From an operational perspective, these power variations may lead to undesirable fluctuations in the power injected into the grid, potentially affecting power quality if not properly mitigated. Therefore, although the DE-tuned PIDA controller improves the response in Region 3, its behavior in Region 2 indicates poorer generalization under real wind conditions.

The observed performance degradation in Region 2 can be further interpreted by analyzing Table 5. Two controllers are particularly noteworthy: DE/best/1/bin, which achieved the lowest J value under the considered Region 2 configuration, and the closely following DE/rand/1/bin. These controllers are of special interest because their gain values differ significantly from those of the remaining controllers, which tend to cluster around a similar local minimum of the objective function J.

This behavior may indicate a degree of overfitting to the optimization scenario. This explains why these controllers achieve strong performance under the tuning conditions. In particular, the best performance is obtained with DE/best/1/bin. However, they exhibit degraded performance when evaluated under different operating scenarios.

Given this scenario, it becomes necessary to select a controller that avoids the previously observed signs of overfitting. According to Table 5, among the PIDA controllers that did not exhibit these signs, DE/rand-to-best/1/bin achieved the lowest objective-function value, with $J=863.06$. However, as shown in Table 6, the PID controller tuned using DE/best/1/bin achieves a lower value of $J=858.60$, outperforming the PIDA controllers that maintained satisfactory generalization behavior under the measured wind profile. Therefore, the PID controller is selected for operation in Region 2, while the DE-tuned PIDA controller selected for Region 3 is retained for above-rated operation. The generator-speed response and electrical-power output obtained with this combined control strategy are presented in Figure 12a and Figure 12b, respectively. The corresponding electromagnetic-torque and pitch-angle signals obtained with the selected DE-tuned controllers are shown in Figure 12c and Figure 12d, respectively. Compared with the previous cases, both control signals exhibit noticeably smoother behavior, which contributes to reducing the oscillations observed in the generator-speed and active-power responses.

With this combination of controllers, both Region 2 and Region 3 outperform the baseline controller, achieving improved reference tracking with reduced error. These results can be quantitatively assessed using integral performance indices, whose values are reported in Table 10.

With the appropriate selection of controllers, the performance evaluated using real data corresponding to a full day of operation shows a significant improvement for the DE-based controllers. Specifically, these controllers achieve reductions of 96.7 %, 95.8 %, 98.1 %, and 99.9 % in the IAE, ISE, ITAE, and ITSE indices, respectively.

These results show that the benefits of the optimal DE-tuned controllers are not limited to improved tracking performance. The smoother control signals also contribute to a more stable active-power response and may reduce unnecessary mechanical stress on the drivetrain and pitch actuators. This highlights the importance of appropriate controller tuning in wind turbine applications, demonstrating that performance does not depend solely on using more complex controller structures, but also on selecting gain parameters that produce a balanced dynamic response.

In contrast to the results obtained under tuning conditions, the most suitable option for the partial-load region under a real wind profile, with the MPPT strategy active, is a well-tuned PID controller, owing to the relative simplicity of the system dynamics in this operating regime. On the other hand, Region 3 presents more complex dynamics, as discussed in Section 6.2.2, where the PIDA architecture offers a clear advantage over the classical PID controller.

Accurate control of the generator angular speed results in a more stable active power output, since electrical power depends directly on generator speed. Accordingly, the active power dynamics closely follow those of the generator speed, as shown in Figure 10b, Figure 11b and Figure 12b, highlighting the role of precise speed regulation in ensuring stable power generation under real operating conditions.

It should be noted that negative electrical power values appear during the initial transient, as the turbine temporarily draws energy from the grid to overcome inertial effects and initiate drivetrain rotation.

7. Conclusions

The results obtained under the considered experimental scenarios indicate that the DE-tuned controllers outperform the MATLAB PID Tuner baseline controllers across the evaluated performance indices. This advantage is particularly pronounced in Region 3, where the stronger nonlinear aerodynamic dynamics and the disturbance-rejection requirements make controller tuning more demanding. Under these conditions, the additional degrees of freedom provided by the PIDA structure result in a substantial performance improvement over the conventional PID controller. In Region 2, however, the performance difference between PID and PIDA structures is considerably smaller, and the subsequent evaluation under a measured wind profile indicates that a well-tuned PID controller provides a more suitable response. These findings indicate that controller performance depends not only on the selected architecture, but also on its compatibility with the dynamic characteristics and control objectives of each wind turbine operating region.

The results also show that the controller achieving the lowest objective-function value under the tuning conditions is not necessarily the most suitable when evaluated under a different wind scenario. In Region 2, some controllers that achieved favorable tuning-stage results exhibited oscillatory behavior when subsequently evaluated in simulation using a measured SCADA wind-speed profile, indicating reduced generalization performance and signs of overfitting. This observation highlights the importance of evaluating tuned controllers beyond the specific scenario used during optimization. Consequently, controller selection should not rely exclusively on the objective-function value obtained during tuning, but should also consider the dynamic response under operating conditions that differ from those used in the optimization process.

A central methodological contribution of this work is the adoption of a controlled and reproducible benchmarking protocol for comparing classical DE variants. Although Differential Evolution is inherently stochastic, the use of a fixed random seed, a predefined random-number generator, common initial populations, and consistent simulation conditions within each operating region enables the reported optimization experiments to be reproduced under the same computational environment. This procedure reduces initialization- and implementation-related inconsistencies and provides a common experimental basis for examining differences in convergence behavior and resulting controller performance. Nevertheless, the reported results correspond to the considered fixed-seed stochastic realization and operating scenarios and should therefore be interpreted as a reproducible within-scenario comparison rather than as a statistical characterization of the DE variants across independent random realizations.

The proposed framework provides a structured basis that can be extended to other controller architectures, optimization variants, and nonlinear control problems following the same experimental methodology. Its current implementation, however, relies on an accurate mathematical representation of the controlled system because the complete tuning and validation process is performed in simulation. Furthermore, validation using measured wind-speed data does not reproduce all effects associated with physical wind turbine operation, such as measurement noise, communication delays, parameter uncertainty, actuator nonlinearities, or unmodeled dynamics. These aspects represent important limitations when considering direct implementation in real systems.

Future work should therefore extend the framework toward Hardware-in-the-Loop (HIL) and real-time simulation environments, where effects such as delays, measurement noise, parameter uncertainty, and nonlinear actuator behavior can be explicitly evaluated. The present fixed-seed benchmark can also be complemented with multiple independent random seeds and statistical analyses to quantify performance variability and evaluate the stability of the observed DE-variant rankings across stochastic realizations. Another relevant extension is the incorporation of multi-objective optimization, including criteria such as control effort, actuator activity, and electrical power performance, which would provide a broader assessment of controller behavior and further reduce the gap between simulation-based tuning and realistic wind turbine operation.