Optimal Controller Design for Speed Governors of Hydroelectric Power Plant

Speed governors have critical importance on hydroelectric power plants, which are adjusted to the rotating speed of hydroelectric generation based on load demand of the grid. The rotating speed is the main factor to balance power generation and load demand. The well-designed controller is needed to control speed governors with high accuracy. A well-defined model is needed to obtain desired control structure. Therefore, in this study, initially, the mathematical model of a hydroelectric power plant is obtained by using physical characteristics of a real-world. Then by using this model and corresponding real-world data, a set of controller parameters is designed by using tuning methodologies based on heuristic optimization algorithms, and their performances are compared with each other and with a classical tuning methodology. Evolutionarybased and nature-inspired-based heuristic optimization algorithms are selected as the tuning algorithms not only to compare the performance of these algorithms with a classical method but also with different origins. The performance of the optimized controller improves the performance of the overall system and helps to get desired performance. The results also indicate that as long as the desired performance criteria are defined as accurate as possible, the performance of the optimization algorithms is acceptable.


Introduction
A hydroelectric power plant is composed of ball valve, shunt valve, induction manifold, volute, adjusting fan, gear, speed regulator, auxiliary systems (break, cooler, fan, greasing etc.), generator, voltage regulator, warning system, synchronization system, protecting system and switch-yard equipment (breaker, dispersive, transformer, etc.). The operation of hydropower plant can be summarized as: at first, the initialization signal is reached to plant, and then the auxiliary systems are engaged. A typical auxiliary system is formed from cooler liquid pump, speed regulator oil pump, bearing lubrication pump, rubber sleeve valve and carbon dust discharge valve. After that, mechanical breaks are deactivated, and emergency valves become online. Then speed regulator is initialized, and wicket gate slightly opens. After some time, turbine reaches its nominal speed. At this phase, it is desired to be synchronized with the grid. Hence the other conditions for synchronization are fulfilled following the nominal speed. For synchronization with the grid, voltage-level of regulation is reached to its nominal value. The synchronization is assumed to be succeeded when the generator voltage level and frequency meet with grid voltage and frequency levels respectively; and then unit circuit breaker closes. Finally; load is applied to the system after synchronization by opening speed regulator wicket gates. It is clear that the control of speed governors has a direct effect on both generation and the synchronization with the grid. In a general manner, the aim of the speed governors is to set the nominal value by using the information from measurement equipment in the field. The controller is utilized for setting the system to make it reach the nominal value. Speed governor control functions can be summarized as: • Speed Control Mode (Transient Response): The speed of the governor is controlled until the synchronization with the grid succeed. Speed regulator slightly opens the wicket gate, by this way turbine begins to rotate and speed-up. Next, wicket gate moves to a lower position and unwanted increase at the water input is prevented. This control action goes on until the synchronization with the grid is achieved.
• Active Power Control Mode (Steady-State Response): After the speed control mechanism is terminated, in other words, the synchronization between plant and grid is succeeded, the active power control mode runs to support the balance the grid frequency by primary frequency control.
In this study, the transient response of the controlled speed governors at hydroelectric power plant is investigated since the transient response of the generator has the most critical mission and it is expected to give essential response to the imbalances of the grid as fast and as accurate as possible. In another respect, it is desired to achieve synchronization between power plant and grid rapidly. In other words, the error between turbine speed, and the frequency of the grid must be zero (the steady-state error must be zero). For all of these reasons, the controller design has an importance at the hydroelectric power plants. Therefore, a classical controller algorithm is applied to the problem. The parameters of the controller have direct effect on the transient response characteristics of the system. In this study, parameter tuning methodologies with respect to the optimization algorithms are investigated by using the mathematical model of the power plant, which is obtained from the field study on a real-world plant.

