An Easily Understandable Grey Wolf Optimizer and Its Application to Fuzzy Controller Tuning

: This paper proposes an easily understandable Grey Wolf Optimizer (GWO) applied to the optimal tuning of the parameters of Takagi-Sugeno proportional-integral fuzzy controllers (T-S PI-FCs). GWO is employed for solving optimization problems focused on the minimization of discrete-time objective functions deﬁned as the weighted sum of the absolute value of the control error and of the squared output sensitivity function, and the vector variable consists of the tuning parameters of the T-S PI-FCs. Since the sensitivity functions are introduced with respect to the parametric variations of the process, solving these optimization problems is important as it leads to fuzzy control systems with a reduced process parametric sensitivity obtained by a GWO-based fuzzy controller tuning approach. GWO algorithms applied with this regard are formulated in easily understandable terms for both vector and scalar operations, and discussions on stability, convergence, and parameter settings are offered. The controlled processes referred to in the course of this paper belong to a family of nonlinear servo systems, which are modeled by second order dynamics plus a saturation and dead zone static nonlinearity. Experimental results concerning the angular position control of a laboratory servo system are included for validating the proposed method.


Introduction
Fuzzy control has proven efficiency in coping with various simple and complex processes.The combinations with nature-inspired optimization algorithms, also known as metaheuristics, have a positive impact on the performance of fuzzy control systems as they contribute to the systematic design and tuning.Some recent applications of metaheuristics to the parameter tuning of fuzzy controllers dedicated to servo systems will be outlined as follows.Representative overviews on such industrial applications are offered in [1][2][3].Proportional-integral-derivative (PID)-fuzzy controllers are tuned in [4,5] using Particle Swarm Optimization (PSO)-based approaches and are tested by experiments or simulations on an electrical direct current (DC) drive benchmark and a laboratory micro air vehicle controller, respectively.A hybrid PSO and pattern search optimized PI-fuzzy controller is applied in [6] to the automatic generation control of multi-area power systems and is validated by simulation.The Gravitational Search Algorithm (GSA) is employed in [7,8] in the optimal tuning of PI-and PID-fuzzy controllers for DC servo systems and load frequency control in power systems and is validated by simulation and experimental results.Charged System Search algorithms are applied in [9] to the optimal tuning of PI-fuzzy controllers for DC servo systems.Ant Colony Optimization is applied in [10,11] to the optimal tuning of fuzzy controllers for robots and ball and bam systems.
The Grey Wolf Optimizer (GWO) algorithm is proposed in [12] on the basis of modeling grey wolf social hierarchy and hunting habits towards finding prey, represented by the solution to the optimization problem.The social hierarchy is simulated by categorizing the population of search agents into four types of individuals, i.e., alpha, beta, delta, and omega, based on their fitness.The search process is modeled with the aim of mimicking the hunting behavior of grey wolfs making use of three stages, searching, encircling, and attacking the prey.The first two stages are dedicated to exploration and the last one covers the exploitation.The reduced number of search parameters is an important advantage of GWO algorithms reflected in various applications, which include blackout risk prevention in smart grids [13], training multi-layer perceptrons [14], optimization of reactive power dispatch [15], solutions to benchmarks generally used to test optimization algorithms [16], hyperspectral band selection [17], maximum power point tracking [18], and optimal tuning of PI-and PID-fuzzy controllers [19][20][21].
This work builds upon the authors' previous results on the optimal tuning of Takagi-Sugeno proportional-integral fuzzy controllers (T-S PI-FCs) and fuzzy models by means of various standard nature-inspired algorithms like: Gravitational Search Algorithm [7,22], Charged System Search [9], adaptive variations [23], or applications of multiple optimization algorithms [3,24].This paper proposes an easily understandable presentation of GWO.Based on this, GWO is later applied in an innovative method for the optimal tuning of T-S PI-FCs to ensure fuzzy control systems with a reduced process parametric sensitivity by the inclusion of sensitivity functions in the objective functions, which are minimized by the GWO in appropriately defined optimization problems.The presentation is focused on a family of servo systems considered as controlled processes, which are modeled by second-order dynamics with an integral component, variable parameters, along with a saturation and dead zone static nonlinearity.The results are confirmed through comparisons with other metaheuristics that have the same unsupervised characteristics.
This paper is organized as follows: the optimization problem is set in the next section.GWO is presented and analyzed in Section 3 along with the tuning approach for T-S PI-FCs dedicated to the considered family of servo systems.The case study that targets the GWO-based tuning of T-S PI-FCs for the angular position control of a nonlinear DC servo system is treated in Section 4. The conclusions are summarized in Section 5.

