A Mixed-Strategy-Based Whale Optimization Algorithm for Parameter Identiﬁcation of Hydraulic Turbine Governing Systems with a Delayed Water Hammer Effect

: For solving the parameter optimization problem of a hydraulic turbine governing system (HTGS) with a delayed water hammer (DWH) effect, a Mixed-Strategy-based Whale Optimization Algorithm (MSWOA) is proposed in this paper, in which three improved strategies are designed and integrated to promote the optimization ability. Firstly, the movement strategies of WOA have been improved to balance the exploration and exploitation. In the improved movement strategies, a dynamic ratio based on improved JAYA algorithm is applied on the strategy of searching for prey and a chaotic dynamic weight is designed for improving the strategies of bubble-net attacking and encircling prey. Secondly, a guidance of the elite’s memory inspired by Particle swarm optimization (PSO) is proposed to lead the movement of the population to accelerate the convergence speed. Thirdly, the mutation strategy based on the sinusoidal chaotic map is employed to avoid prematurity and local optimum points. The proposed MSWOA are compared with six popular meta-heuristic optimization algorithms on 23 benchmark functions in numerical experiments and the results show that the MSWOA has achieved signiﬁcantly better performance than others. Finally, the MSWOA is applied on parameter identiﬁcation problem of HTGS with a DWH effect, and the comparative results conﬁrm the effectiveness and identiﬁcation accuracy of the proposed method.


Introduction
The hydraulic turbine governing system (HTGS) is the core control system of hydroelectric generating units (HGUs), which undertakes the tasks of load and frequency adjustment. The HTGS is a non-linear and non-minimum phase system [1][2][3][4][5][6][7]. Traditional identification methods, such as least squares method [8], input response method [9], and maximum likelihood estimate (MLE) [10], have been applied in parameter identification of HTGS in the past years, but there are too many restrictions on these methods. For example, the least squares method demands enough system input, and MLE easily gets trapped in local optima points. In addition, most of those methods are not global optimization methods and not suitable for parameter identification of nonlinear or complicated systems. Therefore, it is difficult to apply these methods in HTGS parameter identification if nonlinearities have been considered. To conquer these problems, identification methods based on meta-heuristic algorithms have been developed in recent years, which treat the problem of parameter identification as an optimization problem [11]. Because meta-heuristic algorithms are global optimization methods, In this phase, whales in the population update their positions using Equations (1) and (2). In Equations (1) and (2), A and C are calculated with Equations (3) and (4) when |A| is greater than 1. → r is a vector which is composed of a random number in [0,1], and a is evaluated according to Equation (5) and it varies within [0,2]: In this phase, each whale of the population updates its position with Equations (6) and (7). A and C are calculated with Equations (3) and (4): 2.1.3. Modeling the Course of Bubble-Net Attacking (Getting) Prey In this phase, each whale of the population updates its position using Equations (8) and (9). The flowchart of WOA is shown in Figure 1. The term b is a constant to which a positive number is assigned in this paper, l is a random number in [−1.1].
Energies 2018, 11 In this phase, each whale of the population updates its position with Equations (6) and (7). A and C are calculated with Equations (3) and (4): X ⃗ (t + 1) = X * ⃗ (t) − A ⃗ • D ⃗ 2.1.3. Modeling the Course of Bubble-Net Attacking (Getting) Prey In this phase, each whale of the population updates its position using Equations (8) and (9). The flowchart of WOA is shown in Figure 1. The term b is a constant to which a positive number is assigned in this paper, l is a random number in [−1.1].

The Mixed-Strategy Based WOA
The mixed strategy, which is composed of three improvements, not only enhances the algorithm comprehensively but also balances exploration and exploitation. Compared with the standard WOA in Figure 1, how the improvements strategies enhance MSWOA is illustrated in Figure 2. The details are described as follows:

The Mixed-Strategy Based WOA
The mixed strategy, which is composed of three improvements, not only enhances the algorithm comprehensively but also balances exploration and exploitation. Compared with the standard WOA in Figure 1, how the improvements strategies enhance MSWOA is illustrated in Figure 2. The details are described as follows:

Improvement 1: Hybrid Movement Strategy
The basic movement strategy of WOA is composed of three strategies, namely the strategy of encircling prey, the strategy of searching for prey and the strategy of bubble-net attacking the prey. To enhance MSWOA comprehensively, three different improvements are applied on the three movement strategies of WOA, respectively: (1) The dynamic ratio for strategy of searching for prey.
In standard WOA, a random probability value named p was employed to allow WOA to effectively transit between exploration and exploitation. In this paper, besides p, another random probability value named p1 and a dynamic ratio named q1 which is defined as Equation (10), are applied on strategy of searching for prey: where c1 is a coefficient in [0,1]; in this paper c1 = 0.5 and q1 varies in [0,0.5].

