Optimal Fuzzy Controller Based on Chaotic Invasive Weed Optimization for Damping Power System Oscillation

Abstract This paper proposed an optimal fuzzy controller for damping the power system oscillation in multi-machine environment. In this strategy, the proposed controller is optimized by new chaotic invasive weed optimization algorithm. Furthermore, a new objective function has been considered to test the proposed controller in different load conditions which increase the stability of system after disturbances. For this purpose, the damping factor, damping ratio, and a combination of the damping factor and damping ratio were analyzed and compared with the proposed objective function. The effectiveness of the proposed strategy has been applied in two multi-machine power system test cases as three machine 9-bus IEEE standard and 10-machine 39-bus New England power systems. The eigenvalue analysis and nonlinear time-domain simulation results proof the effectiveness of the proposed strategy. Schematic seed production in a colony of weeds

Recently, the research for control method based on fuzzy logic controllers (FLC) as PSS has greatly improved the dynamic characteristics of power system [8]. Fuzzy rules and membership functions shape should be adjusted to obtain the best control performance in FLC. Conventionally, the adjustment is done by the experience of experts or trial and error methods. Therefore, it is difficult to determine the suitable membership functions and rule base without the knowledge of the system. These problems make the design process more difficult [9]. In [10] fuzzy logic-based controller has been proposed for solving the stability problem in multi-machine power system environment. In the design of fuzzy logic controllers, unlike most conventional methods, a mathematical model is not required to describe the system under study. In fuzzy control, it is synthesized from a collection of fuzzy If-Then rules which describe the behavior of the unknown plant. It has been applied to the design of PSSs in a number of publications [11].
In this paper to tackle the mentioned problem in controllers and meta-heuristic algorithms, we proposed a new optimization algorithm which is applied over fuzzy controller. So, we can claim the contributions of this manuscript in literature as: (1) A new meta-heuristic algorithm has been proposed in this paper which applied over fuzzy controller in multi-machine environment. (2) The proposed controller has been tested over 2 case studies of 3-machine 9-bus IEEE power system and 10 machine 39 bus New England power system with a new objective functions. In this method, feedback signal has been considered from a remote machine and local and remote input signal ratios for each machine in a multi-machine power system under various operating conditions. (3) In proposed fuzzy controller, the self-tuning mechanism has been proposed which decreases the equivalent integral control component of the fuzzy controller gradually with the system response process time.
The effectiveness of the proposed algorithm has been tested over two case studies in different load conditions and faults. The remaining parts of the paper are organized as follows. The problem statement has been presented in Section 2. The proposed fuzzy controller is presented in Section 3. Section 4 presents the proposed new meta-heuristic algorithm with the self-tuning mechanism over fuzzy controller. The numerical results obtained from the proposed strategy are presented in Section 5. Finally, Section 6 concludes the paper.

Introduction
Instability problems are caused mainly by the low-frequency oscillations developed in the interconnected power system network. The major causes for this poor low-frequency oscillations are sudden load changes, switching surges, both internal and external faults, transient disturbances, etc., Damping of these oscillations in interconnected power system is essential for secure and stable operation of the system [2]. High gain excitation system in addition with automatic voltage regulators (AVR) will damp out the oscillations of low frequency (0.2-2.5 Hz), but considerably it is favor only for the improvement of steady-state stability rather than dynamic stability. Thus, an auxiliary signal in addition with the excitation system and AVR will greatly improve the system dynamic stability which is termed as power system stabilizer (PSS). PSSs are auxiliary control devices on synchronous generators, used in conjunction with their excitation systems to provide control signals toward enhancing the system damping and extending power transfer limits.
Several techniques have been proposed by researchers for tuning the PSS parameters to enhance the stability in power systems. Conventional power system stabilizers are designed based on eigenvalue analysis which utilizes two basic tuning techniques phase compensation and root locus [1]. Phase compensation is widely used and compensates for the phase lags by providing a damping torque component. Root locus involves shifting of eigenvalues related to the power system modes of oscillation by shifting the poles and zeros of the stabilizer [1]. In [2] application of conventional lead-lag PSS (CPSS) has been proposed which is usually designed with a fixed gain, with an aim to stabilize at the nominal operating condition. However, the inherent nonlinearity and multiple operating points of a power system degrade the performance of such a fixed gains CPSS. In [3], H ∞ optimization techniques have been presented for PSS tuning problem. However, the additive and/or multiplicative uncertainty representation cannot treat situations where a nominal stable system becomes unstable after being perturbed. On the other hand, the order of the H ∞ -based stabilizer is as high as that of the plant. This gives rise to complex structure of such stabilizers and reduces their applicability. Particle swarm optimization (PSO) has been proposed in [4,5] and genetic algorithm (GA) in [6,7] for the mentioned problem. Although GA is very sufficient in finding global or near global optimal solution of the problem, it requires a very long-run time that may be several minutes or even several hours depending on the size of the system under study (Appendix 1).