Problem Setting
The control system structure is presented in Figure 1, where: P-the process, which belongs to a family of nonlinear servo systems, FC-the fuzzy controller, F-the reference input filter, r-the reference input, r 1 -the reference input filtered through F, y-the controlled output, u-the control signal, and e-the control error.The load-type disturbance inputs are not applied in Figure 1 as the integral component of the controller copes with them.
power systems and is validated by simulation and experimental results.Charged System Search algorithms are applied in [9] to the optimal tuning of PI-fuzzy controllers for DC servo systems.Ant Colony Optimization is applied in [10,11] to the optimal tuning of fuzzy controllers for robots and ball and bam systems.
The Grey Wolf Optimizer (GWO) algorithm is proposed in [12] on the basis of modeling grey wolf social hierarchy and hunting habits towards finding prey, represented by the solution to the optimization problem.The social hierarchy is simulated by categorizing the population of search agents into four types of individuals, i.e., alpha, beta, delta, and omega, based on their fitness.The search process is modeled with the aim of mimicking the hunting behavior of grey wolfs making use of three stages, searching, encircling, and attacking the prey.The first two stages are dedicated to exploration and the last one covers the exploitation.The reduced number of search parameters is an important advantage of GWO algorithms reflected in various applications, which include blackout risk prevention in smart grids [13], training multi-layer perceptrons [14], optimization of reactive power dispatch [15], solutions to benchmarks generally used to test optimization algorithms [16], hyperspectral band selection [17], maximum power point tracking [18], and optimal tuning of PI-and PID-fuzzy controllers [19][20][21].
This work builds upon the authors' previous results on the optimal tuning of Takagi-Sugeno proportional-integral fuzzy controllers (T-S PI-FCs) and fuzzy models by means of various standard nature-inspired algorithms like: Gravitational Search Algorithm [7,22], Charged System Search [9], adaptive variations [23], or applications of multiple optimization algorithms [3,24].This paper proposes an easily understandable presentation of GWO.Based on this, GWO is later applied in an innovative method for the optimal tuning of T-S PI-FCs to ensure fuzzy control systems with a reduced process parametric sensitivity by the inclusion of sensitivity functions in the objective functions, which are minimized by the GWO in appropriately defined optimization problems.The presentation is focused on a family of servo systems considered as controlled processes, which are modeled by second-order dynamics with an integral component, variable parameters, along with a saturation and dead zone static nonlinearity.The results are confirmed through comparisons with other metaheuristics that have the same unsupervised characteristics.
This paper is organized as follows: the optimization problem is set in the next section.GWO is presented and analyzed in Section 3 along with the tuning approach for T-S PI-FCs dedicated to the considered family of servo systems.The case study that targets the GWO-based tuning of T-S PI-FCs for the angular position control of a nonlinear DC servo system is treated in Section 4. The conclusions are summarized in Section 5.