Improvement 1: Hybrid Movement Strategy
The basic movement strategy of WOA is composed of three strategies, namely the strategy of encircling prey, the strategy of searching for prey and the strategy of bubble-net attacking the prey. To enhance MSWOA comprehensively, three different improvements are applied on the three movement strategies of WOA, respectively: (1) The dynamic ratio for strategy of searching for prey.
In standard WOA, a random probability value named p was employed to allow WOA to effectively transit between exploration and exploitation. In this paper, besides p, another random probability value named p 1 and a dynamic ratio named q 1 which is defined as Equation (10), are applied on strategy of searching for prey: where c 1 is a coefficient in [0,1]; in this paper c 1 = 0.5 and q 1 varies in [0,0.5].
In Figure 2, q 1 is compared with p 1 in each iteration. If q 1 < p 1 , the improved JAYA algorithm is employed to update the agent randomly selected, otherwise we use a strategy of randomly selecting an agent. The usage of q 1 and p 1 is helpful to search the space thoroughly.
(2) Strategy of searching for prey based on an improved JAYA algorithm Because the dynamic ratio q 1 and p 1 are applied on strategy of searching for prey, it makes MSWOA transit between the improved JAYA algorithm and the strategy of randomly selecting an agent, which are described as Equations (11) and (12), respectively, in the phase of searching for prey. Compared with the JAYA algorithm [39], Equation (11) is proposed, which is named the improved JAYA algorithm. In Equation (11), c 2 and c 3 are employed to enlarge the search space. This is useful to not only enhance global search capability but also to avoid prematurity.
In each iteration, q 1 is compared with p 1 . If q 1 < p 1 , − −− → X rand is upgraded according to Equation (11), otherwise − −− → X rand is upgraded according to Equation (12). The step of q 1 < p 1 originates from the acceptance/rejection step in the Markov Chain Monte Carlo Techniques which is helpful to remove some samples by some mechanism [39,40].
(3) Strategies of encircling prey and bubble-net attacking based on chaotic dynamic weight In this section, the logistic chaotic map is applied on the strategies of encircling prey and bubble-net attacking. The logistic chaotic map, which is ergodic, sensitive, non-repetitive and helps MSWOA to avoid local optimum points or prematurity. ω(t), which denotes value of the logistic chaotic dynamic weight in the t-th iteration, is evaluated using Equation (13) [41]. The initial value, ω(1), is a random number in [0,1]: Each position of the whale agents is updated by Equation (14) in the strategy of encircling prey and is updated with Equation (15) in the strategy of bubble-net attacking while each whale agent position is updated with Equations (7) and (8) in the standard WOA. In Equations (14) and (15), ω(t) is the chaotic dynamic weight. The global best position so far − −− → X gbest and the best position in current iteration → X * are both taken into account, respectively. The memory of previous iterations is helpful for MSWOA to enhance the global optimum finding capability and avoid local optimum points and prematurity: A guidance of the elite's memory inspired by PSO is applied on the movement of the whale agents of the population. In each iteration, the movement strategy is introduced to update each position of the whale agents of the population after searching space by using the strategies proposed above. The strategy is described as Equations (16) and (17). This is helpful to improve search capability: where ω 1 is a constant in [0,1], c 11 and c 12 are coefficients in [0,2],r 8 , r 9 are random numbers in [0,1], − −− → X pbest denotes the global best position so far, − −− → X gbest denotes the personal best position so far.

Improvement 3: Mutation Operator Based on The Sinusoidal Chaotic Mutation
The chaotic mutation operator, which can overcome the shortcoming of local convergence or prematurity, is one of the best approaches to find the best solution through thoroughly evaluating the search space. In this section, the sinusoidal chaotic map is selected as the chaotic mutation operation. The mathematical function of sinusoidal chaotic map can be written as Equation (18) [41]: where a1 is chaotic constant and a1 = 2.3, N denotes the total number of members of the population. At each iteration, each position of the whale agents of population is brought into Equation (18) to achieve a new position. The values of objection function of each position and its new position achieved by chaotic mutation operation are calculated. Between any position and its new position achieved by Equation (18), the one whose value of objection function is better remains in the population while the other is abandoned.

Procedure of Mixed-Strategy WOA
The pseudo-code of the algorithm is shown as follows (Algorithm 1): Energies 2018, 11, x FOR PEER REVIEW 9 of 30 Algorithm 1. Pseudo-Code. 1. Initialize the agents population X, and x (i = 1,2…N，N is total number of all agents) denotes position of the i-th agent of X; 2. While (t < tmax) For x to x (1) Evaluate objective function value of each agent of the population and select the worst solution and the best solution and update X * , X , X in the current iteration; (2) Update A, C, a, l, p and calculate the chaotic weight ω in the current iteration; (3) if 1 p < 0.5 (p and p1are random in [0,1], q1 = c1·(1 − t/tmax)) if 2 |A| ≥ 1 if 3 q1 < p1 select X in Equation (11) and update the position of each agent of population as in Equations (1) and (2) else if 3 q1 > p1 pdate X in Equation (12) and update the position of each agent of population as in Equations (1)  Check whether the position of any agent of population is out of the boundaries. t = t + 1 end while 3. end

Experiments and Result Discussion
In this section, PSO [10], GSA [13], ALO [15], WOA [19], the Enhanced Whale Optimization Algorithm (EWOA) [34] and CWOA [36] are selected to be compared with the MSWOA proposed above on 23 benchmark functions which are depicted in detail in Tables 2-4. In all tests on 23 benchmark functions, the population size is 30, and the total number of iterations is 500 for all algorithms proposed in Table 5. To prove that MSWOA outperforms other algorithms in Table 5, the

Experiments and Result Discussion
In this section, PSO [10], GSA [13], ALO [15], WOA [19], the Enhanced Whale Optimization Algorithm (EWOA) [34] and CWOA [36] are selected to be compared with the MSWOA proposed above on 23 benchmark functions which are depicted in detail in Tables 2-4. In all tests on 23 benchmark functions, the population size is 30, and the total number of iterations is 500 for all algorithms proposed in Table 5. To prove that MSWOA outperforms other algorithms in Table 5, the experimental results of MSWOA are compared with those of other algorithms by two methods, namely the Wilcoxon' test and the box and whisker method.

Function
Dimension Range fmin

Experiments Setting and Benchmark Function
Each algorithm of the seven algorithms proposed in Table 5 is tested on 23 benchmark functions which can be divided into three groups in Tables 2-4, namely unimodal test functions (F1-F7), multimodal test functions (F8-F13) and multimodal test functions with fixed dimensions (F14-F23). Generally, F1-F7 are employed to calculate the exploitation ability of algorithms benchmarked by benchmark functions and F8-F23 are used to calculate exploration ability. In Tables 2-4, the first column is the expression of function, the second column is dimension of function, the third column is the domain of variable and the last column is the standard minimum value of the function.