Literature Review
In the approaches to the design of controllers for power plants, common property is noticed, namely that the controllers are designed based on the various sections of the different types of power plants. Therefore, many approaches are proposed for each power plant. Furukowa et al. investigated a blackout scenario where hydroelectric plants are initially started up after blackout to support the grid [1]. Two control structures are applied to the system; conventional Proportional (P) control, and lag compensator. The authors showed that if the same control parameters are applied to the proposed control system, the overall system becomes unstable. To improve the performance of the controller, the authors suggested adding an Integrator (I) part to the proportional controller instead of P controller and PI controller is applied to the problem [1]. It is also indicated in the paper that speed regulator stability is the important aspect for a start-up hydroelectric power plant. Rastrepo and Galiana investigated reducing the total production cost in power plants (a similar study is applied for wind turbines [2]) [3]. That study showed that speed control action has a direct impact on the production costs, and a welldesigned controller help to reduce it. In the study by Nanda et al., ideal I controller and PI controller are compared on thermal power plants, and it is concluded that almost the same performance is obtained from both algorithms [4]. Abdolmaleki et al. suggest an adaptive PI controller based on fuzzy logic for hydropower plants [5]. The results indicated that the transient response becomes faster, and causes large stability margin. Also, the robustness increases with the proposed algorithm. The study showed that conventional PI controller with a proper parameter selection greatly helps to improve the overall performance of the system. These results are also supported by a different study from the same authors [6]. In that study, the best-possible controller parameters are obtained from classical control design approaches (root locus), and the results showed that performance of the system is improved with a properly adjusted set of optimal parameters. Also, from all these researches, it is learned that almost same performance can be obtained from another set of parameters and controllers.
As another methodology, intelligent controllers, such as neural network controller, are applied to the problem. A neural network controller is applied to this problem. Ram and Jha suggested a controller framework where PI controller is supported by a soft computing approaches -Fuzzy Logic, Neural Network and Neuro-Fuzzy System [7]. The performance of these joint systems is almost the same with the rise time and settling time reduced. Nanda et al. [8] discussed the performance of I, PI, PID and ID controllers on hydrothermal power plant. Even these controllers present almost the same performance, Differentiation (D) term increases the noise of the controller. Khuntia ve Panda [9] are discussed the Neuro-Fuzzy System (ANFIS) and compared with ideal I compensator. The results also suggest that an intelligent controller gives more robust and fast response against the compensator. Laghardi et al. present a general discussion about conventional control algorithms (PID, lag and lead compensators) [10]. It is indicated that PID and PI controllers give the almost same overshoot performances, and the I term occurs as the important term for reducing (elimination) steady-state error. Jagatheesan and Anand [11] present a paper and compare I and PI controllers. The performance of the overall system indicates that with I (ideal integral compensator) controller the overshoot becomes larger than in case of PI controller.
As a joint study that merges intelligent optimization algorithm and conventional controllers, Singh et al. give an optimal parameter selection/determination for PI controller by using the Genetic Algorithm (GA) [12]. The result of this paper also supports the previous discussions that the performance of the controller significantly depends on the appropriate parameter selection. The idea of tuning the parameters of PID (or other algorithms) is an efficient methodology to get desired performance. Currently, even the manyobjective optimization algorithms can be applied to the tuning process [13], single objective optimization algorithms are still successfully applied to the tuning problem [14], where overshoot of the output of the inverted pendulum system decreases and better stabilizes [14]. As a power system application for nuclear reactors, GA and PSO is applied for PID controller [15]. The performance index of ITSE is selected as cost function, and the results support PSO against GA algorithm. The time varying property of the PID controller tuned by PSO is investigated in [16]. Also, in that study, CMA-ES algorithm is also applied to the benchmark problems of optimization algorithms. As an application of the tuning methodology, 10 kW air heater system is evaluated as the problem [17]. The results of this study, which is based on comparing PSO-PID tuning with the classical tuning algorithm like ZN, shows the heuristic-based tuning algorithms performance. Among these applications MIMO systems were also investigated with a PID variant named as 2-DOF PID controller and PSO is applied and compared with other heuristics [18]. Also, as a cost function minimizing Integral of Squared Error (ISE) and Integral of Time multiplied by the Squared Error (ITSE) are selected. It is possible to implement these controllers tuned by heuristic algorithms as the engineering problem [19]. The coding is one of the fundamental problems of convolution optimization algorithms. The correct and intelligent coding will increase the optimization performance. Therefore, a new coding is also proposed, and impact of the coding is investigated [20]. Other optimization algorithms like Firefly [21], Chaotic Firefly Algorithm [22], self-adaptive firefly algorithm [23], gravitational search algorithm [24], and Bat-Inspired Algorithm [25] are applied to controller tuning problem. The performance of the tuning algorithm is demonstrated on all these cases.
By evaluating the previous studies explained above, the following statements which are the direction of this study are extracted: • Integrator is an important control parameter term which eliminates steady-state error.
• Derivative term increases the noise in the system and may cause instability.
• PI gives lower overshoot and faster response when compared to ideal integral compensator.
• Intelligent methods improve the performance of the controller with respect to the steady-state and transient responses.
• The optimal parameter selection not only improves the performance but also increases robustness of the system. Therefore, in this study, conventional PI controller parameters are determined by using heuristic optimization algorithms named as Differential Evolution, Particle Swarm Optimization, Firefly Algorithm and Evolutionary Strategy. Then, the optimized controller parameters are applied to the system model, which is extracted using the real-world hydropower plant. This paper is organized into five sections. The sub-section of introduction gives the survey of the similar papers in the literature. Also, in this sub-section with the aid of previous papers, a direction of this paper is emphasized. Section 2. gives the mathematical model of the power plants. Section 3. explains the optimization algorithms. Section 4.
presents the results of the implementations, and the last section is the conclusion of this study.