Power System Modeling
In this paper, a multi-machine power system is considered as a test case where the third-order model is presented in [12,13]. Actually, the proposed power system consists of 10 generators and the electrical and mechanical part of ith generator is modeled as follow: Also, it can be presented after mathematical transformers where, ΔP e , Δw, and ΔV t are quantities [12,13]; where ΔP ei is the state deviation in generator electromagnetic power for the ith subsystem; Δw i is the state deviation in rotor angular velocity for the ith subsystem; ΔV ti is the state deviation in the terminal voltage of the generator for the ith subsystem: The linearized rotor motion equation for synchronous generator can be described as: where, J, D : the rotor inertia and the damping factor; T ′ doi : the direct axis transient time constant; x di : the direct axis reactance; x ′ di : the direct axis transient reactance; V ti : the terminal voltage; x dsi : the is the total reactance which takes into account x d ; V si : the voltage at the infinite bus; x si : the line impedance; δ i : the power angle of the ith generator, in rad; w i : the relative speed of the ith generator, in rad/s; E fi : the equivalent EMF in the excitation coil; E qi : the transient EMF in the quadrature axis; E qi : the EMF in the quadrature axis; P m : the mechanical input power; P ei : the active power delivered by the ith generator; Q ei : the reactive power of the ith generator; I di , I qi : the direct and the quadrature axis stator currents; I fi : the excitation current of the ith generator; x adi : the mutual reactance between the excitation coil and the stator coil; G ij + B j : the ith row and jth column element of nodal transient admittance matrix after eliminating all physical busbars; G ii : the self-conductance of the ith node; B ii : the self-admittance of the ith node; ΔT m is the mechanical input torque, ΔT e (11) Which can be presented as; When the function value of the above equation is equal to or less than zero, the response is that the maximal damping factors in each operating conditions are exactly or less than the excepted value σ 0 [15]. According to the Figure 1, all damping factors are limited in a rectangular region that is not more than σ 0 .
(2) For the second objective function with the damping ratio: The minimal damping ratio in all eigenvalues can be defined by; max 1≤q≤n q q , where ζ q is the damping ratio of the qth eigenvalue. The expected damping ratio constant ζ 0 (positive value) in this equation can be presented as; where, the difference values of all operating conditions can be presented as; Which can be presented as; A fan-shaped region with the tip at the origin can be formulated as; ≡ − √ 2 + 2 × 100%, where, as it is described in part (b) of figure one, all damping ratios are no less than ζ 0 . Also, the slope of the straight line with inverse ratio by ζ 0 can be presented as; is the electromagnetic torque and ΔT e = K 1 Δδ + K 2 Δδ, By neglecting the K 2 ΔE ′ q , the formulation can be described as; ΔT e = K 1 Δδ + ΔT D , D is the natural damping constant. Accordingly, the above equation after Laplace transformer and Δw = sΔδ/w 0 can be described as; Which can be described as; OR Accordingly, we can achieve the following equation from above equations; where ξ n is the damping factor; w n is the undamped mechanical oscillation frequency.