Algorithm Parameter Settings
PSO12 c 1 = 2, c 2 = 2, inertia weight ω = 1, inertia weight damping ratio ω2 = 0.99 GSA13 G0 = 100, α = 20 ALO16 All parameter settings are recommended in [15] WOA20 a decrease linearly from 2 to 0 CWOA37 a decrease linearly from 2 to 0, other settings are same to [19] EWOA35 a decrease linearly from 2 to 0, c = 0.3, other settings are same to [19] MSWOA a decrease linearly from 2 to 0, q 1 = 0.5(1 − t/t max ), inertia weight ω = 1, inertia weight damping ratio ω2 = 0.99, other settings are same to WOA [19] Because the meta-heuristic algorithms proposed in Table 5 are stochastic, the experiment of any kind of benchmark function is repeated 20 times. Meanwhile the mean value and the standard deviation value of 20 repeated experiments are saved as test results. The population size of each algorithm is set at 30 and total number of iterations is 500, other parameter settings are listed in Table 5. In addition, all benchmark functions experiments are run on MATLAB R2016a (R2016a, Math Works, Natick, MA, USA).

Results Analysis Based on Statistical Test Methods
All meta-heuristic algorithms presented in Table 5 are applied on the 23 benchmark functions in comparative experiments. For a fair competition, the experiments are repeated 20 times, while the average values and standard deviation are calculated and presented in Table 6. To analyze these results reasonably, two statistical test methods, namely the Wilcoxon's and box and whisker tests, are adopted. The specific results are discussed as follows.

The Wilcoxon's Test
The Wilcoxon's test, a nonparametric statistical test, is applied to test which is the more significant between two algorithms. The Wilcoxon's test is divided into the single-problem statistical analysis and multiple-problem statistical analysis.
In the single-problem statistical analysis, any two algorithms mentioned in Table 5 are selected to be compared with each other on the same benchmark function, F1 to F23. In Table 6, "Ave" represents a mean value of a group of experiments for 20 times on one benchmark function. Hence the "Ave" is used as the test sample for the Wilcoxon's test. If algorithm A obtains a smaller "Ave" than algorithm B, algorithm A is considered to have a better significance level than algorithm B. In other words, algorithm A is superior to algorithm B by the Wilcoxon signed-rank test at alfa = 0.05. In Table 7, "W/T/L" represents three relations, namely Win, Tie and Lose. "W" means MSWOA obtains a smaller "Ave" than another algorithm, "L" means MSWOA obtains a bigger "Ave" than another algorithm, "T" means MSWOA obtains an equal "Ave" to another algorithm. In the test, the total number of "W/T/L" is counted in the last row of Table 7.  Tie  Tie  Tie  Tie  Tie  Tie  F17  Tie  Tie  Tie  Tie  Tie  Tie  F18  Tie  Tie  Tie  MSWOA  Tie  Tie  F19  MSWOA  GSA  MSWOA  MSWOA  MSWOA  MSWOA  F20 MSWOA The multiple-problem statistical analysis is used to compare two algorithms in the several similar benchmark functions while the single-problem statistical analysis is used to compare two algorithms in the same one benchmark function. In Table 6, the results of each meta-heuristic algorithm tested on 23 benchmark functions are listed. All mean values of the results of each algorithm tested on the 23 benchmark functions are treated as the input vectors of the Wilcoxon test.
In Table 8, the test results of multiple-problem statistical analysis by the Wilcoxon signed-rank test are listed. The R+, R− and the p-values can be evaluated by the Wilcoxon signed-rank test applied in multiple-problem statistical analysis. If MSWOA obtained the higher R+ than R− values in any one pair of comparison and the p-value is less than 0.05, it means that MSWOA is statistically significantly better than another. From Tables 7 and 8, we can conclude that MSWOA has a better performance in a statistical manner than the other six algorithms, especially for WOA, CWOA and EWOA. Based on WOA, MS-WOA is improved greatly on the capability of stability, exploitation and exploration and avoidance of local optimum points and prematurity.

Box and Whisker
The stability is a concept used to evaluate an algorithm by checking the randomness of solutions. In addition to the standard deviation statistical index, the box and whisker plot is also effective in estimation of algorithm stability. The box and whisker methods used in this section to evaluate the distribution of results of the repeated 20 runs on benchmark functions F1 to F23. A box and whisker plot can illustrate the variation in a set of data and provide more information about the data. It provides five statistics, namely, the minimum value, the second quartile, the median value, the third quartile and the maximum value. The second quartile is the value below which the lower 25% of the data are contained. Third quartile is the value above which the upper 25% of the data are contained.
In Figure 3, the distributions of the objective function values of the optimal solutions of multiple runs are revealed by box and whisker plots, while the proposed algorithm is compared with the other six algorithms. According to Figure 3, it is obvious that box and whisker of MSWOA are all nearly horizontal lines except for (h) and (t). It is shown that the variation of the values of the best solution obtained from MSWOA after 20 times are very small and the stability of MSWOA is superior to that of the other algorithms. Compared with homogenous algorithms, like WOA, CWOA and EWOA, the stability of MSWOA has been greatly improved. The multiple-problem statistical analysis is used to compare two algorithms in the several similar benchmark functions while the single-problem statistical analysis is used to compare two algorithms in the same one benchmark function. In Table 6, the results of each meta-heuristic algorithm tested on 23 benchmark functions are listed. All mean values of the results of each algorithm tested on the 23 benchmark functions are treated as the input vectors of the Wilcoxon test.
In Table 8, the test results of multiple-problem statistical analysis by the Wilcoxon signed-rank test are listed. The R+, R− and the p-values can be evaluated by the Wilcoxon signed-rank test applied in multiple-problem statistical analysis. If MSWOA obtained the higher R+ than R− values in any one pair of comparison and the p-value is less than 0.05, it means that MSWOA is statistically significantly better than another. From Tables 7 and 8, we can conclude that MSWOA has a better performance in a statistical manner than the other six algorithms, especially for WOA, CWOA and EWOA. Based on WOA, MS-WOA is improved greatly on the capability of stability, exploitation and exploration and avoidance of local optimum points and prematurity.