Problem Definition and Mathematical Modelling
Mathematical model of hydropower plant is obtained by extraction of mathematical relations inside hydropower plant [26]. The relation between tunnel dynamics is given (Fig. 1).
where the following parameters are used: U tu -speed of water in tunnel (pu); H r -level of reservoir (pu); H d -level of surge tank (pu); H l2 -lost in tunnel (pu); t t -time to reach water to tunnel (s). The relation with respect to the surge tank is given below: where C d is the constant (s) and U d is water speed at surge tank (pu). To add penstock dynamics to the model the Eq. (3), Eq. (4) and Eq. (5) are used.  where the variables are: H tr -turbine level (pu); H l -loss in penstock (pu); z p -stability impedance in penstock; t c -penstock time (s); U tr -speed of water to penstock (pu); G, wicket gate opening (pu). First, the no load condition is extracted for calculation of necessary power to move turbine. The Equation (6), Eq. (7) and Eq. (8) are used for this purpose.
where the used parameters mean: A tr -turbine gain (pu); U N L -no load flow (pu); U F L -full load flow (pu); w m -angular velocity (rad·s −1 ); T mek -torque (Nm); P mek -power (W). To find the acceleration torque, the difference between mechanical and electromechanical torque is calculated as: where T elek is the electromagnetic torque (Nm) and J is the joint moment between generator and turbine (kg · m 2 ): where H is a constant (kWs/kVA); w m0 -nominal angular velocity (rad·s −1 ); V A -nominal power (VA); where L c is the length of penstock (m), and α is wave speed (m · s −1 ). Hydraulic impedance is calculated from: where t c2 -time for water to reach penstock (s); A c -penstock cross sectional area (m 2 ); Q -nominal flow rate (m 3 ); H -nominal head (m). Similarly, the time for water to reach the tunnel is calculated from: where L t -length of tunnel (m), A t -tunnel cross section (m 2 ) and C d -penstock constant, which is calcu- . In this study, by using the real hydropower plant (Kadincik -I power plant) parameters and the nonlinear equations, the control algorithms are designed and compared with each other. The next section explains the control structure and algorithms used in this paper.

The Controller Structure
The PID controller is one of the basic but frequently preferred control algorithms due to its structure, and already discussed nature and necessary tools needed for tuning the parameters of the controller. The general frequency domain representation of the PID controller is given below: The PID controller in Eq. (15) has three terms and each of these terms has an impact on transient response. The proportional term (K P ) decreases rise time and steady-state error but increases overshoot. The integral term (K I ), gives similar performances like proportional term but increases settling time. But with the aid of integral term, steady-state error can be eliminated. The derivative term (K D ) has no positive effect on steady-state error and slightly changes rise time. However, this parameter decreases overshoot and settling time with a proper parameter value. The study [27] showed that like conventional tuning algorithms, heuristic algorithms are applied successfully to the tuning problem. In contrary to the results of conventional methods, heuristics gives much better results for systems of various order and even systems with time delay [28].