Objective Function Problem
Actually, in multi-machine power system environment, damping ratio over PSS controller is the main problem. So, in these kinds of problems the objective function sets to enhance the stability in power system [14][15][16][17][18][19][20][21]. Where, when a disturbance occurs in a power system, the duration of an oscillation in the time domain is determined by several eigenvalues at the most right sides on the plane. For this purpose, three different kinds of objective functions have been proposed in this paper which has been used in several research manuscript as: (1) The first objective function with the damping factor may be derived as follows: In a power system, if the number of eigenvalues is n q under the yth operating condition, then the maximal damping factor in all eigenvalues is represented as; max 1≤q≤n q q , where σ q is the damping factor of the qth eigenvalue. According to this equation, the subtracts expected for damping factor constant σ 0 (negative value) to obtain a difference value as; max 1≤q≤n q q − 0 . If the total operating conditions of the test system are n y , the difference values of all operating conditions are summed as: where χ q is the damping scale of the qth eigenvalue, and χ is a constant value of the expected damping scale. This function region with the tip at σ 0 is presented in section (d) of Figure 1. Also, if we assume the σ 0 as 0, Slope 0 is equal to Slope 0 ; otherwise, σ 0 is a negative value and the magnitude of is smaller than that of Slope 0 . If σ < σ 0 and ζ < ζ 0 we can write:

Fuzzy logic controller
Recently, the research for control method based on fuzzy logic controllers (FLC) as PSS has greatly improved the dynamic characteristics of power system [10]. Fuzzy rules and membership functions shape should be adjusted to obtain the best control performance in FLC. Conventionally, the adjustment is done by the experience of experts or trial and error methods. Therefore, it is difficult to determine the suitable membership functions and rule base without the knowledge of the system. These problems make the design process more difficult [11]. Because of the complexity and multi-variable conditions of the power system, conventional control methods may (27) Slope The third objective function for damping ratio or damping factor can be evaluated as: where ϖ is the weight for combining both damping factor and damping ratio. This equation forms a D-shaped region, as presented n section (c) of Figure 1, which all damping factors and damping ratios are no more than σ 0 and no less than ζ 0 , respectively. In this paper, we applied a new objective function which is obtained by subtracting σ 0 from the eigenvalue, and was then derived as the damping ratio equation by considering the damping scale (χ) as: Thus, the equivalent control components of the PID-type FLC are obtained as follows: Proportional gain = αK e P + βK d D, Integral gain = βK e P, Derivative gain = αK d D

PID Self-tuning Mechanisms
In this self-tuning method [22], the PID-type FC integral and derivative components updating are achieved based on scaling factors β and Kd, using the information on system's error as: where, β 0 and K d0 are the initial values of β and K d′ , respectively, Φ(.) and Ψ(.) are the empirical tuner functions defined, respectively, by: In these equations, the parameters to be tuned φ 1 , φ 2 , ψ 1 , and ψ 2 are all positive. The empirical function related to integral component decreases as the error decreases while the function related to derivative factor increases. Indeed, the objective of the function is to decrease the parameter with the change of error. However, the function has an inverse objective to make constant the proportional effect.