Box and Whisker
The stability is a concept used to evaluate an algorithm by checking the randomness of solutions. In addition to the standard deviation statistical index, the box and whisker plot is also effective in estimation of algorithm stability. The box and whisker methods used in this section to evaluate the distribution of results of the repeated 20 runs on benchmark functions F1 to F23. A box and whisker plot can illustrate the variation in a set of data and provide more information about the data. It provides five statistics, namely, the minimum value, the second quartile, the median value, the third quartile and the maximum value. The second quartile is the value below which the lower 25% of the data are contained. Third quartile is the value above which the upper 25% of the data are contained.
In Figure 3, the distributions of the objective function values of the optimal solutions of multiple runs are revealed by box and whisker plots, while the proposed algorithm is compared with the other six algorithms. According to Figure 3, it is obvious that box and whisker of MSWOA are all nearly horizontal lines except for (h) and (t). It is shown that the variation of the values of the best solution obtained from MSWOA after 20 times are very small and the stability of MSWOA is superior to that of the other algorithms. Compared with homogenous algorithms, like WOA, CWOA and EWOA, the stability of MSWOA has been greatly improved.

Performance Comparison
According to the test results in Tables 6 to 8, a conclusion is drawn that the MSWOA algorithm proposed in the paper has been improved significantly on its performance compared to the other algorithms in Table 5. In this subsection, the performance of MSWOA will be analyzed further based on the results of 23 benchmark functions.

Analysis of Test Results of F1-F7
Based on character of benchmark functions, unimodal test functions (F1-F7) are applied to test the capacity of exploitation of algorithms while the multimodal test function and multimodal test functions with fixed dimension (F8-F23) are employed to test the capacity of exploration of algorithms. In Table 6, its manifest that the MSWOA gives a superior test result than the other algorithms in all unimodal test functions except F5, while EWOA gives a better test result than the others on F5, MSWOA outperforms other algorithms on both the final test results and convergence rate. In Figure 4a-w, the convergence curves of PSO, GSA, ALO, WOA, CWOA, EWOA and MSWOA are the average curves with mean values tested on the 23 benchmark functions 20 times, respectively. In Figure 4a-g, it is clear that the MSWOA has the fastest convergence rate to search for the global optimization point. In the unimodal test functions test, the convergence rates of the algorithms are more vital than the final results. As for Figure 3a-g, we can find that MSWOA outperforms the others in the box and whisker test results. That means MSWOA is more stable than PSO, GSA, ALO, WOA, CWOA and EWOA for unimodal test functions (F1-F7).

Performance Comparison
According to the test results in Tables 6-8, a conclusion is drawn that the MSWOA algorithm proposed in the paper has been improved significantly on its performance compared to the other algorithms in Table 5. In this subsection, the performance of MSWOA will be analyzed further based on the results of 23 benchmark functions.

Analysis of Test Results of F1-F7
Based on character of benchmark functions, unimodal test functions (F1-F7) are applied to test the capacity of exploitation of algorithms while the multimodal test function and multimodal test functions with fixed dimension (F8-F23) are employed to test the capacity of exploration of algorithms. In Table 6, its manifest that the MSWOA gives a superior test result than the other algorithms in all unimodal test functions except F5, while EWOA gives a better test result than the others on F5, MSWOA outperforms other algorithms on both the final test results and convergence rate. In Figure 4a-w, the convergence curves of PSO, GSA, ALO, WOA, CWOA, EWOA and MSWOA are the average curves with mean values tested on the 23 benchmark functions 20 times, respectively. In Figure 4a-g, it is clear that the MSWOA has the fastest convergence rate to search for the global optimization point. In the unimodal test functions test, the convergence rates of the algorithms are more vital than the final results. As for Figure 3a-g, we can find that MSWOA outperforms the others in the box and whisker test results. That means MSWOA is more stable than PSO, GSA, ALO, WOA, CWOA and EWOA for unimodal test functions (F1-F7).

Performance Comparison
According to the test results in Tables 6 to 8, a conclusion is drawn that the MSWOA algorithm proposed in the paper has been improved significantly on its performance compared to the other algorithms in Table 5. In this subsection, the performance of MSWOA will be analyzed further based on the results of 23 benchmark functions.

Analysis of Test Results of F1-F7
Based on character of benchmark functions, unimodal test functions (F1-F7) are applied to test the capacity of exploitation of algorithms while the multimodal test function and multimodal test functions with fixed dimension (F8-F23) are employed to test the capacity of exploration of algorithms. In Table 6, its manifest that the MSWOA gives a superior test result than the other algorithms in all unimodal test functions except F5, while EWOA gives a better test result than the others on F5, MSWOA outperforms other algorithms on both the final test results and convergence rate. In Figure 4a-w, the convergence curves of PSO, GSA, ALO, WOA, CWOA, EWOA and MSWOA are the average curves with mean values tested on the 23 benchmark functions 20 times, respectively. In Figure 4a-g, it is clear that the MSWOA has the fastest convergence rate to search for the global optimization point. In the unimodal test functions test, the convergence rates of the algorithms are more vital than the final results. As for Figure 3a-g, we can find that MSWOA outperforms the others in the box and whisker test results. That means MSWOA is more stable than PSO, GSA, ALO, WOA, CWOA and EWOA for unimodal test functions (F1-F7).