Performance Metric
The tuning of PID parameters by using optimization algorithms is an objective function needed to obtain a better transient response. Therefore, in this study, Mean Square Error (MSE) is selected as the main component of the objective function at the entire optimization algorithm. The previous study [29] showed that MSE could help to improve the steady-state error without critical change of the other transient response characteristics. Also, the difference between reference input and output of the system is selected as the objective function parameter with the steady-state error.

Differential Evolution (DE)
Differential Evolution (DE) algorithm is a populationbased evolutionary algorithm like genetic algorithm, which is composed of the mutation, crossover, and selection operators. The first introduction of the algorithm was proposed by Storn and Price in 1995 [30], [31] and [32]. In their first paper [30], the performance of their new algorithm by comparing with improved Simulated Annealing (SA) algorithm and Genetic Algorithm (GA) is presented in detail. (This paper is the reason for that SA and GA are not included as an optimization algorithm since it is shown that DE outperforms SA and GA). The organization of DE is very similar to GA. The algorithm begins with the random initialization of the solution candidates on the search space. In a general manner, the uniformly distributed random number generator is preferred. Instead of the last operator in GA, the mutation operator is evaluated as the first operator to obtain a new population set. The equation given below mathematically describes the mutation operator.
where x -member of the initial population (X), v -member of the produced population (V ), F -control parameter selected as 0.8, and r1, r2, and r3 -randomly selected indexes inside X.
As the second operator, the crossover is applied to the population V and population U is formed. Two approaches, proposed for crossover, are called as exponential and binomial operator. Since the performances of these approaches are almost the same [31], generally, binomial operator is preferred due to the low complexity. The binomial crossover operator is given below: where C R -second control operator of DE and the value is selected as 0.8 [30] and [32]. As the final operator; the best members among X and U are selected and survived to the next generation.

Particle Swarm Optimization (PSO)
Particle Swarm Optimization (PSO) algorithm is a nature-inspired heuristic population-based algorithm proposed by Kennedy and Eberhart in 1995 [34]. The algorithm is formed by investigating the behaviours of the members in an animal swarm [34], [35], [36] and [37], mainly their motion. It is observed from the animal swarm, that the members are deciding their movements by following the swarm leader and their instinct to return their most benefit position. Therefore, the best member at every iteration of the algorithm is determined with the personal best position. In summary, the algorithm begins with the randomly distribution of the members on search space. Then for each member, the objective value is calculated. By using these objective values, the member with smaller objective value (for minimization problem) is found. Also, from these objective values, each member updates its personal best position. By using all these calculated val-ues, the position and velocity are updated as given by equations: where x k (d) -position of the particle in the d th dimension at the k th iteration, similarly v -velocity component, p best -best position among the memory of the member, g best -best position of the whole swarm, w -inertia weight which is initialized as 0.9 and decreased linearly to 0.4, c 1 and c 2 -control parameters and usually chosen to be 1.494, rand -pseudo-random number generated from a uniform distribution, and ∆t -time step = 1.