Self-tuning mechanism by Relative Rate Observer Method
In this self-tuning method [23,24], the PID controller integral and derivative components updating are achieved as: where δ k is the output of the fuzzy relative rate observer (RRO), K f is the output scaling factor for δ k and K fd is the additional parameter that affects only the derivative factor of the controller. The rule base for δ k , [23], is considered for the fuzzy RRO. This fuzzy RRO block has as inputs the absolute values of error |e k | and the variable r k , defined subsequently, as shown in Table 1.
The linguistic levels assigned to the input |e k | and the output variable δ k are as follows: L (Large), M (Medium), give satisfactory solutions. On the other hand, their robustness and reliability make fuzzy controllers useful for solving a wide range of control problems in power systems. The controller block is formed by the fuzzification of (∆ω i ), the interface mechanism, and the defuzzification. Therefore, Ui is a control signal that applies to the excitation system in each generator. By taking ∆ω i as the system output, the control vector for the conventional PID controller is given by: The parameters, K Ii , K di , and K pi are determined by a set of fuzzy rules of the form: , and E i are fuzzy sets on the corresponding supporting sets.
Basically, structural parameters include input/output variables to fuzzy inference, fuzzy linguistic sets, membership functions, fuzzy rules, inference mechanism, and defuzzification mechanism. Tuning parameters include input/output scaling factors and parameters of membership functions. Usually, the structural parameters are determined during offline design while the tuning parameters can be calculated during online adjustments of the controller to enhance the process performance, as well as to accommodate the adaptive capability to system uncertainty and process disturbance. In this study, we will deal with fuzzy PID-type controllers formed using one PD-type FLC with an integrator at the output.
The new method adjusts the input scaling factor corresponding to the derivative coefficient and the output scaling factor corresponding to the integral coefficient of the PID-type FLC in a similar fashion to the methods mentioned above. The adjustments are done through a fuzzy inference mechanism in an online manner.
The output of the fuzzy PID-type FLC is given by; u = U + ∫ Udt. Where, U is the output of the FLC and the value of parameters are presented in paper. Moreover, the fuzzy controllers with product-sum inference method, center of gravity defuzzification method and triangular uniformly distributed membership functions for the inputs and a crisp output, the relation between the input and the output variables of the FLC is given by; where, E = K e e and Ė = K dĖ . Therefore, we can write the output of controller as; (28) u = A + At + K e Pe + K d De + K e P ∫ edt + K d Ḋe like tuning of robot controller [25], optimal positioning of piezoelectric actuators [26], development of recommender system [27], antenna configuration optimization [28], design of E-shaped MIMO antenna [29], and DNA computing [30]. IWO is a meta-heuristic algorithm which mimics the colonizing behavior of weeds. In invasive weed optimization algorithm, the process begins with initializing a population. It means that a population of initial solutions is randomly generated over the problem space. Then members of the population produce seeds depending on their relative fitness in the population. In other words, the number of seeds for each member is beginning with the value of S min for the worst member and increases linearly to S max for the best member [25]. This technique can be summarized as: A. Initialization In this step, a finite number of weeds are initialized at the same element position of the conventional array which has a uniform spacing of 'γ/2' between neighboring elements.

B. Reproduction
The individuals, after growing, are allowed to reproduce new seeds linearly depending on their own, the lowest, and the highest fitness of the colony (all of plants). The maximum (S max ) and minimum (S min ) number of seeds are predefined parameters of the algorithm and adjusted according to structure of problem. The schematic seed production in a colony of weeds is presented in Figure 3. In this figure, the best fitness function is the lower one [27].

C. Spatial distribution
The generated seeds are being randomly distributed over the d-dimensional search space by normally distributed random numbers with mean equal to zero; but varying variance. This step ensures that the produced seeds will be generated around the parent weed, leading to a local search around each plant. However, the standard deviation (SD) of the random function is made to decrease over the iterations, which is defined as: SD max and SD min are the maximum and minimum standard deviation, respectively. And pow is the real no. In this equation, the 'iter' define the iteration in proposed algorithm. This step ensures that the probability of dropping a seed in a distant area decreases nonlinearly with iterations, which result in grouping fitter plants and elimination of inappropriate plants.

D. Competitive exclusion
When the maximum number of population in a colony is reached (Pmax), each weed can produce seeds and The variable r k , called normalized acceleration, gives 'relative rate' information about the fastness or slowness of the system response as shown in Table 2, as: where Δe k and Δ(Δe k ) are the incremental change in error and the so-called acceleration in error given, respectively, by: where, Δe * is chosen as: For this self-tuning approach, the uniformly distributed triangular and the symmetrical membership functions are assigned for the fuzzy inputs r k and |e k |, and fuzzy output variable δ k . The view of the above fuzzy rule base is illustrated in Figure 2