Analysis of Test Results of F8-F13
The ability which makes an algorithm escape from local optimum points and locate a global optimization point, can be tested by multimodal functions (F8-F13). When an algorithm is run on multimodal functions (F8-F13), the test results can reflect the ability better than the convergence rate.
In Table 6, it is clear that MSWOA outperforms the other algorithms in F9, F10, F11 and F12 in the test results. MSWOA has the second best performance in F13. Although GSA has a better test result than others in F13, MSWOA outperforms WOA, CWOA and EWOA on both the final test results and convergence rate. In Figure 3h-m, it is found that MSWOA outperforms the others in the box and whisker test results, which means the stability of MSWOA is superior to those of PSO, GSA, ALO, WOA, CWOA and EWOA for multimodal functions (F8-F13).

Analysis of Test Results of F14-F23
The ability of an algorithm to avoid a small quantity of local optima can be tested by multimodal functions with fixed dimensions (F14-F23). In Table 6, it is clear that MSWOA outperforms other the algorithms in F14 to F23 apart from F20 in the test results. MSWOA achieves the second best performance in F20. Though ALO has a better test result than the others in F20, MSWOA outperforms WOA, CWOA and EWOA on both the final test results and convergence rate. In Figure 3n-w, it is found that MSWOA outperforms others in test results in F14 to F23 apart from F20. That means the stability of MSWOA is superior to those of PSO, GSA, ALO, WOA, CWOA and EWOA for multimodal functions (F14-F23) apart from F20. In Figure 3t, the stability of MSWOA is only inferior

Analysis of Test Results of F8-F13
The ability which makes an algorithm escape from local optimum points and locate a global optimization point, can be tested by multimodal functions (F8-F13). When an algorithm is run on multimodal functions (F8-F13), the test results can reflect the ability better than the convergence rate.
In Table 6, it is clear that MSWOA outperforms the other algorithms in F9, F10, F11 and F12 in the test results. MSWOA has the second best performance in F13. Although GSA has a better test result than others in F13, MSWOA outperforms WOA, CWOA and EWOA on both the final test results and convergence rate. In Figure 3h-m, it is found that MSWOA outperforms the others in the box and whisker test results, which means the stability of MSWOA is superior to those of PSO, GSA, ALO, WOA, CWOA and EWOA for multimodal functions (F8-F13).

Analysis of Test Results of F14-F23
The ability of an algorithm to avoid a small quantity of local optima can be tested by multimodal functions with fixed dimensions (F14-F23). In Table 6, it is clear that MSWOA outperforms other the algorithms in F14 to F23 apart from F20 in the test results. MSWOA achieves the second best performance in F20. Though ALO has a better test result than the others in F20, MSWOA outperforms WOA, CWOA and EWOA on both the final test results and convergence rate. In Figure 3n-w, it is found that MSWOA outperforms others in test results in F14 to F23 apart from F20. That means the stability of MSWOA is superior to those of PSO, GSA, ALO, WOA, CWOA and EWOA for multimodal functions (F14-F23) apart from F20. In Figure 3t, the stability of MSWOA is only inferior to that of GSA which gains the best value while the stability of MSWOA is superior to others including PSO, ALO, WOA, CWOA and EWOA.

Parameter Identification of the Hydraulic Turbine Governing System
In this section, the MSWOA is applied to solve the parameter identification problem of a complicated HTGS with a delayed water hammer effect. The mathematical model of the HTGS is established, and then the parameter identification methodology as well as simulation experiments are conducted.

Model of Hydraulic Turbine Governing System with a Delay Water Hammer Effect
A HTGS is composed of speed governor, penstock, hydraulic turbine and generator [1,4,25,26], as shown in Figure 5. In Figure 5, X denotes the speed of the generator, Y denotes the position of the main servomotor, H denotes the water head, Q denotes the flow and M denotes the moment. In this paper, the penstock system model is a hyperbolic tangent function which is different from the rigid water hammer equation and elastic water hammer equation. The hyperbolic tangent function can describe the dynamic process of the penstock system in the most precise way. In Figure 6, each part of the HTGS model is described. to that of GSA which gains the best value while the stability of MSWOA is superior to others including PSO, ALO, WOA, CWOA and EWOA.

Parameter Identification of the Hydraulic Turbine Governing System
In this section, the MSWOA is applied to solve the parameter identification problem of a complicated HTGS with a delayed water hammer effect. The mathematical model of the HTGS is established, and then the parameter identification methodology as well as simulation experiments are conducted.

Model of Hydraulic Turbine Governing System with a Delay Water Hammer Effect
A HTGS is composed of speed governor, penstock, hydraulic turbine and generator [1,4,25,26], as shown in Figure 5. In Figure 5, X denotes the speed of the generator, Y denotes the position of the main servomotor, H denotes the water head, Q denotes the flow and M denotes the moment. In this paper, the penstock system model is a hyperbolic tangent function which is different from the rigid water hammer equation and elastic water hammer equation. The hyperbolic tangent function can describe the dynamic process of the penstock system in the most precise way. In Figure 6, each part of the HTGS model is described.

PID Controller Model
The PID controller is one part of the hydraulic turbine governor. The PID controller model is described as follows: to that of GSA which gains the best value while the stability of MSWOA is superior to others including PSO, ALO, WOA, CWOA and EWOA.

Parameter Identification of the Hydraulic Turbine Governing System
In this section, the MSWOA is applied to solve the parameter identification problem of a complicated HTGS with a delayed water hammer effect. The mathematical model of the HTGS is established, and then the parameter identification methodology as well as simulation experiments are conducted.