Firefly Algorithm (FA)
Firefly Algorithm (FA) is a population-based optimization algorithm, which mimics the behaviour of the firefly flocks in their nest. The algorithm proposed by Yang [38] and [39] in 2008 is based on the light intensity and attractiveness between fireflies. Fireflies use their flashing light properties as a warning mechanism against predators and to attract for mating. The amount, rate and time duration of the light is a communication mechanism among the fireflies. The light intensity is changing with respect to the absorbing the light source and it is related to the distance between fireflies. Therefore, by using the light intensity and corresponding attractiveness between fireflies, an algorithm called as Firefly Algorithm is proposed for solving optimization problems. To form a relation between optimization algorithm and firefly behaviours, it is assumed that a) All the members can be attracted from each other. In other worlds all the members are unisex. b) Attractiveness is only depended on the light intensity. Therefore, less bright member moves to brighter one. If there is no brighter member due to the distance, then it moves randomly. c) The search space is corresponding to the light distribution. Therefore, brightness corresponds to the objective value. The performance of the FA among other heuristic algorithm on same well-known benchmark problems are discussed in detail [40], and results showed that for some of the benchmark/test problems, FA gives same or better performances against well-known optimization algorithms. Then, the improved version of this algorithm also discussed and converted to solve multimodal problems [41]. Next, since it can present good performances, after some time for obtain mature FA, it is applied to some real-world problems [39] and [42] and produced satisfactory performances. The members of the populations (firefly) assumed to have a position (x) at the same dimension of the search space. Based on the position and distance, the light intensity is changed. For a uniform medium, which has a constant light absorption coefficient (γ), initial light intensity (I 0 ) at the distance (r), and light intensity can be calculated as given: Similarly to the light intensity, as the attractiveness is directly proportional to the light intensity, the attractiveness is calculated from the following equation.
The attractiveness is mutual for both fireflies. Therefore, the distance between them is a parameter at the formulation. The distance is calculated from Euclidean function (Eq. (22)).
In the algorithm, the position of the members corresponds to the solution of the optimization problem. Therefore, the position update rule is the main mechanism of the algorithm. The following Eq. (23) is the position update rule of the FA algorithm. The update mechanism has three parts connected with addition. The first term of the equation is the previous position, the second term is the attractiveness, and the last term is the randomness of movement, since the distance or light intensity between members is the same and the second term is zero.

Covariance Matrix Adaptation Evolution Strategy (CMA-ES)
Initially, Evolution Strategy (ES), as the name expressed, is designed as a set of rules based on the evolution process in nature for using at automatic design of the engineering problems [43]. Then, the single solution-based ES has only one parent and one offspring per generation. This new form of ES is called as (1+1)-ES. Next, in the master thesis of Schwefer [44], it is found that the algorithm is stuck on local optimum frequently. Finally, multi-member (populationbased) ES is proposed and called as (µ + 1)-ES, and then (µ+λ)-ES is proposed. In this algorithm weighted average of λ is selected from µ offspring which is generated at each generation. After years, the algorithm becomes a population-based evolutionary algorithm with the operators of marriage, recombination and mutation. Since the local optimum problem is the reason of the improvements on ES algorithm, in 2001, Hansen proposed an adaptive covariance approach to ES algorithm [45], and this CMA-ES algorithm is evaluated [46] and discussed in the literature in detail [47] and [48].
The CMA-ES has structure with five stages: • initialization, • update of covariance and step size, • generation of offspring, • recombination.
Unlike other optimization algorithms discussed in this paper, for this population-based method, the solution candidates are selected and updated randomly by using multivariate Gaussian distribution with Different Mean (M ), Variance (σ) and Covariance (C) values.
where g -index of generation, k -index of members in population. Like the other optimization algorithms, x corresponds to the vector of elements and d to the dimension of the problem. At the first step, the population is randomly distributed on the search space. Then, two fundamental matrices are updated from the Eq. (25), Eq. (26), Eq. (27) and Eq. (28). These vectors are covariance matrix and step size, where n -number of decision variables, µ -best individual of the population, c cov = [0, 1] -control parameter (determines the rate of change of the covariance matrix) for step size and covariance matrix [49].
Similarly, step size σ -related to the parameter called conjugate evolution path P After all these steps are completed, the algorithm is repeated until the termination conditions are met.