Modified Invasive Weed Optimization
Invasive weed optimization (IWO) is inspired from weed colonization and motivated by a common phenomenon in agriculture that is colonization of invasive weeds. Weeds have shown very robust and adaptive nature which turns them to undesirable plants in agriculture. Since its advent IWO has found several successful engineering applications Step 1: Generate an initial chaos population randomly for chaotic local search.
where, the chaos variable can be generating as follows: Step 2: Determine the chaotic variables.
where, N chaos is the number of individuals for chaotic local search. cx N g i is the ith chaotic variable. Rand (.) generate a random number between [0-1].
Step 3: Mapping the decision variables. (40) Ng 0 cx j 0 = X j ch,0 − P j,min P j,max − P j,min , j = 1, 2, … , N g (42) spread them according to the mechanism mentioned in step 3.2 and step 3.3, respectively. Then, new seeds with their parents are ranked together with respect to their fitness. This method is known as competitive exclusion and is also a selection procedure of IWO. Next, weeds with lower fitness are eliminated to reach the maximum allowable population size in a colony. This mechanism using the 'survival of the fittest' idea [29] gives a chance to plants with lower fitness to reproduce, and if their offsprings have good fitness, they can survive in their offspring's existence.

E. Termination condition
The whole process continues until the maximum number of iterations has been reached, and we hope that the plant with the best fitness is the closest one to the optimal solution.

F. Chaotic local searc
Because of the properties of unique ergodicity, inherent stochastic property, irregularity of chaos, a chaotic search can overpass every state in a certain space, and every state is visited only once, which is useful to prevent being trapped in local optima. Accordingly, this search technique is applied over IWO to modify its local search capabilities. For chaotic local search, there are two types of  are evaluated by Figure of demerit (FD) criteria as defined in [2]. In the case of the PID-type FC structure with the Empirical Functions Tuner Method self-tuning mechanism, the scaling factors to be optimized are ϕ1, ϕ2, ψ1, and ψ2. The formulated optimization problem is defined as follows:

Numerical Results
In this section, two case studies have been presented by applying the proposed control strategy. For this purpose, we present each case with the provided results in each subsections in the following:

Three-machine Nine-bus Power System
In this study, the three-machine nine-bus power system shown in Figure 4 is considered. Furthermore, operating of system is tasted in three different conditions as nominal, lightly, and heavily loading conditions. Details of the system data and operating condition are given in Ref. [32].
In this case, we applied the proposed controller over three generators where all of them have consists one controller. The convergence trend of proposed algorithm has been presented in Figure 5. Also, the obtained fuzzy rule base for the output u fz and its figure is presented in Table  3 and Figure 6, respectively. Also, for more analysis, the proposed controller has been tested in different load condition with comparison by genetic algorithm (GA) and robust PSS with PSO (RPSSPSO) in [32].
In this section, we applied a 6-cycle 3-phase fault at t = 1 s, on bus 7 at the end of line 5-7. The fault is cleared Step 4: Convert the chaotic variables to the decision variables.
Step 5: Evaluate the new solution with decision variables.

Proposed CIWO-Based Controller
The problem of scaling factors tuning is presented in this section for all defined PID-type FLC structures by proposed meta-heuristic algorithm.