Model of Hydraulic Turbine Governing System with a Delay Water Hammer Effect
A HTGS is composed of speed governor, penstock, hydraulic turbine and generator [1,4,25,26], as shown in Figure 5. In Figure 5, X denotes the speed of the generator, Y denotes the position of the main servomotor, H denotes the water head, Q denotes the flow and M denotes the moment. In this paper, the penstock system model is a hyperbolic tangent function which is different from the rigid water hammer equation and elastic water hammer equation. The hyperbolic tangent function can describe the dynamic process of the penstock system in the most precise way. In Figure 6, each part of the HTGS model is described.

PID Controller Model
The PID controller is one part of the hydraulic turbine governor. The PID controller model is described as follows:

PID Controller Model
The PID controller is one part of the hydraulic turbine governor. The PID controller model is described as follows: where σ(s) is the Laplace transform of output of the PID controller; c(s) is the Laplace transform of a given speed; x(s) is the Laplace transform of the speed; b p is a permanent transition coefficient; y 1 (s) is the position of the auxiliary servomotor; T d is the differential time; K p is the proportional gain, K i is the integral gain and K d is the differential gain.

Servomechanism Model
The servomechanism is another part of the hydraulic turbine governor. The model is described as follows: where y(s) is the position of the main servomotor; T y1 is the response time constant of the auxiliary servomotor; T y is the response time constant of the main servomotor, y 0 is the initial value of y(s).

Hydraulic Turbine Model
When the HTGS is working, the dynamic process of the hydraulic turbine system can't be obtained by experiments or model tests, but if the speed of the hydraulic turbine fluctuates in a small range, a linear model of the hydraulic turbine can be given by Equation (21), which includes the express moment and flow characteristics, and is applied to depict the characteristics of a hydraulic turbine: m t (s) = e x ·x(s) + e y ·y(s) + e h ·h(s) q(s) = e qx ·x(s) + e qy ·y(s) + e qh ·h(s) (21) where m t (s) is the Laplace transform of the moment of the hydraulic turbine, q(s) is the Laplace transform of the water flow of the hydraulic turbine, h(s) is the water head of the hydraulic turbine, e x , e y , e h , e qx , e qy , e qh are the transfer coefficients of the hydraulic turbine which are obtained from the comprehensive characteristic curve of the hydraulic turbine. Usually we can obtain the transfer coefficients of some working point to build the hydraulic turbine model. The detailed instructions of Equation (21) as well as the calculation of those transfer coefficients may be found in [42].

Penstock System Model
The elastic water hammer and rigid water hammer models are main expressions of the penstock system models and are frequently used to describe the characteristics of a penstock system. In recent studies, many novel expressions of elastic water hammer [27,28] and rigid water hammer [2] models are proposed. Compared with a hyperbolic tangent function water hammer model as given by Equation (22), those models only keep parts of a hyperbolic tangent function water hammer model. This inevitably reduces the accuracy of the penstock system model. In fact, the hyperbolic tangent function water hammer model is the most proximate to the real penstock system and can describe the dynamic process of the penstock system in the most precise way. Therefore, the hyperbolic tangent function water hammer model given by Equation (22) is employed as the penstock system model in this paper. The transfer function of the penstock system can be expressed as: where h w is the pipeline characteristic coefficient, T r is the reflection time of the water hammer pressure wave. The dynamic equation of the synchronous generator is described as Equation (23) in which m g (s) means the load disturbance. When the HTGS is working under load conditions, m g (s) = 0 and e g = 0. When HTGS is working under no-load conditions, m g (s) = 0 and e g = 0.
x(s) m t (s) − m g (s) = 1 T a s + e n (23) where T a is the inertial time constant of the generator, e g is the adjustment coefficient of the generator and e n = e g − e x .

Identification Strategy Based on MSWOA
The proposed MSWOA is employed to identify the HTGS parameters. In the MSWOA identification strategy, firstly, an input is employed to activate the real system and the system to be identified. The outputs of the real system and the system to be identified are input into Equation (24) to evaluate the value of C OF (θ).
where θ is six-dimensional vector and θ = T y1 T y h w T r T a e g , z = x y m t is the output of the real system,ẑ = xŷm t is the output of the identified system. N is the total number of samples, M is the dimension of the outputs. After that, in the MSWOA-based optimizer, the parameters to be identified are identified by minimizing C OF (θ). As the optimization loop in the MSWOA-based optimizer goes on, the parameters to be identified accurately match the real values. The process of identification is illustrated in Figure 7, where θ is a parameter vector of which each element is a real value whileθ is a parameter vector of which each element is to be identified.

Generator System Model
The dynamic equation of the synchronous generator is described as Equation (23) in which m (s) means the load disturbance. When the HTGS is working under load conditions, m (s) ≠ 0 and e ≠ 0. When HTGS is working under no-load conditions, m (s) = 0 and e = 0.
where T is the inertial time constant of the generator, e is the adjustment coefficient of the generator and e = e − e .

Identification Strategy Based on MSWOA
The proposed MSWOA is employed to identify the HTGS parameters. In the MSWOA identification strategy, firstly, an input is employed to activate the real system and the system to be identified. The outputs of the real system and the system to be identified are input into Equation (24) to evaluate the value of C θ .
where θ is six-dimensional vector and θ = T T h T T e , z = x y m is the output of the real system, z = x y m is the output of the identified system. N is the total number of samples, M is the dimension of the outputs. After that, in the MSWOA-based optimizer, the parameters to be identified are identified by minimizing C θ . As the optimization loop in the MSWOA-based optimizer goes on, the parameters to be identified accurately match the real values. The process of identification is illustrated in Figure 7, where θ is a parameter vector of which each element is a real value while θ is a parameter vector of which each element is to be identified. PE (parameter error) and APE (the average parameter error) which are described by Equations (25) and (26) are employed to evaluate the accuracy of the parameters identified by MSWOA: where is the parameters of the real system, is the parameters of the identified system, m is the total number of parameters of while k is the k-th parameter of . PE (parameter error) and APE (the average parameter error) which are described by Equations (25) and (26) are employed to evaluate the accuracy of the parameters identified by MSWOA: where θ k is the parameters of the real system,θ k is the parameters of the identified system, m is the total number of parameters of θ while k is the k-th parameter of θ.