Implementation Procedure and Results
In the study, the aim is to present a better control parameter selection methodology for hydropower plants. Therefore, as an initial step, the control action is selected. Based on the previously published papers (discussed in Sec. 1.1. ), it is observed that the PI controller gives the sufficient/better performance with respect to the desired criteria. The reason is that the integral term is applied to reduce/cancel the steady-state error, even to make it zero. If a Differentiation term (D) is added to the controller, then a noisy output signal is obtained and none of the papers for hydropower plants from our research is suggested to use D term. Also, from the previous studies, it is possible to make an inference that the soft computing approaches help to increase the performance of the classical controller. Similarly, when a good controller parameter is selected, not only the performance of the controller increases, but also robustness and overshoot improve. Therefore, in this study, the optimal parameters of PI controller are obtained using heuristic optimization algorithm. In this study, DE, PSO, FA and CMA-ES algorithms are applied to the problem. These algorithms can be categorized into two groups. One includes the natureinspired optimization algorithms like PSO and FA; the other group is for the evolutionary algorithms like DE and CMA-ES. The reason is to compare the performance of the algorithms of different categories. Therefore, two relatively new and relatively widely preferred algorithms are selected due to their acceptable performance.
To show the performance of these tuning algorithms, a classical tuning method; Ziegler-Nichols is initially applied to the problem. Ziegler-Nichols (ZN) tuning method is based on determination of the critical value at marginal stable condition. In other words first it is assumed that the system has only proportional part and it increases until the closed-loop poles are almost at the jw-axis, which means the almost marginal stability of the system. The ZN tuning algorithm is based on determination of the PI parameter with respect to the critical gain and critical period. The proportional value for marginal stability is named as critical gain K cr , and the period of oscillation is called critical period T cr . The PI controller parameters are selected as K P = 0.45K cr and T i = T cr · 1.2 −1 with respect to the critical gain and critical period. Therefore, initially critical gain and critical period are determined such that the output exhibits an oscillation (marginal stability), and the period of the oscillation is the critical period.
In order to determine the parameters of the control method with heuristic optimization algorithms, simula-tion studies are carried out based on the predetermined model. The purpose function in simulation studies is defined by the following equation.
where J -cost function of the problem, t b and t f -beginning and final time of the simulation, K p -position constant (steady state error), and e(t) -error function between reference input (desired speed) and the output of the system. The minimum value of the variable 'J' defined in the equation expected from the optimization algorithm is found. As different goal values can be selected here, in this study, only the smallest value of the error signal given by 'e' was requested (Fig. 3).
The reason for choosing this way is not that the desired system reaches the desired value quickly, unlike other mechanical and electrical systems. Instead, it is neither too slow nor too fast and minimizes the steadystate error. The examined system consists of many sub-systems in connection with each other and these systems must be put into operation in order to work together. For this reason, the system must first become a steady state. After this stage, commissioning should be done slowly in order not to put too much load on other mechanical systems. The purpose function determined in the equation is therefore chosen.
Initially, to compare the performance of the heuristic algorithms, the system is tuned by using ZN, where K cr = 1.56 and T cr = 7.2. The controller parameters are calculated as K P = 0.7 and T i = 6. Note that, instead of T i , K I can be used as controller parameter which is equal to K I = K P · T −1 i = 0.116. In the second part of the implementation, heuristic optimization algorithms are applied as a tuning strategy to the problem. For this purpose, all these algorithms have executed in the same conditions such as population size 60, number of iteration/generations 100. The solution candidates are selected inside 2]. Controller parameters obtained after the completed iteration are reported in Tab. 2, and their transient responses are presented graphically in Fig. 2   of these controllers have the integral term, the steadystate errors are converged to approximately zero. In case of overshoot values, PSO and CMA-ES give the worst performances with respect to the overshoot, such that almost 9 % overshoot is observed from the transient response of the CMA-ES parameter system. At the same time, CMA-ES gives the oscillation output at the stability border, which is expected to reach steady state at time approaching infinity. For peak time values of ZN, DE and FA are larger than settling time, which means that there is no overshot and systems reach final value. Therefore, for these three systems, the settling time is the important indicator of the speed of the systems. Therefore, as seen in Tab. 1, DE gives the best performance among all of the tuning algorithms. Then, ZN gives the smoothest and fast performance, and finally, FA gives the third best performance among all tuning algorithms.

Conclusion
The hydropower plant composes of different mechanical and electrical systems where all of these subsystems run to produce electricity. The control of speed governors has an important duty for synchronization with the grid. Therefore, a well-designed controller is needed for this duty. In this study, initially a model of hydropower plant is constructed with the aid of realworld power plant parameters. Then, PI controller is selected as the main control algorithm and the parameters are optimized using ZN, DE, PSO, FA and CMA-ES algorithms. The results showed that heuristic optimization methodology improves the performance with respect to the settling time. However, for this relatively complex problem, it is not obtained best performance from all of the heuristic tuning algorithms. The CMA-ES gives almost unstable performance, PSO gives large overshoot but relatively fast response, but even the overshoot is not a good decision for synchronization, FA presents good but slow response, ZN gives acceptable average performance and DE gives the best performance among all the tuning methodologies.