PID-type FC Tuning Problem Formulation
The structure of fuzzy controller in most of the times set by trial and error. For this purpose, the optimization of these scaling factors is proposed like a promising solution. This tuning problem can be formulated as the following constrained optimization problem: x min ≤ x ≤ x max } the initial search space, which is supposed containing the desired design parameters, and g l : R m → R the problem's constraints.
The optimization-based tuning problem consists in finding the optimal decision variables x * = (x * 1 , x * 2 , … , x * m ) T , representing the scaling factors of a given PID-type FC structure, which minimize the defined cost function, chosen as the performance criteria. These cost functions are minimized, using the proposed meta-heuristic algorithm, under various time-domain control constraints such as overshoot D, steady-state error Ess, rise time tr, and settling time ts of the system's step response.
Hence, in the case of the PID-type FC structure without self-tuning mechanisms, the scaling factors to be optimized are Ke, Kd, α, and β. The formulated optimization problem is defined as follows: where D max , E max ss , t max r and t max s are the specified overshoot, steady state, rise, and settling times, respectively, that constraint the step response of the IWO-tuned PID-type FC controlled system, and can define some time-domain templates. Furthermore, in this paper the objective function assess the effectiveness and robustness of the proposed method over a wide range of loading conditions, three different cases designated as nominal, light, and heavy loading are considered. Details of the system data and operating condition are given in Ref. [32]. The G 5 , G 7 , and G 9 are indicated for the place of CPSS. The proposed model has been compared with GA and no PSS techniques which are reported in [32]. For this purpose, performance of the proposed strategy under transient conditions is verified by applying a 3-cycle 3-phase fault at t = 1 s, on bus 25 at the end of line 25-26. The fault is cleared by permanent tripping the faulted line. For obtained results from simulation, Figure 9, is presented to show the speed deviations of the generators G 5 , G 7 , G 9 , and G 11 under heavy load condition. The information of load for this power system is presented in [32]. Also, the convergence trend of proposed algorithm is presented in Figure 10. The obtained fuzzy rule base for the output u fz is presented in Table 5.
The oscillation frequencies, damping factors, damping ratios, and damping scales over four objective functions by permanent tripping the faulted line. Speed deviations of the generators G 1 to G 3 under light load condition are shown in Figure 7. It can be seen that the overshoot, undershoot, settling time, and speed deviations of all machines are greatly reduced by applying the proposed method in comparison with other techniques. Obtained results for this case study over the fourth objective function is presented in Table 4.

New England Power System
In this study, the 10-machine 39-bus power system shown in Figure 8 is considered. This system consists of 10 generators and each generator is assumed to be provided with governors, AVR, and IEEE ST1A-type static exciter. To   • When the system response is slow, the derivative effect of the PID-type fuzzy logic controller must decrease. • When the error is small and the system response is fast, the derivative effect of the PID-type FLC must increase.
The rule base that has been given in Table 1 is proposed using the above meta-rules. The inputs of the fuzzy rule and three operating conditions, which demonstrate the worst three values of damping factors, damping ratios, and damping scales have been presented in Table 6. This table presents that all damping factors of F 1 , F 3 , and F 4 are smaller than or equal to the value of σ 0 , and all damping ratios of F 2 and F 3 are more than or equal to ζ 0 . The damping scales of F 4 are more than χ. These results are also shown in Figure 11. The results show that F 1 , F 3 , and F 4 not only relocated unstable or lightly damped modes, but also shifted other modes to the left side in the plane.
In proposed controller model, the meta-rules for determining δ k can be summarized as follows:   The parameters that will be kept constant for all of the controllers used throughout the simulations have been selected as α = 0.2, β = 1, K e = 0.8, K d = 0.25.
base are the absolute values of error and the normalized acceleration.
The symmetrical rule base given in Table 1 is shown in Table 7. This rule base is used as the main fuzzy logic controller that is common to the controller of all methods.   needs more time to provide the appropriate results. So, we leave this section as a future work in this paper.

Conclusion
In this research, a new strategy based on fuzzy controller is proposed to solve the power system stability problem. In this strategy the proposed controller has been optimized by new chaotic invasive weed optimization (CIWO) algorithm. Also, four types of objective function applied by considering the damping factor, damping ratio, and a combination of them to analysis the power system behavior in different conditions. By considering the mentioned objective functions, it has been obtained that the fourth function had good characteristics to test the power system constants, are searched using the proposed algorithm within the range [0 5].
The optimum values of these parameters are found as follows: φ 1 = 0.68; φ 2 = 0.18; ψ 1 = 3.57; ψ 2 = 0.73. Actually, in this paper we have implemented our proposed method over power system by simulation through MATLAB software. All run times reported in this table are measured on the simple hardware set of 64-bit computer with 8 GB of RAM and Intel Core i7-4500U CPU@ 2.4 GHz by means of the computer timer. By this software, we have tried to provide a real-world engineering problem based on its simulation in detail.
Moreover, implementation of proposed controller is hart in experimental environment which is so costly and