Experiments and Analysis of Parameters Identification
In this section, MSWOA is employed to identify the parameters of a mathematical model of HTGS which is simulated in MATLAB R2016a. T y1 , T y ,T r , h w , T a , e g are the parameters to be identified. In identification experiments, two working conditions of HTGS, namely the no-load condition and the load condition, are both taken into account. The step disturbance of frequency and load are utilized to excite the system, respectively. Under no-load conditions, the amplitude of the step disturbance of the frequency is 0.1 p.u and under no-load conditions the amplitude of the step disturbance of the load is 0.1 p.u, respectively. Other parameters are set as follows: the total time of the simulation experiments is 30 s and the sampling time is 0.01 s. K p = 5.59, K i = 1.06, K d = 3.3, T d = 0.28, b p = 0.04, the initial value of the parameter vector θ is θ = 0.1 0.3 1.5 0.5 12 0.5 , the transfer coefficients of HTGS are selected referring to [2]. The transfer coefficients are listed in Table 9. To prove that the MSWOA is superior to other algorithms in parameter identification, PSO, ALO, WOA, EWOA and CWOA were employed to identify the parameter of HTGS as a comparison. The identification experiments of any algorithm are independently repeated 20 times and the average values of the identified parameters of each algorithm have been obtained. The size of the population and the total number of iterations are set as 30 and 100 for each algorithm. The parameter setting of the algorithms proposed in this paper is described in details in Table 5.

Comparison of Different Identification Methods under No-Load Condition
In Figure 6, c and mg are the frequency disturbance and the load disturbance signals, respectively. Under no-load conditions, c is added to HTGS without loads, namely c, and mg = 0. Moreover, because the value of mg is zero, e g = 0. Therefore θ has five elements, namely T y1 , T y , h w , T r and T a , that need to be identified. In addition, because a linear model of the hydraulic turbine is used in the identification experiments, the transfer coefficients are given values.
Mean values of the identified parameters and mean best cost for 20 repetitions of PSO, ALO, WOA, EWOA, CWOA and MSWOA are listed in Tables 10 and 11, respectively. It is obvious that the mean best cost and mean APE of MSWOA are smaller than those of the others. Furthermore, the mean values of the identified parameters except h w of MSWOA are more approximate to the real values than the others. In short MSWOA can get a comprehensively better accuracy in parameter identification than other algorithms. In Figure 8, the ability of exploration and exploitation of WOA and EWOA are worse than those of the others. PSO, which is easily trapped in local optimum points, achieves a better ability of exploration and exploitation than those of WOA and EWOA. Though CWOA and ALO have a similar iteration process, ALO achieves a better exploration and exploitation ability than CWOA. MSWOA achieves a better ability of escaping from of local optimum points than others. Therefore, MSWOA has a better ability of exploration, exploitation and escaping from local optimum points than other algorithms.  In Figure 8, the ability of exploration and exploitation of WOA and EWOA are worse than those of the others. PSO, which is easily trapped in local optimum points, achieves a better ability of exploration and exploitation than those of WOA and EWOA. Though CWOA and ALO have a similar iteration process, ALO achieves a better exploration and exploitation ability than CWOA. MSWOA achieves a better ability of escaping from of local optimum points than others. Therefore, MSWOA has a better ability of exploration, exploitation and escaping from local optimum points than other algorithms. The outputs of HTGS with real values are compared with those of HTGS with the identified parameters in Figure 9. The compared outputs are the turbine speed, guide vane opening and turbine torque, respectively. The curves of the identified system are approximately the same as those of the real system. This verifies that MSWOA is effective for HTGS parameter identification with a delayed water hammer effect under no-load working conditions.
(a) The outputs of HTGS with real values are compared with those of HTGS with the identified parameters in Figure 9. The compared outputs are the turbine speed, guide vane opening and turbine torque, respectively. The curves of the identified system are approximately the same as those of the real system. This verifies that MSWOA is effective for HTGS parameter identification with a delayed water hammer effect under no-load working conditions.  In Figure 8, the ability of exploration and exploitation of WOA and EWOA are worse than those of the others. PSO, which is easily trapped in local optimum points, achieves a better ability of exploration and exploitation than those of WOA and EWOA. Though CWOA and ALO have a similar iteration process, ALO achieves a better exploration and exploitation ability than CWOA. MSWOA achieves a better ability of escaping from of local optimum points than others. Therefore, MSWOA has a better ability of exploration, exploitation and escaping from local optimum points than other algorithms. The outputs of HTGS with real values are compared with those of HTGS with the identified parameters in Figure 9. The compared outputs are the turbine speed, guide vane opening and turbine torque, respectively. The curves of the identified system are approximately the same as those of the real system. This verifies that MSWOA is effective for HTGS parameter identification with a delayed water hammer effect under no-load working conditions. (a)

Comparison of Different Identification Methods under Load Conditions
Under load condition, only a load disturbance is added to HTGS. Therefore e ≠ 0 and the parameter vector θ have six elements, namely T , T , h , T , T and e , that need to be identified.
In Table 12, it is obvious that the four parameters identified by MSWOA, namely h , T , T and e , more accurately match the real values than the other algorithms while T and T are not the most accurate match with the real values, but T and T are more accurately matched to the real values than those of WOA and EWOA. Furthermore, the mean best cost and mean APE achieved by MSWOA are both better than those of the other algorithms in Table 13.

