Optimal neuro-fuzzy model and PID controller of a hydro-generator unit from the identification of its LTI model

ORIGINAL ARTICLE Optimal Neuro-fuzzy model and PID controller of a Hydro-generator unit from the identification of its LTI model Jesu ´s Rodrı´guez-Flores a , Victor Herrera-Perez b,* , Mayra Pacheco-Cunduri c , Jorge Herna ´ndez-Ambato c , Alejandro Paredes-Camacho d , Miguel Delgado-Prieto d a Universidad Regional Auto ´noma de Los Andes (UNIANDES), Ambato 180101, Ecuador b Universidad San Francisco de Quito – USFQ, Institute for Energy and Materials, Quito, 170901, Ecuador c Escuela Superior Polite ´cnica de Chimborazo (ESPOCH), Facultad de Informa ´tica y Electro ´nica, Riobamba, 060101, Ecuador d Universitat Polite `cnica de Catalunya (UPC), Electronic Engineering Department, MCIA Research Group, Terrassa, 08222, Spain Received 25 October 2022; revised 23 February 2023; accepted 14 March 2023 Available online 28 March 2023 KEYWORDS LTI model; Neuro-fuzzy control; PID control; PI + D control; Hydro-generator unit Abstract This paper describes a methodology to obtain a neuro-fuzzy model of a hydro-generator unit (HGU) and its PID controller from the identification of its linear time-invariant (LTI) model. The study performs a cause-effect record, which allows a continuous time identification of the LTI type, then obtaining a classic PID controller capable of complying with a performance specification given by a pole of a second order system, which enables the training of the proposed neuro-fuzzy model. Applying the cost function of the root mean square error and the root of the percentage relative mean square error, the system parameters were adjusted using the decreasing gradient method. By means of linear models, the initialization of the singletons of the neuro-fuzzy models was done in two stages using the decreasing gradient and a cost function. The first stage was carried out without dynamics and the second stage with the dynamics of the simulated system. A case study of the Hydro-Agoya ´n HGU was selected and the results showed that the development of the LTI model allowed the development of the neuro-fuzzy model able to represent the behavior of the power plant and its response under variations of the power setpoint of 5, 10, 15 and 20%. Ó2023 THE AUTHORS. Published by Elsevier BV on behalf of Faculty of Engineering, Alexandria University. This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/ licenses/by-nc-nd/4.0/). 1. Introduction Hydroelectric energy is one of the main sources of energy in the world given the availability of water resources around the planet. It is generated, transmitted, and distributed in the form of electricity in order to satisfy all kinds of power demands. Hydroelectric energy provides certain advantages: *Corresponding author. E-mail addresses: [email protected]du.ec (J. Rodrı´guez-Flores), [email protected] (V. Herrera-Perez), mayra.pacheco@espoch. edu.ec (M. Pacheco-Cunduri), jorge.hernand[email protected] (J. Herna ´ndez-Ambato), alejandro.pa[email protected] (A. Paredes-Cama- cho), [email protected]u (M. Delgado-Prieto). Peer review under responsibility of Faculty of Engineering, Alexandria University. Alexandria Engineering Journal (2023) 71, 309–337 HOSTED BY Alexandria University Alexandria Engineering Journal www.elsevier.com/locate/aej www.sciencedirect.com https://doi.org/10.1016/j.aej.2023.03.039 1110-0168 Ó2023 THE AUTHORS. Published by Elsevier BV on behalf of Faculty of Engineering, Alexandria University. This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/licenses/by-nc-nd/4.0/). large load capacity, minimal fluctuation of electric power, and reduced environmental impact compared to fossil fuels [1].In fact, the importance of a constant supply of energy in all socio-economic sectors cannot be exaggerated, since it is directly linked to social stability and development [2]. In this regard, the living standards of people can be analyzed based on their per capita energy consumption. In 2000, the global electricity net consumption was 13,277 billion kilowatt-hours (kWh), a number that increased to 22,347 billion kWh in 2017 [3], and to 25,300 billion kWh in 2021, in which hydropower accounts a share of 17 % approximately of total electricity production [4]. The aim of each hydropower plant is to supply the necessary power to meet the energy demand following economic, sustainable and quality metrics [5]. Apart from the fact that the electrical energy demand is growing, as the demand follows non-stationary patterns a control system is required to automatically adjust the power generated in each generating unit, thus, trying to keep the operation parameters of the unit within certain limits [6]. In general, in an electrical power system, relevant parameters of power utility operation are variables such as frequency, voltage, turning speed and power, among others. In this regard, the frequency is one of the main variables, which should remain practically constant during the operation. An adequate control of the frequency ensures, for example, a constant speed of related synchronous/induction machines, which is essential for the correct operation of a generating unit since the auxiliary equipment used for controlling fuel, supplying water and combustion are based on this machinery. Nomenclature r2 yQuadratic Multivariable Factor i;jIndexers NNumber of samples ySamples bySample Prediction y Mean value of samples Recm Root Mean Square Error Recr%Root Squared Relative Percent Error sLaplace operator GsðÞ Transfer function aiCoefficients of the polynomial in sof the numerator of the transfer function biCoefficients of the polynomial in sof the denominator of the transfer function @Partial Differential Operator L1Inverse Laplace transform dsConvolution differential JCost function hParameter identifier of a system Cii-th singleton of the i-th fuzzy membership function lifuzzy i-th membership function bfxðÞ Fuzzy prediction of the function ffor the variable x aLearning factor biFuzzy i-th membership function Wij Weight or ij-th Singleton of the fuzzy ij-th rule DDiscrete differential Ksm Servo motor gain Tsm Servo motor time constant TWTurbine time constant (water hammer timer) KgGenerator efficiency TgGenerator time constant eError uInput yOutput xii-th state of the system YPID PID controller output d WG Prediction of the position of the Servomotor or Guide Vanes (Wicket Gate) min Minimum value max Maximum value b PmPrediction of Mechanical Power WG Position of the Servomotor or Guide Vanes (Wick- et Gate) b PeElectric Power Prediction MF Membership Function FIS Fuzzy Inference System PmMechanical power tupwg Rise time of the servomotor trRise time xdSystem oscillation frequency nConstante de amortiguamiento del sistema mp Maximum System Override rInverse of system response time SnPole assigned to the system xnNatural frequency of the system unPhase condition of the Pole assigned to the system ksmi Servo motor magnitude condition usm Servo motor phase condition kTTurbine Magnitude Condition uTTurbine phase condition kGGenerator Magnitude Condition uGGenerator phase condition hController operating frequency k2 dQuadratic Condition of Frequency of operation of the controller kcController parameterization response time KpController Proportional Gain TiController Integral Time TdController Derivative Time Gsm sðÞ Servo Motor Transfer Function GTsðÞ Turbine Transfer function GGsðÞ Generator Transfer Function b Gsm sðÞ Servo Motor Transfer Function Prediction b GTsðÞ Turbine Transfer Function Prediction b GGsðÞ Generator Transfer Function 310 J. Rodrı´guez-Flores et al. However, related control strategies are still inadequate facing undesired scenarios when turbines and generators problems arise, mainly due to high-demanding operation periods, parts wear and tear, or even improper usage by the unit’s operator managing the governor’s load limit device. These problems force the control area to maintain constant monitoring of the behavior of the units to take corrective actions that improve their performance [7]. It is in this regard when enhanced robust control methodologies, such as the one proposed in this study based on a Neuro-fuzzy and Proportional Integral Derivative Control (PID) model, are required. A reference review of the state of the art related to hydraulic turbine control can be found in [8], characterized by digital governors, mostly employing Proportional-Integral (PI) controllers and a few modern control techniques. On the one hand, referring to the classic PI control, stability issues and tuning have been studied mostly. In this regard, some reference works are the one presented by S. Hagihara et al.,related to the effect analysis of governor parameters on the stability boundaries of a hydraulic generation turbine unit through the use of the root locus method [9], or the one presented by D. T. Phi et al., in which an analysis to interpret stability limits as a function of speed regulations, turbine loading and system damping, was carried out [10]. Later, other significant works were published, such as [11], where an extensive and complete overview about the operation of governed hydraulic turbines is presented, [12] where a predictive approach was presented to proportional and integral gains tuning, or [13], where speed dynamics neutralization in control strategies is discussed in order to achieve improved power tracking capabilities. Aligned with this research line, one of the most outstanding contributions was made by the IEEE Working Group [14], establishing the necessary dynamic models to study failures conditions due the electrical system or the hydro-generator unit. The study revealed that the turbine firstorder model is feasible to regulation purposes when the penstock is relatively short, otherwise, in front of relatively long pen-stock or power oscillations in the hydrogenation process, second or higher order models (using Pade ´approximants) are necessary. Vournas et al.[15], and Choo et al.[16–17], incorporated the computational tools to ratify and improve some aspects of such aforementioned studies. On the other hand, referring to modern controls, several studies have covered nonlinear, optimal, and robust controls using state variable representations. In this regard, a fuzzy model control, for a nonlinear model of a hydro generator unit, was introduced in [18]. Also, a self-tuning power system stabilizer based on minimization of the quadratic performance index was studied in [19]. A genetic algorithm, for a proportional plus integral control parameters optimization applied to a hydro generator plant, was presented in [20]. Neural Network was successfully considered in an ANN-NARMA-L2 controller and compensating scheme applied to load frequency control of a two-area power system in [21]. The robust control of hydraulic turbines through frequency–response was presented in [22]. Optimal linear time-invariant (LTI) approximation of discrete-time nonlinear systems is studied using a quasistationary signal description in [23]. A significant advancement in LTI control systems design was presented in [24] trough the proposal of the so-called H/sub/spl infin//optimal control theory, applied to a hydroelectric generating speed system based on state-space solution techniques [24]. It is also worth mentioning the work carried out in [25] related with a mothflame optimization control algorithm employed for automatic generation control mechanism of a three-area hydropower system model. Recently, although without the predictive capability of a neural network, a fuzzy systems approach applied to generation power systems has been studied in [26]. Summarizing, and based on the literature reviewed, including classic and modern control proposals, the conventionalbased speed rotation regulation represents still an open problem that needs accurate and novel statements. One of the most challenging issues is related with the slow response of the governor system when returning to the nominal value after a frequency disturbance occurred in the plant. It is in this regard that ANFIS is positioned as one of the most appropriate techniques by combining learning capacities of artificial neural networks and fuzzy inference [27]. The neural capabilities of such hybrid technique make it able to provide accurate predictions based on input–output data collected from the plant, while the data representation is fuzzyified to improve learning and optimization capabilities. An ANFIS based scheme has potential to be used to control the response of the two areas as well as that of the single-area hydropower plant by combining ANFIS-based governorcontrolling techniques [28]. Thereby, the contribution of this work lies in the proposal of a novel methodology to be followed in order to generate a neuro-fuzzy model of a hydro-generation unit (HGU) based on a study carried out using the data of a real hydroelectric plant in Ecuador [29]. The proposed neuro-fuzzy modelling methodology stablishes a common framework for a control strategy applied to Francis-type hydroelectric plants. Unlike most of lookup table-based methods, a novel approach is proposed for PID controller tuning based on the identification process of the system. One of the main advantages of the proposed method is that it requires just the control variation and the electrical power register and, also considering the pole assignment, the programmed maximum opening time for the servomotor. The proposed methodology considers, initially, to perform a cause-effect recording, which allows a continuous time identification of the corresponding LTI type, then, resulting in a classic PID controller capable of complying with the performance specification given by a second order system pole. This enables the training of the proposed neuro-fuzzy model. Indeed, the relevance of using an LTI model lies in the fact that it facilitates the study and analysis of complex systems that can be represented by a mathematical model following linearity and time invariance conditions. According to the carried-out research, because of the redundancy of neuro-fuzzy controller parameters configured in the Proportional-Integral plus Derivative (PI + D) control form, the performance of a classic LTI PI + D controller is improved. Thus, it can be confirmed that any adjustment over the neuro-fuzzy model of the hydrogenerator trained with experimental tests will approximate better the real behavior of the HGU. As it has been aforementioned, this work constitutes a common framework of a functional neuro-fuzzy modelling methodology valid for all subsystems of the HGU operating in a closed loop, which can later be trained to incorporate certain non-linearities due, mainly, to the turbine and the controller. Therefore, the contributions of this paper can be summarized as follows: Optimal Neuro-fuzzy model and PID controller of a Hydro-generator unit from the identification of its LTI model 311 The proposal of detailed methodology for the identification of an LTI-type model for a system. The proposal of the detailed process for the development and training of a neuro-fuzzy model and its PID controller based on the LTI-type system identified. The application of the proposed methodology in a case study (HGU) from the recording of signals from the real systems, system parametrization, stability analysis of the proposed control and validation under disturbances. Fig. 1 illustrates the proposed methodology (with its corresponding sections) for the identification of the LTI model of the HGU and the development and parameterization of the PID control based on neuro-fuzzy logic, which is finally validated by testing it under different disturbances. 2. Modeling of the Hydro-generator unit This section briefly addresses the theoretical aspects necessary to identify the modeling of a Francis-type HGU in continuous time. The modeling aims to determine the main parameters for the identification process to define the PID controller. The parameterized model of the system, together with the PID controller, allow to obtain the neuro-fuzzy models, which aim to reproduce the LTI (Linear and Time-Invariant) behavior of the HGU. Fig. 2 describes the model developed in detail in this section. The black elements correspond to the plant modeling where SM denotes the Servomotor, Tis the turbine and Gis the electric generator. The red path denotes the predictive model of the system which, once properly adjusted, should provide the same electric power output as considering the real plant model. Besides, urecord corresponds to the data recorder in the output of the PID controller in the generator and yrecord is the output of the system (electric power). The signals uand yare recorded input and output signals of the real plant. Being uthe output of the controller, which has been manipulated using limiters, and gradually allows the generation of an increase signal in the position of the servomotor, producing a peak variation of 0.1 in p.u. in a time of 20 s. With this variation, an input signal is achieved for the identification of the system in a frequency range of 0 to 0.45 Hz that allows capturing the energy components of the signals convolved with the system. The result of the disturbance of the real plant is the yrecord, and this signal corresponds to the electrical power of the HGU, which is used to evaluate the adjustment of the parameters of the models located on the red path (Fig. 1) and whose signal is the prediction of electrical power. The use of the PID controller is proposed and a set of rules, applicable to Francis turbines, is generated. These rules are implemented by means of the neuro-fuzzy controller and structured in a PI-D type configuration to avoid control bumps during setpoint changes. 2.1. Statistics used to identify models and systems performance In order to carry out a valuation of the system in temporal samples, it is necessary to apply statistics for sequential data. In this case, the Lennart Ljung criteria are taken as a reference to carry out these comparative studies [14]. The main indicator is the multivariable quadratic correlation factor (r2 y), to indicate that the mean trend level of the signal to be predicted has been suppressed, as observed in (1). r2 y¼PN i¼1b yy   2 PN i¼1yy   2ð1Þ Other statistical criteria used are the root of the mean square error equivalent to standard deviation, and the root of the relative mean square error percentage, represented in (4) and (8) respectively. Consider the i-th deviation between the i-th recorded data and the i-th prediction as shown in (2). ri¼yib yið2Þ The i-th squared deviation allows obtaining a global value based on obtaining its mean value (3). r2¼1 NX N i¼1 r2 ið3Þ Consequently, the square root of the resulting mean square deviation is as follows: Recm ¼ffiffiffiffiffi r2 p¼ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1 NX N i¼1 yib yi ðÞ 2 v u u tð4Þ Considering the i-th recorded value yias realistic, and the ith prediction b yias a prediction, the percentage relative error can be calculated for the case yi–0 and for any integer value of iwith the following equation: eri%¼yib yi yi 100 ¼ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi yib yi yi  2 s100 ð5Þ Therefore, the i-th relative percentage error squared is represented by (6). e2 ri%¼yib yi yi 100  2 ¼yib yi yi  2 1002 ¼yib yi yi 100  2 ð6Þ The relative mean square percentage error can be determined with the following expression: Fig. 1 Proposed methodology for the hydro-generator LTI identification and control. 312 J. Rodrı´guez-Flores et al. e2 r%¼1 NX N i¼1 yib yi yi 100  2 ð7Þ Therefore, an indicator that penalizes the system identification process is the root of the mean squared relative percentage error, which can be calculated as follows: Recr%¼ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1 NX N i¼1 yib yi yi 100  2 v u u tð8Þ 2.2. Decreasing gradient applied to the identification of dynamic LTI systems Research works focused on machine learning as the ones presented in [16–18] use the Nealder-Mead hill descent technique as optimization method, which does not require derivative calculations, and depends on a basic geometric figure that can be formed in a dimension nand with nþ1 sides. Even more intuitive, and not dependent on any derivative process, is the hill climb technique, or valley descent technique, which is the basis of genetic algorithms [18 19]. This technique is intended to produce a constant variation in one direction until the improvement in the cost function to be optimized is achieved, which leads to the variation of the next variable, repeating the cycle until the error condition wanted is met. On the other hand, decreasing gradient method or Newton- Raphson is an optimization method which starts from a random point located on the cost surface and move until they find a minimum which is expected to be the global minimum. Fig. 3 illustrates this optimization search process. As in the methods that not require derivate process, in Newton-Raphson one, it is necessary to establish an algorithm to perform the optimization search. For this case, the criteria to produce the variation of the parameters is bounded between a maximum and minimum step, being necessary to determine the first derivative of the cost function, for the case of the decreasing gradient, and the second derivative for the Newton-Raphson case. The process of obtaining the derivatives of the cost function has been generalized for case studies when the candidate model represents an LTI system. Therefore, a transfer function of order nis assumed, as shown in (9). GsðÞ¼a0snþa1sn1þþan1s1þans0 snþb1sn1þþbn1s1þbns0ð9Þ From (9), the parameters of the numerator and denominator of the transfer function can be generalized and indexed as shown in (10) and (11), respectively. ai)0inð10Þ bi)1inð11Þ Considering (9) and the generalization of aiand bi, the derivative process can be expressed as follows: @b y @ai ¼Zt 0 @ @aiutsðÞL1GsðÞ fg  dsð12Þ @b y @bi ¼Zt 0 @ @biutsðÞL1GsðÞfg  dsð13Þ If the cost function Jis defined as the integral of the square of the error, and the candidate model is considered as LTI- type, the cost function is as shown in (14). J¼1 2Ztf 0 ytðÞRt 0utsðÞL1GsðÞfgds  2dt ð14Þ In (14),utðÞand ytðÞrepresent the input and output variables of the system to be identified, respectively. Thus, the definition of a parameter in a system, when the decreasing gradient method is used, is as shown in (15). On the other hand, when the modified Newton-Raphson method is used, it is as shown in (16), which presents the partial derivative of the cost function with respect to each of the parameters. The matrix representation of (16) is shown in (17) and (18). hinþ1¼hina@J @hi ð15Þ hnþ1¼hnaJ@J @h  1@J @hi  ð16Þ @J @hi¼ Rtf 0ytðÞRt 0uts ðÞ L1Gs ðÞ fgds   Rt 0 @ @hiutsðÞL1GsðÞ fg  ds  dt ð17Þ @J @h¼ @J @h1 @J @h2 . . . @J @hi1 @J @hi 2 6 6 6 6 6 6 6 6 4 3 7 7 7 7 7 7 7 7 5 ð18Þ Fig. 2 Proposed model of the Hydro-generator unit. Optimal Neuro-fuzzy model and PID controller of a Hydro-generator unit from the identification of its LTI model 313 The Jacobian matrix formed with the second partial derivative of the cost function (14) with respect to each of the parameters is shown in (19) and (20). However, if the model is LTI- type, the parameters are adjusted from (17) and (18), for the case of the modified Newton-Raphson method. @ @hj @J @hi  ¼Rtf 0Rt 0 @ @hjutsðÞL1GsðÞ fg  ds   Rt 0 @ @hiutsðÞL1GsðÞ fg  ds  dt Rtf 0ytðÞRt 0utsðÞL1GsðÞfgds   Rt 0 @ @hj @ @hiutsðÞL1GsðÞ fg  ds  dt ð19Þ J@J @h¼ @ @h1 @J @h1  @ @h2 @J @h1  @ @hj1 @J @h1  @ @hj @J @h1  @ @h1 @J @h2  @ @h2 @J @h2  @ @hj1 @J @h2  @ @hj @J @h2  . . .. . ... .. . .. . . @ @h1 @J @hj1 @ @h2 @J @hj1  @ @hj1 @J @hj1 @ @hj @J @hj1  @ @h1 @J @hj  @ @h2 @J @hj  @ @hj1 @J @hj  @ @hj @J @hj  2 6 6 6 6 6 6 6 6 6 6 6 4 3 7 7 7 7 7 7 7 7 7 7 7 5ð20Þ Once defined the expressions to measure the minimization of the cost function, a stop conditions for the optimization process is needed. This condition could be a minimum mean square error, a specific number of cycles, or the impossibility of producing a variation in the parameters that leads to a better value of the cost function. When neuro-fuzzy functions are used, the membership functions (MF) are the relationships proposed to involve the universe of discourse of each of the variables. These MF are related among then by the rules that are generated under the concept of T-norms [20] and the assigned weights or singletons under the Takagi-Sugeno-Kang concept of zero order [1]. Equations (21) and (22) summarize the concept of fuzzy inference and training using the decreasing gradient to each of the singletons for a variable xwith lMF. b fxðÞ¼ X N i¼1 Cilið21Þ Ci¼CiaX N i¼1 liX N i¼1 Cili ! fxðÞ !ð22Þ As example, in the case of two variables (x;y), the case depends on two universes of discourse and a set of membership functions land b. By assigning an indexer jfor l, and an indexer kfor b, the i-th element is expressed as shown in (23), the fuzzy prediction function is as shown in (24), and the training function for the weights or singletons is shown in (25). i¼j1ðÞNbþkð23Þ b fx;yðÞ¼ X Nb Nl j¼1 k¼1 Cj1ðÞNbþk  ljbk  ð24Þ Cj1ðÞNbþk¼Cj1ðÞNbþkaX Nb Nl j¼1 k¼1 ljbk  X Nb Nl j¼1 k¼1 Cj1ðÞNbþk  ljbk  0 B B B B B @ 1 C C C C C A fx;yðÞ 0 B B B B B @ 1 C C C C C A ð25Þ When training a neuro-fuzzy system, it is necessary to consider which rules are activated and which are not, by defining the concept of clustering. Clustering allows developing an adaptive learning system for the algorithm which must determine which rules are activated in order to adapt the algorithm and train only the weights or singletons of these rules. Fig. 4 illustrates the case, in (a) if no rules are activated, the algorithm should not train any weight, on the contrary, in (b) only the weights or singletons of the activated rules should be trained, otherwise will produce a numeric error. Fig. 3 Search process for a minimum point in a cost surface. 314 J. Rodrı´guez-Flores et al. 2.3. Parameterization of the HGU This process considers the dynamic training presented in Fig. 5 for the servomotor, turbine and generator to adjust the operation parameters by applying the decreasing gradient method and considering initial values close to the standards. The numerical convergence towards a global minimum will depend on the initial conditions of the parameters and the selection of aand Dh.ais a learning factor, being much greater than unity when the differential approximation of the cost function is used, and less than unity when using a numerical approximation of the derivative of the cost function. Dhis the differential of the parameter at the time of carrying out the approximate numerical study for both the differential or the derivative of the cost function. The differential behavior of the mechanical power supplied by the turbine in the event of a variation in the position of the servomotor is represented with the differential equation (26). The time constant Twis a value that depends on the operating point of the turbine and is defined as the time it takes for the water circulating through the forced pipe to reach its nominal velocity starting from rest. TW 2 dPmtðÞ dt þPmtðÞ¼WG tðÞTW dWG tðÞ dt ð26Þ The differential behavior of the generator, considering a simplified first-order model, is represented by the differential equation (27). In this equation, the equivalent time Tgis determined based on the time that the generator takes to reach a stable electrical power in the face of a sustained variation of the mechanical power supplied by the turbine. This equation presents the parameter Kgrepresenting the efficiency in the transformation of mechanical power into electrical power in the generator. Tg dPetðÞ dt þPetðÞ¼KgPmtðÞ ð27Þ The behavior of the servomotor is reproduced with the differential equation (28). There are two parameters to be consider: the time of the pilot valve that allows commanding the three-way valve that moves the servomotor (Tsm), and the other parameter (Ksm) is associated with the accumulation of oil in the servomotor which allows, in a mathematical model, convert into an expression that integrates the flow resulting in the position (p.u.). Tsm d2WG t ðÞ dt2þdWG t ðÞ dt þKsm WG tðÞ¼Ksm YPID tðÞ ð28Þ The differential equations (26),(27) and (28) are represented by block diagrams in Fig. 6 for (a) the turbine, (b) generator and (c) servomotor, of the HGU. Fig. 4 Universes of discourse for two variables, with membership functions land b(a) no activated rule (b) excitation of some rules during the performance of the system. Fig. 5 Parameterization of the servomotor, turbine and generator using a register and the decreasing gradient method. Optimal Neuro-fuzzy model and PID controller of a Hydro-generator unit from the identification of its LTI model 315 2.4. Neuro-fuzzy model of the servomotor The model of the servomotor is obtained in three stages. In the first stage, Fig. 7 shows the approach to the fuzzy structure (a) and the differential neuro-fuzzy model of the servomotor (b). In the second stage, Fig. 8 presents the universe of discourse and its parameterized MF (a) which are defined based on the reproduction of the dynamics of the identified model of the servomotor presented in (b). Finally, the third stage showed in Fig. 9, presents the two steps to identify and adjust the neuro-fuzzy differential model: in (a) the fuzzy inference system (FIS) is adjusted as a characterizer, thus creating the initial conditions to adjust the singletons, and in (b) being able to reproduce the LTI dynamic system of the servomotor. Considering that (28) can be represented by (29), by dividing the expression by the time of the servomotor (Tsm), resulting in: d2WG tðÞ dt2þ1 Tsm dWG tðÞ dt þKsm Tsm WG t ðÞ¼Ksm Tsm YPID t ðÞ ð29Þ Equation (29) can be written as shown in equation (30) defining the constants a0,b0and b1. d2WG tðÞ dt2þb1 dWG tðÞ dt þb0WG tðÞ¼a0YPID tðÞ ð30Þ Based on (30) the behavior of the servomotor can be reproduced, with the model presented in Fig. 8(b) where a1and a2 are zero. This LTI model for the servomotor allows to estimate the universe of discussion for the fuzzy model. In this research, the proposed fuzzy systems establish the universe of discussion based on an LTI model. The output of the system yis associated with the state x1, and _ yis associated with the state x2. Therefore, when determining the fundamental harmonic C0and x0for y(or for YPID)thereare possible ranges for the universes of x1and x2:Rang YPID ðÞ¼ max YPID ðÞ,Rang x1 ðÞ¼max yðÞ,Rang x2 ðÞ¼max yðÞx0. It is important to clarify that these are suggested values and must be adjusted depending on the characteristics of the system. During the modeling process, inference systems with membership functions of the triangular type were implemented due to the following reasons: Fig. 6 Block diagrams of the differential models of (a) the turbine, (b) generator and (c) servomotor. Fig. 7 (a) Arquitecture of the FIS to predict the second derivative of the position of the servomotor (b) Differential neuro-fuzzy model of the servomotor. 316 J. Rodrı´guez-Flores et al. They allow the application of the overlapping rule, such that, if only two membership functions of the universe of discourse are activated, their degrees of truth when added always result in unity. The application of triangular membership functions decreases the computational cost and allows a direct defuzzification based on the sum of the product of the singletons by the degree of truth of the activated membership function. Triangular functions allow the implementation of fuzzy systems that behave like piecewise linear functions. Intelligent systems support random initialization, however, as the complexity of the systems increases and they incorporate dynamics, they do not allow random initialization, requiring pre-training (not involving dynamics) as shown in Fig. 9(a). Subsequently, training is carried out to refine the dynamic behavior as shown in Fig. 9(b). 2.5. Neuro-fuzzy model of the turbine Similarly to the servomotor, the model of the turbine is obtained in three stages. Fig. 10 shows the approach to the fuzzy arquitecture in (a) and (b), and the neuro-fuzzy differential model of the turbine in (c). Fig. 11 in (a) presents the universe of discourse and its MF which are parameterized with the reproduction of the dynamics of the identified model of the turbine presented in (b). Fig. 12 presents the two stages to identify and adjust the neuro-fuzzy differential model: in (a) the FIS is adjusted as a characterizer, thus creating the initial conditions to adjust the singletons in (b) and be able to reproduce the dynamic LTI system of the turbine. Considering variations with the maximum amplitude and the maximum gradient of change (in a periodic excitation) the maximum value of the WG input is determined, and its fundamental frequency x0. These parameters make it possible Fig. 8 (a) Membership functions of the servomotor, (b) LTI model to estimate the universe of discourse of the servomotor. Fig. 9 (a) Initialization and (b) dynamic adjustment, of singletons of the neuro-fuzzy differential model of the servomotor. Optimal Neuro-fuzzy model and PID controller of a Hydro-generator unit from the identification of its LTI model 317 To obtain the initial value of Tw(56), the gravity at the location of the hydroelectric plant was determined using the geodesic equation presented in (57) and the gps data for altitude hand latitude l. Data from the hydroelectric plant, such as the gross head, the length and diameters of the penstock and the flow rate measured at full load of the HGU are presented in Table 2. TW¼4QL pgHD2ð56Þ g¼9:780327 1þ5:3024 103sin lðÞðÞ 25:8 106sin 2lðÞðÞ 2  3:086 106hð57Þ The Kgparameter is related to the efficiency of the HGU, considering the losses at both the turbine and generator level and for the proposed case study it was initially defined at 10 %. The value of Tgis an estimate for the previously defined values of Kgand Tw, and a disturbance is evaluated by registering the position of the servomotor WG, and the electric power of the generator Pe. With these determined values, Tgis estimated using equations (46),(47) and (48). ci¼ þTwPei þPei12Ts3Tw ðÞ þPei23Tw4Ts ðÞ þPei32TsTw ðÞ 8 > > > < > > > : 9 > > > = > > > ; ð58Þ fi¼ TSTW2KgWGi1þPei1  þ2TSKgTSþ2TW ðÞWGi2TSTW ðÞPei2  TS2KgTSþTW ðÞWGi32TSTW ðÞPe  8 > < > :ð59Þ Tg¼PN i¼1cifi PN i¼1c2 i ð60Þ Fig. 20 (a) Evolution of the cost function of the RPRMSE of the identification of the Hydro-Agoya ´n HGU, (b) Record of the control signal (yPIDÞ, the electrical power (Pe) and its prediction (Pep). Fig. 21 MF of the universe of discourse of: (a) state x1, (b) position WG, (c) derivative of state Dx1. 324 J. Rodrı´guez-Flores et al. During maintenance operations, smooth perturbations were made to the servomotor and the YPID (system input) and Ysm (system output) signals were recorded. Equations (61)-(65) allow calculating a value to initialize the servomotor parameters to proceed with the joint identification of the model. Table 3 Singletons of the turbine FIS: input fuzzification. UNIVERSE OF WG DISCOURSE lWG;1ðÞ lWG;2ðÞ lWG;3ðÞ UNIVERSE OF THE DISCOURSE OF X1 lX1;1 ðÞ 8.8817e-16 6.5041 13.0083 lX1;2ðÞ 6.5041 0 6.5041 lX1;3ðÞ 13.0083 6.5041 8.8817e-16 Fig. 22 (a) Evolution of the RPRMSE in the process of adjusting the singletons to predict the behavior of the mechanical power of the turbine. Evolution in the adjustment of weights of the singletons in the neuro-fuzzy model of the turbine: (b) input fuzzification, (c) output fuzzification. Optimal Neuro-fuzzy model and PID controller of a Hydro-generator unit from the identification of its LTI model 325 ci¼TsYsmi1Ysmi2 ðÞ ð61Þ f1i¼T2 sYPIDi1Ysmi2 ðÞ ð62Þ f2i¼1ðÞYsmi2Ysmi1þYsmi2 ðÞð63Þ Ksm ¼PN i¼1cif2i  PN i¼1f1if2i  PN i¼1cif1i  PN i¼1f2 2i  PN i¼1f1if2i  2PN i¼1f2 1i  PN i¼1f2 2i ð64Þ Tsm ¼PN i¼1cif1i  PN i¼1f1if2i  PN i¼1cif2i  PN i¼1f2 1i  PN i¼1f1if2i  2PN i¼1f2 1i  PN i¼1f2 2i ð65Þ In Fig. 20(a) the evolution of the cost function of the RPRMSE is shown, and (b) presents: the order of the PID controller (green line), the power signal of the HGU (red line) and the prediction of the fitted model (blue line). Fig. 20(b) shows the validation of the proposed methodology for the identification of the system, where it is evident that it manages to reproduce the response dynamics of the real system. 4.2. Result of the parameterization of the neuro-fuzzy model of the turbine Once defined the LTI model of the turbine, the universe of discourse of the input variables of the neuro-fuzzy system has been determined, as shown in Fig. 21. Once the neuro-fuzzy systems of the turbine have been initialized, it is trained to adjust the weights. Fig. 22 (a) shows the evolution of the cost function of the RPRMSE, where (b) and (c) are the evolution of the singletons of both neuro-fuzzy systems. Table 3 and Table 4 show the values of the singletons of the turbine in the input and output of the neuro-fuzzy system, respectively. Fig. 23 shows how the neuro-fuzzy model of the turbine fits the LTI model, being the statistical data of the training the one presented in Table 5. Table 4 Singletons of the turbine FIS: output fuzzification. UNIVERSE OF DX1 DISCOURSE lWG;1ðÞ lWG;2ðÞ lWG;3ðÞ UNIVERSE OF THE DISCOURSE OF X1 lX1;1ðÞ 0.9199 48.9199 lX1;2ðÞ 4.9199 0 4.9199 lX1;3ðÞ 8.9199 4 0.9199 Fig. 23 Comparison between an LTI model (Pm) and a neurofuzzy model (Pmp) for the prediction of the mechanical power of the turbine. Table 5 Statistical indicators of adjustments of the neurofuzzy model of the turbine. Root of Initial Percentage Relative Mean Square Error 263.5644 Root of Final Relative Percentage Mean Square Error 1.1356e-05 Quadratic Multivariate Correlation Factor r2 y 1 Fig. 24 Membership functions of the universe of the discourse of the: (a) state x1, (b) electric power of the generator. 326 J. Rodrı´guez-Flores et al. 4.3. Result of the parameterization of the neuro-fuzzy model of the generator Once defined the LTI model of the generator, the universe of discourse of the input variables of the neuro-fuzzy system has been determined, as shown in Fig. 24. Once the neuro-fuzzy system of the generator has been initialized, it is trained to adjust the weights. Fig. 25 shows (a) the evolution of the cost function of the RPRMSE, where (b) is the evolution of the singletons of the neuro-fuzzy system. Table 6 presents the values of the singletons of the neuro-fuzzy system of the generator. Fig. 26 shows how the neuro-fuzzy model of the generator fits the LTI model. Being the statistical data of the training presented in Table 7. 4.4. Result of the parameterization of the neuro-fuzzy model of the servomotor Fig. 27 depicts the universe of discourse and MF of the variables of the neuro-fuzzy system of the servomotor. Fig. 25 (a) Evolution of the RPRMSE in the process of adjusting the singletons to predict the behavior of the electric power of the generator. (b) Evolution in the adjustment of weights of the singletons in the neuro-fuzzy model of the generator. Table 6 Singletons of the generator FIS. UNIVERSE OF Pm DISCOURSE lPm;1ðÞ lPm;2ðÞ lPm;3ðÞ UNIVERS OF Pe DISCUOURSE lPe;1 ðÞ 0.7960 13.2868 25.7777 lPe;2ðÞ 12.4908 0 12.4908 lPe;3ðÞ 25.7777 13.2868 0.7960 Fig. 26 Comparison between an LTI model (Pe) and a neurofuzzy (Pep) model during the prediction of the electrical power of the generator. Table 7 Statistical indicators of adjustments of the neurofuzzy model of the turbine. Root of Initial Percentage Relative Mean Square Error 35.7123 Root of Final Relative Percentage Mean Square Error 0.0126 Quadratic Multivariate Correlation Factor r2 y1 Optimal Neuro-fuzzy model and PID controller of a Hydro-generator unit from the identification of its LTI model 327 Fig. 27 Membership functions of the universe of discourse of: (a) state x1, (b) position WG, (c) derivative of state Dx1. Fig. 28 (a) Evolution of the RPRMSE in the process of adjusting the singletons to predict the position of the servomotor. (b)(c)(d) Evolution of the weight adjustment of the singletons of the neuro-fuzzy model of the servomotor. 328 J. Rodrı´guez-Flores et al. Fig. 29 Three-dimensional table structure and location of singletons. Table 8 Servomotor FIS. UNIVERSE OF DISCOURSE OF Ypid lYPID;1ðÞ UNIVERSE OF DISCOURSE OF Sm lSm;1ðÞ lSm;2ðÞ lSm;3ðÞ UNIVERSE OF DISCOURSE OF DSm lDSm;1ðÞ 79.8730 378.5283 836.9297 lDSm;2ðÞ 5.6843e-14 458.4014 916.8028 lDSm;3ðÞ 79.8730 538.2744 996.6759 Table 9 Servomotor FIS. UNIVERSE OF DISCOURSE OF Ypid lYPID;2ðÞ UNIVERSE OF DISCOURSE OF Sm lSm;1ðÞ lSm;2ðÞ lSm;3ðÞ UNIVERSE OF DISCOURSE OF DSm lDSm;1ðÞ 538.2744 79.8730 378.5283 lDSm;2 ðÞ 458.4014 0 458.4014 lDSm;3ðÞ 378.5283 79.8730 538.2744 Table 10 Servomotor FIS. UNIVERSE OF DISCOURSE OF Ypid lYPID;3ðÞ UNIVERSE OF DISCOURSE OF Sm lSm;1ðÞ lSm;2ðÞ lSm;3ðÞ UNIVERSE OF DISCOURSE OF DSm lDSm;1ðÞ 996.6759 538.2744 79.8730 lDSm;2 ðÞ 916.8028 458.4014 5.6843e-14 lDSm;3ðÞ 836.9297 378.5283 79.8730 Optimal Neuro-fuzzy model and PID controller of a Hydro-generator unit from the identification of its LTI model 329 Once the neuro-fuzzy system of the servomotor has been initialized, it is trained to adjust the weights. Fig. 28 shows (a) the evolution of the cost function of the RPRMSE, where (b) is the evolution of the singletons of the neuro-fuzzy system. Considering that the representation in four dimensions is a complex issue, the membership functions were defined in a three-dimensional space, with three membership functions per dimension (Fig. 29). This representation defines 27 rules or singletons that define a fourth dimension (position of the servomotor). This fourth dimension is represented with values distributed in 9 values by layers, which are presented in 3 tables. In Tables 8-10 the values of the singletons of the neurofuzzy system of the servomotor are presented. Fig. 30 shows how the neuro-fuzzy model of the servomotor fits the LTI model, being the statistical data of the training the one presented in Table 11. 4.5. Result of the parameterization of the neuro-fuzzy model of the PID controller Based on the equations proposed for the linear PID controller, the parameters were initialized to apply the decreasing gradient method and improve performance with respect to the pattern given by the specified second-order pole. The statistical results of this comparison are observed in Table 12 and the comparative temporal response between the standard signal, the initial and final adjustment conditions of the HGU controller, is observed in Fig. 31. With the LTI model of the controller already determined, the universe of the discourse of the neuro-fuzzy controller variables are defined, as shown in Fig. 32. Once the neuro-fuzzy system of the PID controller has been initialized, it proceeds to its training to adjust the weights. Fig. 33 shows (a) the evolution of the RPRMSE, being (b) the evolution of the singletons of the neuro-fuzzy system. Tables 13-15 show the values of the singletons of the neurofuzzy system of the PID controller output. Fig. 34 shows how the neuro-fuzzy model of the PID controller adjusts to the specified pole pattern model, before and after training the weights for the singletons. The statistical data of the training process is presented in Table 16. 4.6. Stability analysis of the control system PID LTI and PID Neuro-Fuzzy The stability study of the HGU, with the PID LTI serial controller, can be studied with the Bode criteria, taking as reference the optimal values determined during the evaluation of the system. For this, 4 equations are determined to establish Fig. 30 Comparison between an LTI model and the neuro-fuzzy model of the prediction of electrical power of the generator. Table 11 Statistical indicators of adjustments of the neurofuzzy model of the servomotor. Root of Initial Percentage Relative Mean Square Error 15.7358 Root of Final Relative Percentage Mean Square Error 0.6552 Quadratic Multivariate Correlation Factor r2 y1 Table 12 Initial and final parameters of the LTI PID controller and statistical values of the system performance. Assigned Pole PID LTI Controller Initial Setting rx dKp 0.0496 Ti 2.9583 Td 0.0018 Final adjustment of the PID LTI controller 0.4 + 0.5457i Kp 0.2031 Ti 2.5502 Td 0.0263 Initial CV%Final CV%Initial R2 yFinal R2 y 17.6411 11.8770 0.5958 0.8167 Fig. 31 Response of the HGU simulation with the initial parameters and after adjusting using the decreasing gradient method. 330 J. Rodrı´guez-Flores et al. Fig. 32 Membership functions of the universe of discourse for: (a) integral error, (b) proportional error, (c) derivative error. Fig. 33 (a) Evolution of the RPRMSE cost function in the process of adjusting the singletons to predict the behavior of the controller output. (b) Evolution in the weight adjustment of the Singletons of the neuro-fuzzy model in the PID controller. Table 13 Singletons of the FIS of the controller output. UNIVERSE OF THE DISCOURSE OF THE DERIVATIVE OF ERROR lDerErr;1ðÞ UNIVERSE OF THE DISCOURSE OF PROPORTIONAL ERROR lProErr;1ðÞ lProErr;2ðÞ lProErr;3ðÞ UNIVERSE OF THE DISCOURSE OF THE INTEGRAL OF ERROR lIntErr;1ðÞ 4.2275 2.9102 1.5929 lIntErr;2ðÞ 1.71365 0.3963 0.9209 lIntErr;3ðÞ 0.8002 2.1175 3.4348 Optimal Neuro-fuzzy model and PID controller of a Hydro-generator unit from the identification of its LTI model 331 the stability range of the parameters Kp,Tiand Td, of the PID LTI controller. Equation (66) represents the transfer function of the servomotor, the turbine and the generator. Equation (67) is the equivalent transfer function for Kpthat allows determining the relative stability for this parameter. Equation (68) is the equivalent transfer function for Tithat allows determining its relative stability. Finally, (69) is the equivalent transfer function for Tdthat allows determining its relative stability. GSTG sðÞ¼ Ksm sT smsþ1ðÞ  1þKsm sT smsþ1ðÞ  1TWsðÞ 1þTW 2s  Kg Tgsþ1  ð66Þ GeqKp sðÞ¼ GSTG sðÞ 1 TisþTds  GSTG sðÞþ1 ð67Þ GeqTi sðÞ¼sK pþTds  GSTG sðÞþs GSTG sðÞ ð68Þ GeqTd sðÞ¼ sGSTG sðÞ Kpþ1 Tis  GSTG s ðÞ þ1 ð69Þ Table 17 shows the values for the stability range for Kp,Ti and Td, of the PID controller connected in series with the HGU. Fuzzy control systems are essentially nonlinear ones. For this reason, it is difficult to obtain general results on the stability analysis of controllers designed under this technique. Furthermore, the knowledge of the dynamic behavior of the process to be controlled is normally poor. Therefore, the robustness of the fuzzy control system must be studied to ensure stability despite variations in the process dynamics. For the case of two linguistic variables, the analysis can focus on a geometric interpretation of the state maps. This technique is based on the study of field vectors and the rules of the control system based on fuzzy knowledge. When the fuzzy system has more linguistic variables, the analysis is based on the field vector defined by the partial derivative of the con- Table 14 Singletons of the FIS of the controller output. UNIVERSE OF THE DISCOURSE OF THE DERIVATIVE OF ERROR lDerErr;2ðÞ UNIVERSE OF THE DISCOURSE OF PROPORTIONAL ERROR lProErr;1ðÞ lProErr;2ðÞ lProErr;3ðÞ UNIVERSE OF THE DISCOURSE OF THE INTEGRAL OF ERROR lIntErr;1ðÞ 3.8312 2.5139 1.1966 lIntErr;2ðÞ 1.3172 0 1.3172 lIntErr;3ðÞ 1.1966 2.5139 3.8312 Table 15 Singletons of the FIS of the controller output. UNIVERSE OF THE DISCOURSE OF THE DERIVATIVE OF ERROR lDerErr;3ðÞ UNIVERSE OF THE DISCOURSE OF PROPORTIONAL ERROR lProErr;1ðÞ lProErr;2ðÞ lProErr;3ðÞ UNIVERSE OF THE DISCOURSE OF THE INTEGRAL OF ERROR lIntErr;1 ðÞ 3.4348 2.1175 0.8002 lIntErr;2ðÞ 0.9209 0.3963 1.7136 lIntErr;3ðÞ 1.5929 2.9102 4.2275 Fig. 34 Comparison of the pattern system with the specified pole vs the HGU with neuro-fuzzy control. Table 16 Statistic performance comparison before and after training the weights of the neuro-fuzzy controller. Initial CV%Final CV%Initial R2 yFinal R2 y 17.6411 10.9291 0.5958 0.8448 332 J. Rodrı´guez-Flores et al. trolled variable with respect to the deviation or error of the control system. Said field vector can be approximated as shown in (70). @dPe dt @YE ¼ L1PesðÞ s 5Tssþ1  2  L1YesðÞ s 5Tssþ1 no ð70Þ If the error of the controlled system is considered as the difference between the power setpoint and the system output (in this case the generator power), this can be determined using (71). Ye¼SetPoint Peð71Þ Finally, in an approximate way, the derivative of the power can be determined with (72). dPe dt ¼L1PesðÞ s 0:25Tssþ1  ð72Þ The behavior of the derivative of the power with respect to time versus the control error, and the direction and sense of the field vector, allow determining the absolute stability condition of the system. This behavior study is carried out in the operating range and considering the model of the system under test in a symmetrical manner considering the setpoint values of 0.25, 0.5, 1 and 0.25, 0.5, 1. Under this analysis, Fig. 35(a) shows the behavior of the power change gradient with respect to the setpoint error, with its punctual representations of the field vector. On the other hand, Fig. 35(b) shows an approach to the point of convergence, where the field vector points towards equilibrium at the origin of the pair (0,0). 4.7. Behavior of the neuro-fuzzy model of the Hydro-Agoya ´n HGU during power setpoint variations The model proposed in this paper and characterized with the data of the case study allows to reproduce the behavior of the Hydro-Agoya ´n HGU, under the concept of a nominal neuro-fuzzy reduced-order model. To verify its ability to reproduce the real behavior, it is initialized at 50 % (0.5 [p. u.]) and several variations in the power setpoint of ± 10 % (0.1 [p.u.]) are applied. The neuro-fuzzy model of the HGU is presented in Fig. 36(a) and the simulation result is observed in Fig. 36(b). The models presented in this paper have been parameterized for the Hydro-Agoya ´n HGU and the main results have been presented under these plant considerations. However the methodology for the parameterization and the structures of the neuro-fuzzy system are of a general nature and applicable both to Francis-type and Kaplan-type turbines, although their application can be extrapolated, with a little more limitations, to Pelton turbines. The neuro-fuzzy models have been proposed to reproduce the behavior of LTI models in order to follow a training procedure and incorporate some nonlinearities that can be detected with this structure. Additionally, disturbances in the power setpoint of 5 %, 15 % and 20 % have been considered, the results of which are shown in Fig. 37. Nevertheless, it is important to mention that the consideration of other conditions in the models, not only does it change the fuzzy structure but also it requires the consideration of new universes of discourses and the incorporation of a large number of additional rules. Therefore, the selection of the main conditions to be included in the modeling of this kind of systems must be evaluated in terms of cost versus benefit (in terms Fig. 35 Stability analysis: (a) behavior of the change gradient and (b) approach to the point of convergence. Table 17 Stability Range for the parameters of the PID LTI. Parameter Range Kp:½ 00:6205½ Tis½ 1:2156 1½ Tds½ 00:1450½ Optimal Neuro-fuzzy model and PID controller of a Hydro-generator unit from the identification of its LTI model 333