Problem Setting
The control system structure is presented in Figure 1, where: P-the process, which belongs to a family of nonlinear servo systems, FC-the fuzzy controller, F-the reference input filter, r -the reference input, 1 r -the reference input filtered through F, y -the controlled output, u -the control signal, and e -the control error.The load-type disturbance inputs are not applied in Figure 1 as the integral component of the controller copes with them.The family of servo systems is characterized by the state-space model where: t ≥ 0-the continuous time, k P > 0-the servo system gain, T Σ > 0-the small time constant, u(t)-a pulse width modulated signal, x 1 (t) = α(t)-the angular position, x 2 (t) = ω(t)-the angular speed, m(t)-the output of the saturation and dead zone static nonlinearity with the parameters that fulfill 0 < u a < u b , 0 < u c < u b , and T-matrix transposition.
The nonlinearity in Equation ( 1) is neglected in the transfer function of the simplified model of the process where k EP is the equivalent process gain, with the expression The transfer function P(s) is used in linear and fuzzy controller design, and PI controllers are recommended in [25,26].The PI controllers are used in the block FC in Figure 1 and their transfer function is where k C > 0 (or k c > 0) is the controller gain and T i > 0 is the integral time constant.
Since the simplified model P(s) is used in the controller design, and it is characterized by only two parameters, it is justified to consider that the process parameters k P and T Σ are variable and the other ones are constant.The process parameter vector is which leads to the definition of the state sensitivity functions λ α j υ , υ = 1...n, and the output sensitivity function σ: where the subscript 0 is inserted to outline the nominal value of the process parameter α j , j ∈ {1, 2}, and n is the number of state variables of the fuzzy control system.The Extended Symmetrical Optimum (ESO) method [25,26] is recommended to tune the two PI controller parameters because it ensures a compromise to the control system performance indices (percent overshoot, settling time, rise time, phase margin, etc.) of the linear control system by means of a single design parameter β within the recommended domain 1 < β ≤ 20.The PI tuning conditions specific to the ESO method are and the transfer function of the reference input filter F that contributes to performance enhancement is The T-S PI-FCs ensure further performance enhancement by design and tuning starting with the linear PI controllers.The typical structure and input membership functions of a T-S PI-FC are presented in Figure 2, where q −1 is the backward shift operator, TISO-FC is the Two Inputs-Single Output fuzzy controller block modeled by a nonlinear input-output static map, ∆e(t d ) is the increment of the control error, and ∆u(t d ) is the increment of the control signal.
The T-S PI-FCs ensure further performance enhancement by design and tuning starting with the linear PI controllers.The typical structure and input membership functions of a T-S PI-FC are presented in Figure 2, where 1 − q is the backward shift operator, TISO-FC is the Two Inputs-Single Output fuzzy controller block modeled by a nonlinear input-output static map, is the increment of the control error, and is the increment of the control signal.These two increments are obtained by discretizing the continuous-time PI controller by Tustin's method, which leads to the recurrent equation of the incremental discrete-time PI controller and its parameters where is the discrete time index and 0 > s T is the sampling period.
The TISO-FC block in the T-S PI-FC structure uses the weighted average method for defuzzification, and the SUM and PROD operators in the inference engine.The complete rule base of this block consists of the following two rules [21]:  (10) where the parameter η, , contributes to the reduction of the control system overshoot.
The application of the modal equivalence principle [30] leads to the tuning equation .] [ The optimization problems that ensure the sensitivity reduction with respect to the modifications of the process parameter j α are defined as , are the weighting parameters, with the subscript indicating the process parameters indicated in Equation ( 5), * ρ is the optimal parameter vector, i.e., the optimal value of These two increments are obtained by discretizing the continuous-time PI controller by Tustin's method, which leads to the recurrent equation of the incremental discrete-time PI controller and its parameters K P and µ [27-29] where t d ∈ Z, t d ≥ 0 is the discrete time index and T s > 0 is the sampling period.
The TISO-FC block in the T-S PI-FC structure uses the weighted average method for defuzzification, and the SUM and PROD operators in the inference engine.The complete rule base of this block consists of the following two rules [21] where the parameter η, 0 < η < 1, contributes to the reduction of the control system overshoot.The application of the modal equivalence principle [30] leads to the tuning equation and the parameters of the T-S PI-FC are grouped in the parameter vector ρ: The optimization problems that ensure the sensitivity reduction with respect to the modifications of the process parameter α j are defined as where γ α j , j ∈ {1, 2}, are the weighting parameters, with the subscript indicating the process parameters indicated in Equation ( 5), ρ * is the optimal parameter vector, i.e., the optimal value of the vector ρ, and D ρ is the feasible domain of ρ.The stability of the fuzzy control is recommended to be taken into consideration in order to set the domain D ρ , and useful approaches are reported in [23,[27][28][29][31][32][33][34][35][36].
The objective function J α j (ρ) is referred to as the weighted sum of the absolute value of the control error and of the squared output sensitivity function, but the state sensitivity functions can be included as well.Other objective functions can also be considered, inspired from several applications with classical and modern algorithms [37][38][39][40][41][42].
The state sensitivity models of the fuzzy control system with respect to α j are derived using Equation ( 6), with n = 4 for T-S PI-FCs.Accepting that the control signal u is changing at the discrete sampling intervals and the zero-order hold is included in the fuzzy control system structure for digital control, the number of four state variables results from the fact that the state variables x 3 and x 4 of the T-S PI-FC are defined in terms of [7].The expression of the state sensitivity model of the fuzzy control system with respect to the parameter k P of the family of nonlinear servo systems is [7] and the expression of the state sensitivity model of the fuzzy control system with respect to the parameter T Σ is [7] λ where the state variables of the T-S PI-FC are defined in terms of f TISO−FC is the nonlinear input-output map of TISO-FC, and the dynamics of the reference input filter are not accounted for.The subscript 0 in Equations ( 14) and ( 15) outlines the nominal value of the process parameter (as already shown in relation with Equation ( 6)) and also the nominal trajectory of the fuzzy control system (i.e., the nominal state variables) and the nominal expression of the input-output map of TISO-FC.

Grey Wolf Optimizer, Analysis, and Optimal Tuning Approach
This section presents the general description of GWO with the aim of introducing it later into the optimization problems described in the previous section.Based on this information a novel optimal tuning approach is defined in the last part of the section.

Grey Wolf Optimizer and Analysis
The operating mechanism of GWO according to its standard formulation given in [12] starts with the random initialization of the agents that comprise the wolf pack.A total number of N agents (i.e., grey wolves) is used, and each agent is assigned to a position vector X i (k) where is the position of ith agent in the f th dimension, f = 1...q, k is the index of the current iteration, k = 1...k max , and k max is the maximum number of iterations.
The search process specific to GWO continues with the exploration stage, which models the search for the prey.During this stage the positions of the top three agents, namely the alpha (α), beta (β), and delta (δ) agents, dictate the search pattern by diverging from other agents and converging to the prey, representing the optimal solution.
The exploitation stage models the attack on the prey.The top three agents constrain the other agents (the omega (ω) agents) to update their positions according to theirs.The following notations are used for the top three agent position vectors, i.e., the first three best solutions obtained at each iteration (or the alpha, beta, and delta solutions): the three vector solutions X α (k), X β (k), and X δ (k) are obtained by the following selection process: and they fulfill the condition A set of search coefficients is then defined: and the coefficients a f (k) are linearly decreased from 2 to 0 during the search process: The approximate distances between the current solution and the alpha, beta, and delta solutions, i.e., d i f The components (agents) of the updated alpha, beta, and delta solutions are computed as and they lead to the updated expressions of the agents' positions obtained as the arithmetic mean of the updated alpha, beta, and delta agents: The vector counterpart of Equation ( 25) is the formula that gives the update vector solution X i (k + 1): The GWO consists of the following steps that are formulated by the revision of the steps presented in [20,21]: Step 1.The initial random grey wolf population, represented by N agents' positions in the q-dimensional search space, is generated.The iteration index is initialized to k = 0 and the maximum number of iterations is set to k max .
Step 2. The performance of each member of the population of agents is evaluated by simulations and/or experiments conducted on the fuzzy control system.The evaluation leads to the objective function value by mapping the GWO onto the optimization problems using Step 3. The first three best solutions obtained so far, i.e., X α (k), X β (k), and X δ (k), are identified using Equation (19).
Step 5.The agents are moved to their new positions by computing X i (k + 1) in terms of Equations ( 24), (25), and (26).
Step 6.The updated vector solution X i (k + 1) ∈ D ρ is validated by checking the steady-state condition for the fuzzy control system with the T-S PI-FC parameter vector ρ = X i (k + 1), so far [7] where t d f is the final time moment.Theoretically t d f → ∞ as pointed out in Equation ( 13), but t d f takes practically a finite value in order to capture the transients in the fuzzy control system responses.The stability analysis of the fuzzy control system can be also employed with this regard.
Step 7. The iteration index k is incremented and the algorithm continues with step 2 until k max is reached.
Step 8.The algorithm is stopped and the final solution obtained so far is actually the solution to the optimization problems defined in Equation ( 13): The substitution of d i f l (k) taken from Equation (23) in Equation ( 24) leads to The computation of the arithmetic mean in Equation ( 30) for l ∈ {α, β, δ} and then the application of Equation (25) at the iterations k + 1 and k are organized as the following state-space equation: which indicates a nonlinear dynamic system.Such expressions are used in [43][44][45] for other metaheuristics as well (PSO and GSA) and are then associated with stability analyses.However, the nonlinearity in Equation ( 31) makes the GWO analysis more difficult in comparison with [43][44][45].
Equation ( 31) is then expressed in the equivalent form and since the right-hand term in Equation ( 32) is negative, this indicates that Using k = k max in Equation ( 22) and then Equation ( 21), the coefficients obtain the expressions Using Equation (34) in Equation ( 32) for k = k max , this gives the condition for the final values of the agents' positions In addition, Equation ( 35) used in combination with Equation (33) gives the following condition related to the dynamics of the agents' positions: Concluding, Equation (36) indicates that {x f i (k)} k=1...k max is a monotonic strictly decreasing sequence as this will guarantee the convergence of GWO.Therefore, it is recommended to initialize as large as possible values of x f i (1), f = 1...q, i = 1...N.The stability analysis of the dynamics of the agents' positions starts with the definition of the Lyapunov function candidate that fulfils the conditions The expression of the increment of V is Using Equation (36) in Equation (37), the sufficient condition ∆V(x f i (k)) < 0 for the global asymptotic stability of the origin of dynamics of the agents' positions is (40) Since the right-hand term in Equation ( 41) is strictly positive, Equation ( 41) leads to the first sufficient condition for the agents' positions: which also implies that Using Equations ( 42) and ( 43) and the properties of the modulus, the right-hand term in Equation ( 41) is upper bounded in terms of To sum up, in order to fulfill the condition ∆V(x f i (k)) < 0 for the global asymptotic stability of the origin of dynamics of the agents' positions, it is sufficient to fulfill the following condition, as it results from Equations ( 41) and ( 44) and the transitive property of inequalities: which is equivalent to Since Equation s (21) and (22) the left-hand term in Equation ( 46) fulfills the condition Summing up, the stability analysis should be re-organized in order to check the sufficient stability conditions of Equations ( 42) and ( 46) at each iteration.This stability check is recommended to be carried out in step 6 in order to validate the next vector solution by dropping out those agents that do not fulfill the conditions of Equations ( 42) and (46).
Moreover, Equation (42) indicates that GWO is recommended to solve optimization problems where the variables take positive values.Such problems occur in the optimal tuning of the parameters of several linear and fuzzy controllers or of the positive parameters specific to linear and nonlinear models, including fuzzy ones.machines with various configurations, thus parallelizing the computations in order to reduce the execution time.
The weighting parameters presented in Equation ( 13) have been defined to satisfy a ratio of {0.1, 1, 10} of the initial values of the terms presented in the sums.The following values have been used: (γ k P ) 2 ∈ {0.006858, 0.06858, 0.6858}, (γ T Σ ) 2 ∈ {0.0066695, 0.066695, 0.66695}.( The upper limit of the sums in Equation ( 13) has been set to t d f = 2000 instead of t d f → ∞ , which corresponds to the evaluation of the objective functions by experiments conducted on a time horizon of 20 s.The vector variables ρ of the objective functions have been initialized and are then validated to belong to the feasible (search) domain The dynamic regimes characterized by r = r 0 = 40 rad step type modification applied to the reference input and disturbance input d = d 0 = 0 have been used.The integral term in the T-S PI-FC structure carries out the disturbance rejection.
The parameters of GWO have been defined based on the designers' experience and aim to achieve an acceptable balance between convergence and computational resources as N = 20 and k max = 100.
The optimal controller parameters and the associated minimum objective function values of the two optimization problems defined in Equation ( 13) are presented in Tables 1 and 2. The quality of GWO is analyzed through three sets of performance indices that try to quantify the use of allocated resources and asses the search process.The first set of performance indexes used in the analysis of the solution quality is the average value of the objective functions J k P and J T Σ obtained by running several independent runs of GWO and are referred to as Avg(J k P min ) and Avg(J T Σ min ): where J k P min and J T Σ min are the values obtained by running GWO, N best is the number of best values (i.e., the smallest values) obtained for the objective functions, and the superscript j, j = 1...N best , indicates the values of the objective functions obtained by one of the best N best runs of GWO, namely by the run j, j = 1...N best .The value N best = 5 is considered in the context of the current analysis.The values of this performance index have already been presented in Tables 1 and 2 for a fair assessment.
For the second set of performance indices, viz. the convergence speed c s , the number of evaluations of the objective functions until the minimum value is found is looked at.This index is presented as an average value of GWO runs for the first performance index.
The third set of performance indices is introduced with the aim of evaluating the recall of the search process.This set of indices, referred to as the accuracy rate a r , is defined in terms of the percentage of standard deviation of the objective function value with respect to the average value: a r = StDev % (J min ) = 100StDev(J min )/Avg(J min ), StDev(J min ) = GWO has been compared with the other two metaheuristics to solve the same optimization problems defined in Equation (13).The resulting parameters from [7,9] for the two metaheuristics GSA and PSO that were considered based on their unsupervised characteristics show similar values to the ones already presented.
Another set of experimental results is presented in Table 3 as the values of the performance indices c s and a r .These results show that no algorithm has a clear advantage as the values of these sets of performance indices are close.However, GWO is the overall best metaheuristics from the point of view of c s , and PSO is the overall best one as far as a r is concerned.The same conclusion has been drawn in [21] for a different optimization problem.Typical responses for these fuzzy control systems are presented in [21].

Conclusions
This paper has proposed an easily understandable optimization algorithm (GWO) applied to the optimal tuning of the parameters of T-S PI-FCs, showing good results in comparison with the other two metaheuristics in the fuzzy control of a family of nonlinear servo systems.
The convergence and stability analysis of GWO indicate that GWO is mainly dedicated to optimization problems where the variables are strictly positive.However, the operating mechanism of GWO based on the use of the best three agents makes it useful in other optimization problems as well without guaranteeing the convergence.
GWO is advantageous over other metaheuristics because of the reduced number of random parameters and user-selected parameters.This illustrates its potential in solving various optimization problems with a reduced user experience and fair comparison with similar metaheuristics.
The performance obtained by GWO is encouraging and comparable to PSO and GSA.Future research will be focused on the performance improvement of GWO by its inclusion in adaptive and hybrid optimization algorithms.

Figure 1 .
Figure 1.Structure of the fuzzy control system, y r e − = 1 .

Figure 1 .
Figure 1.Structure of the fuzzy control system, e = r 1 − y.

Figure 2 .
Figure 2. Structure and input membership functions of the Takagi-Sugeno proportional-integralfuzzy controller.
the parameters of the T-S PI-FC are grouped in the parameter vector ρ :

Figure 2 .
Figure 2. Structure and input membership functions of the Takagi-Sugeno proportional-integralfuzzy controller.
distributed random numbers within 0

Table 1 .
Weighting parameter and controller parameters for the minimization of J k P .

Table 2 .
Weighting parameter and controller parameters for the minimization of J T Σ .

Table 3 .
Average values of c s and a r .