Comparison of Different Identification Methods under Load Conditions
Under load condition, only a load disturbance is added to HTGS. Therefore e g = 0 and the parameter vector θ have six elements, namely T y1 , T y , h w , T r , T a and e g , that need to be identified.
In Table 12, it is obvious that the four parameters identified by MSWOA, namely h w , T r , T a and e g , more accurately match the real values than the other algorithms while T y and T y1 are not the most accurate match with the real values, but T y and T y1 are more accurately matched to the real values than those of WOA and EWOA. Furthermore, the mean best cost and mean APE achieved by MSWOA are both better than those of the other algorithms in Table 13.  In Figure 10, PSO has a better exploration and exploitation ability in the early stages of the iteration process but gets trapped into prematurity later. MSWOA has a better ability of exploration, exploitation and escaping from local optimum points than the other algorithms. In Figure 10, PSO has a better exploration and exploitation ability in the early stages of the iteration process but gets trapped into prematurity later. MSWOA has a better ability of exploration, exploitation and escaping from local optimum points than the other algorithms. The outputs including the turbine speed, guide vane opening and turbine torque of the real system of HTGS and the identified system acquired by MSWOA are compared in Figure 11. The curves of the identified system almost overlap with those of the real system. Water hammer effects would cause a change of flow parameters which have to be managed by the guide vane system, so the fact that the outputs of the real system match those of identified system shows that MSWOA is effective for HTGS identification parameter with a delayed water hammer effect under load working conditions. Therefore, we can draw the conclusion that MSWOA has better global optimization ability and can get better accuracy in parameter identification than the other algorithms proposed in this paper whether under no-load conditions or load conditions. (a) The outputs including the turbine speed, guide vane opening and turbine torque of the real system of HTGS and the identified system acquired by MSWOA are compared in Figure 11. The curves of the identified system almost overlap with those of the real system. Water hammer effects would cause a change of flow parameters which have to be managed by the guide vane system, so the fact that the outputs of the real system match those of identified system shows that MSWOA is effective for HTGS identification parameter with a delayed water hammer effect under load working conditions. Therefore, we can draw the conclusion that MSWOA has better global optimization ability and can get better accuracy in parameter identification than the other algorithms proposed in this paper whether under no-load conditions or load conditions. In Figure 10, PSO has a better exploration and exploitation ability in the early stages of the iteration process but gets trapped into prematurity later. MSWOA has a better ability of exploration, exploitation and escaping from local optimum points than the other algorithms. The outputs including the turbine speed, guide vane opening and turbine torque of the real system of HTGS and the identified system acquired by MSWOA are compared in Figure 11. The curves of the identified system almost overlap with those of the real system. Water hammer effects would cause a change of flow parameters which have to be managed by the guide vane system, so the fact that the outputs of the real system match those of identified system shows that MSWOA is effective for HTGS identification parameter with a delayed water hammer effect under load working conditions. Therefore, we can draw the conclusion that MSWOA has better global optimization ability and can get better accuracy in parameter identification than the other algorithms proposed in this paper whether under no-load conditions or load conditions. (a)

Conclusions
MSWOA, a novel algorithm-based WOA with a mixed strategy, is proposed in this paper. Compared with the standard WOA, three improvements have been made to enhance the searching ability. Firstly, because the strategies of bubble-net attacking and encircling prey can identify the optimization points in a local region, a hybrid movement strategy is applied on MSWOA, in which a dynamic ratio based on improved JAYA algorithm is applied on the strategy of searching for prey and a chaotic dynamic weight is applied on the strategies of bubble-net attacking and encircling prey. Secondly, a guidance of the elite's memory inspired by PSO is applied on the movement of whale agents of the population. Thirdly, the mutation strategy based on the sinusoidal chaotic map is employed to avoid prematurity and local optimum points. Subsequently the MSWOA is compared with six meta-heuristic algorithms on 23 benchmark functions and the results show that MSWOA achieves remarkably better performance than the others, and the significance has been confirmed by box and whisker's tests. Finally, the proposed MSWOA, together with ALO, PSO, WOA and EWOA and CWOA, are applied in parameter identification of a HTGS with a delayed water hammer effect.
The results reveal that the MSWOA demonstrates satisfactory global searching ability and dramatically promotes the identification accuracy of the complicated system studied, compared with other existing algorithms.
Author Contributions: T.D. and C.L. conceived and designed the experiments; L.C., C.F., N.Z. conceived the study; T.D. proposed the algorithm and wrote the paper; T.D., C.L. and L.C. played important roles in the process of revising the paper.

Conclusions
MSWOA, a novel algorithm-based WOA with a mixed strategy, is proposed in this paper. Compared with the standard WOA, three improvements have been made to enhance the searching ability. Firstly, because the strategies of bubble-net attacking and encircling prey can identify the optimization points in a local region, a hybrid movement strategy is applied on MSWOA, in which a dynamic ratio based on improved JAYA algorithm is applied on the strategy of searching for prey and a chaotic dynamic weight is applied on the strategies of bubble-net attacking and encircling prey. Secondly, a guidance of the elite's memory inspired by PSO is applied on the movement of whale agents of the population. Thirdly, the mutation strategy based on the sinusoidal chaotic map is employed to avoid prematurity and local optimum points. Subsequently the MSWOA is compared with six meta-heuristic algorithms on 23 benchmark functions and the results show that MSWOA achieves remarkably better performance than the others, and the significance has been confirmed by box and whisker's tests. Finally, the proposed MSWOA, together with ALO, PSO, WOA and EWOA and CWOA, are applied in parameter identification of a HTGS with a delayed water hammer effect.
The results reveal that the MSWOA demonstrates satisfactory global searching ability and dramatically promotes the identification accuracy of the complicated system studied, compared with other existing algorithms.
Author Contributions: T.D. and C.L. conceived and designed the experiments; L.C., C.F., N.Z. conceived the study; T.D. proposed the algorithm and wrote the paper; T.D., C.L. and L.C. played important roles in the process of revising the paper.