A Frequency Control Approach for Hybrid Power System Using Multi-Objective Optimization

Mohammed Elsayed Lotfy 1,2,*, Tomonobu Senjyu 2, Mohammed Abdel-Fattah Farahat 1, Amal Farouq Abdel-Gawad 1 and Atsuhi Yona 2 1 Department of Electrical Power and Machines, Zagazig University, Zagazig 44519, Egypt; farahat707@hotmail.com (M.A.-F.F.); amgawad2001@yahoo.com (A.F.A.-G.) 2 Department of Electrical and Electronics Engineering, University of the Ryukyus, Okinawa 903-0213, Japan; b985542@tec.u-ryukyu.ac.jp (T.S.); yona@tec.u-ryukyu.ac.jp (A.Y.) * Correspondence: mohamedabozed@zu.edu.eg


Introduction
In most remote and isolated areas, electric power is often supplied by diesel generators.However, diesel generators cause serious impacts on the environment as every liter of diesel releases about three kilograms of CO 2 [1].Furthermore, diesel is expensive because transportation to remote areas adds extra cost.Due to these environmental and economic influences of the diesel generator, interest in alternative cost-effective, sustainable and clean energy sources has grown significantly.Wind, solar, sea, biomass and geothermal powers are sustainable and clean sources.Wind and solar have attracted much attention nowadays and have become the most widely-utilized renewable energy sources in power systems.Furthermore, the fuel cell (FC) has the ability to be considered as one of the green power sources of the future [2].However, hybrid power systems, especially in isolated systems with renewable energy sources, such as wind turbine generators (WTG) and photovoltaic (PV), face some stability problems because the power supplied by these sources is not constant, diverges quickly and cannot be easily predicted [3].Therefore, these oscillations in the renewable power sources can produce instantaneous mismatch in the vital balance between generation and demand.Consequently, continuous variations in frequency and voltage levels usually appear, which negatively affect the electric power system stability [4].Therefore, a continuous control for the supplied power by these renewable sources is required to ensure robust performance of the hybrid power system.
Nowadays, energy storage systems (ESS) are integrated with the renewable sources to maintain the safe operation of the power system and balance the supply and demand sides.These serve as backup devices and store excess power when the generation is more than demand and release power to the system when the demand is more than generation.This action helps in maintaining a steady flow of power irrespective of the load and generation power levels' fluctuations.As a result, it guarantees acceptable levels of systems frequency deviations [5,6].The flywheel (FW) stores electric energy by way of kinetic energy.The advantages of FW are high stored energy density, high power exchange with the system, high conversion efficiency of 80%-90% and long-lived pollution-free design [7].Furthermore, battery energy storage system (BESS) is growing rapidly in ESS technologies and extends the scope and the future of economic applications in hybrid power systems.Moreover, there are upward trends to install controllable loads, such as electric vehicle (EV) in the customer side of isolated areas, to control system frequency, as discussed in [8,9].
Recently, frequency control techniques for hybrid power systems using fuzzy logic control (FLC) [10][11][12][13], µ synthesis [14], H ∞ and µ-synthesis [15], neuro-fuzzy control [16], FLC with the particle swarm optimization (PSO) algorithm [17], FLC with chaotic PSO [18], PSO with mixed H 2 /H ∞ control [19], the quasi-oppositional harmony search algorithm (QOHSA) [20] and sliding mode control (SMC) [21] have been implemented with encouraging results.In spite of this, the door is still opened to more techniques to improve system frequency and reduce supply error in the face of the fluctuations of renewable sources and random load disturbances; especially, to overcome the drawbacks of H ∞ and FLC control techniques.The weighting functions in H ∞ control design cannot be easily selected, which affect the process of design significantly.Furthermore, the order of the H ∞ controller depends on that of the plant.This leads to a complex structure, which is difficult to implement, especially for large-scale systems.On the other hand, deciding the type and number of input and output membership functions and their associated parameters in the case of FLC is a time-consuming and complex task, which may not lead to satisfying results.Besides that, robust estimation for the supply error of the hybrid power system is required urgently to face the lack of measurements for some states inside the system.Furthermore, the main drawback of the majority of these previous studies is that this can capture the dynamic characteristics of the power system only for a specific operating point and may not have adequate performance under various operating conditions because of using a small-signal simplified transfer functions' model.Hence, in this paper, a new control approach is held to overcome this problem even with using the same small-signal simplified transfer functions.Considering the supply error in the frequency domain and dividing the responsibility to compensate this error between all components of the hybrid power system depending on their speed of response to cover all of the frequency space of the supply error is the key issue to guarantee the ability of the proposed control scheme to withstand a wide range of operating conditions and to defeat the previous studies' problems.Therefore, the performance of the hybrid power system with the inclusion of the proposed control scheme is investigated under various operating conditions, such as a sudden increase/decrease of wind speed, solar radiation, load demand and, also, actual wind speed data, to confirm the effectiveness and robustness of the proposed approach for a wide range of operating conditions.
Therefore, this paper proposes a full-order observer-based frequency control scheme for a hybrid power system in isolated small areas.The system consists of a WTG, PV, aqua-electrolyzer (AE), FC, BESS, FW, diesel engine generator (DEG) and EV.AE is used to resolve water into oxygen and hydrogen.The produced hydrogen is then stored in a hydrogen tank and used as fuel for FC.The full-order observer is utilized for the supply error estimation of the power system.The estimated supply error is considered in the frequency domain.BESS and FW, which have short time constants, are used to suppress the high frequency component of supply error; while FC, EV and DEG, which have long time constants, are controlled to reduce the slowly varying component of supply error.Two PI controllers are applied in the proposed hybrid power system.The task of each one of these controllers is to control the system frequency and mitigate one component of supply error.The epsilon multi-objective genetic algorithm (ε-MOGA) is applied to optimize the controllers' parameters.Numerical simulations are created by the MATLAB ® (Release 2016a, The MathWorks, Inc., Natick, MA, USA) environment to validate the robustness of the proposed control technique.This paper is organized as follows.Section 2 describes the hybrid power system configuration.Section 3 discusses the proposed ε-MOGA-based PI controllers' design approach, including the full-order observer implementation with a brief introduction of ε-MOGA theory.Section 4 presents the simulation results for five case studies of the hybrid power system with their detailed analysis.Specific conclusions are then drawn in Section 5.

Hybrid Power System Configuration
The single-line diagram of the proposed hybrid power system in this study is shown in Figure 1.
The system consists of PV, AE, FC, BESS, FW, DEG, EV and four units of WTG.Natural gas or water can be resolved into hydrogen and oxygen using AE.Then, the generated hydrogen is compressed, stored in hydrogen tanks and transported to FC through pipelines.The detailed block diagram of the system is presented in Figure 2. P d , P pv , P WTG , P AE , P FC , P BESS , P FW , P EV , P L and P C are the output power of DEG, PV, WTG, AE, FC, BESS and FW, the consumed power of EV, the load demand and the combined power, respectively.

Wind Turbine Generator Model
The output power of WTG depends on the wind speed at that instant.The fraction of extracted wind power converted to mechanical power of the rotor, P w , can be calculated a follows: where ρ (=1.25 kg/m 3 ) is the air density, A = (πR 2 ) is the rotor swept area, R (=23.5 m) is the radius of the blades, C p is the aerodynamic power coefficient and V w is the wind speed.C p is a function of both tip speed ratio λ and blade pitch angle β. λ is defined as the ratio of the speed at the blade tip to the wind speed and can be expressed as: where ω (=3.14 rad/s) is the rotational speed of blades.Then, an approximate expression for C p as a function of λ and β can be given as [7]:

Photovoltaic Model
The output power of the PV system under study can be determined as: where η is the conversion efficiency of the PV array ranging from 9% to 12%, S is the measured area of the PV array, φ is the solar radiation and T a is the ambient temperature in degrees Celsius.The value of P pv depends on T a and φ only because S and η are constants.However, in this study, T a is kept constant at 25 • C and P pv is linearly varied with φ only [4].Furthermore, the battery state of charge (SOC) is considered to ensure its safe working.SOC is indicated as a remaining energy level (REL), which expresses dischargeable energy as the percentage of the battery's rated capacity.The REL is obtained by the integral of the BESS output power, and its range is with that of SOC [22].
For the accurate simulation of the dynamic performances of practical WTG, PV, DEG, AE, FC, BESS, FW and EV, one should employ high-order mathematical models with non-linearities.However, for large-scale power systems' simulations, simplified models are generally used to capture the dynamic characteristics of a specific operating point [7].However, here, with applying the proposed control approach, the performance of the hybrid power system is investigated under various operating conditions to confirm the robustness of the proposed control scheme to overcome this problem and to clarify its superiority against the previous studies.Therefore, all of the systems are represented by a first-order lag transfer function in this paper, as shown in Figure 2.

Power and Frequency Deviations
To ensure a robust operation of the hybrid power system, the total power generated must be controlled to meet the required load demand, as the output power of WTG and PV fluctuates with wind speed and solar radiation, respectively.This is controlled based on the supply error ∆P e , which is the difference between the net power generation P S and P L and can be calculated as follows: ∆P e = (P WTG + P PV + P d + P FC ± P BESS ± P FW − P EV − P L ) (5) The system frequency deviates due to the generated power variations.Therefore, The fluctuations of the system frequency ∆F can be expressed as: where K f is the system frequency characteristic constant of the hybrid power system.Practically, due to the time delay between power variation and frequency fluctuation, the transfer function for system frequency variations is modified as: where T, D and M are the frequency characteristic time constant, load damping constant and inertia constant, respectively [23].
The continuous time dynamic performance of the hybrid power system is modeled by a set of state-space differential equations as follows: where X, U and Y are the state, input and output vectors, respectively.A, B, C and D are real constant matrices of the appropriate dimensions associated with the above vectors.X is a 14th order vector whose its elements are clearly indicated in Figure 2; while, U and Y can be given as:

Full-Order Observer Implementation
To improve frequency control in the power system including renewable sources, estimation for ∆P e is necessary to face model uncertainties and the lack of measurements for some states inside the system.P d and P C are estimated using ∆F and the input signal of governor P 1 , which can be calculated from ∆F.The state-space equations have the following form: Y = C X + D P 1 (13) where: C = 1 0 0 0 (16) Four state variables are considered, representing ∆F, P d , output power of the governor and P C , respectively.Y is ∆F only.The main construction of the full-order observer is shown in Figure 3, which can be illustrated by the next equations: substituting from ( 19) in (18): The purpose of this observer is to estimate the true state X.There may exist some error in the estimate at the initial time, but it is wished to decrease with time.The estimation error can be expressed as: The error signal obeys the following differential equation: Therefore: Thus, the state equation for this estimation error is a homogeneous differential equation governed by the (A − LC) matrix.Therefore, the duty is to choose the gain matrix L in a good manner, so that the eigenvalues of (A − LC) are placed in the left-half of the complex plane, which ensures the system stability with decaying estimation error over time.In this study, the poles of the observer are chosen as (γ 1 = −46, γ 2 = −0.25,γ 3 = −9, γ 1 = −65) by a good estimate.The estimated supply error ∆ Pe is estimated by the estimated value of the output power of DEG Pd and the estimated value of combined power PC using the following equation: Then, ∆ Pe is used as the input signal for the PI controllers as discussed later in the next sections.

Epsilon Multi-Objective Genetic Algorithm (ε-MOGA) Theory
ε-MOGA is an elitist multi-objective evolutionary algorithm based on the concept of epsilon-dominance, which is used to control the content of the archive A(t) where the result of the optimization problem is stored.ε-MOGA obtains an ε-Pareto set, Θ * p (which is not unique), that converges toward the Pareto optimal set, Θ p , in a smart distributed manner around the Pareto front µ(Θ p ) with limited memory resources.Moreover, it adjusts the limits of the Pareto front dynamically and prevents the solutions belonging to the ends of the front from being lost [24,25].To reach this goal, the objective space is split into a fixed number of boxes n − box i .Therefore, for each dimension i [1....n], n − box i cells of ε i width are created [26]; where: This grid preserves the diversity of µ(Θ * p ) since each box can be occupied by only one solution.For a solution x solution space, box i (x) can be defined by: A solution x 1 with value µ(x 1 ) ε-dominates the solution x 2 with value µ(x 2 ), denoted by Hence, a set Θ * p ⊆ Θ p is ε-Pareto, if and only if: Therefore, ε-MOGA updates A(t) by saving only ε-dominant solutions that do not share the same box.When two mutually ε-dominant solutions compete, the solution that remains in A(t) is the one closer to the center of the box.Therefore, preventing solutions belonging to adjacent boxes and increasing the diversity of the solution can be achieved.The algorithm is composed of three populations [27]: the main population, P(t), which explores the search space D s during the iterations, and its population size is Nind p ; the auxiliary population G(t), and its size is Nind G , which must be an even number; the last one is the archive, A(t), which stores the solutions Θ * p , and its size is Nind A , which is variable, but bounded by: where The main steps of the proposed algorithm are as follows [28]: Step 1. Begin and create empty A(t).
Step 2. P(0) is initialized with Nind p individuals that have been randomly selected from D s .
Step 3. Calculate the function value of each individual in P(t).
Step 4. Check individuals in P(t) that might be included in A(t), as follows: (1) Non-dominated individuals in P(t) are detected, Θ ND .
(2) Pareto front limits µ max i and µ min i are calculated from µ(x),∀ x Θ ND .(3) Individuals in Θ ND are analyzed, and those that are not ε-dominated by individuals in A(t) are included in A(t).
Step 5. Create G(t) as follows: (1) Two individuals are randomly selected, x p from P(t) and x A from A(t).
(3) If u > P c/m (probability of crossing/mutation), x p and x A are crossed over by means of the extended linear recombination technique.(4) If u < P c/m , x A and x p are mutated using Gaussian distribution and then included in G(t).This procedure is repeated Nind p /2 times until G(t) is filled up.
Step 6. Calculate the function value of each individual in G(t).
Step 7. Checks, one by one, which individuals in G(t) must be included in A(t) on the basis of their location in the objective space.
Step 8. Update P(t) with individuals from G(t).Every individual x G from G(t) replaces an individual x p that is randomly selected from the individuals in P(t) that are dominated by x G .However, x G will not be included in P(t) if there is no individual in P(t) dominated by x G .Finally, individuals from A(t) compose the smart characterization of the Pareto front, Θ * p .

Epsilon Multi-Objective Genetic Algorithm (ε-MOGA)-Based PI Controllers' Scheme
Two PI controllers are designed in this section as indicated in Figure 4.The input of the first controller is the sum of the low-frequency component of e∆P e and e∆F multiplied by a constant factor K, where e∆P e is the difference between the estimated supply error ∆ Pe and the command value of supply error ∆P * e , which is always zero; while e∆F is the difference between frequency deviation ∆F and the command value of frequency deviation ∆F*, which is always zero.The reason for multiplying e∆F by factor K is the fact that ∆F is very small compared to ∆ Pe .The value of K is determined by trial and error, such that ∆F is small.Controller output is the control signal to adjust the output powers of DEG and FC and the charged power of EV, respectively.On the other hand, the input for the second controller is the sum of the high-frequency component of e∆P e and e∆F multiplied by K. Controller output is then fed as the control signal to modify the stored and released powers of BESS and FW.
After that, the parameters of the two PI controllers have been tuned using ε-MOGA to meet the proposed performance.Two objective functions are introduced for each controller using integral absolute error (IAE) criteria as follows: (30) Subject to: The typical range selected for each of K p and K i is [0, 100].The parameters of ε-MOGA algorithm were set to: State-space equations of the hybrid power system including the proposed PI controllers are utilized to calculate the values of J 1 and J 2 every iteration till composing the optimal Pareto front.After simulation running of the ε-MOGA algorithm using the MATLAB environment, the final Pareto front of ε-MOGA is achieved, which consists of nine points.Each point contains one value for each of K p1 , K i1 , K p2 and K i2 , respectively.the final optimal Pareto front obtained is shown in Figure 5.The performance of the small power system with these values of the PI controllers parameters associated with each point is then investigated individually using the SIMULINK ® (Release 2016a, The MathWorks, Inc., Natick, MA, USA) environment to get the best operating point that can achieve the proposed performance and minimize the interactions between the objective functions.Finally, among all points of the Pareto front, the best point that can ensure adequate performance of the small power system and at the same time minimize the conflict between the objective functions has PI controller parameters as follows: K p1 = 0.0174, K i1 = 1.3129,K p2 = 98.8895, and K i2 = 0.5371.The flowchart of the complete design process for the proposed ε-MOGA-based PI controllers is presented in Figure 6.The singular value plots of the control loops with/without the proposed control approach are shown in Figure 7 to validate their stability performance.It is clear from Figure 7a that the singular value plot of the frequency control loop (∆F* → ∆F) has a high gain in the low-frequency domain (below 1 rad/s).Therefore, the frequency deviations of the hybrid system will increase in the low-frequency domain.Furthermore, there is a high gain in the low-frequency domain (below 1 rad/s) clarified in the supply error control loop (∆P * e → ∆P e ) as presented in Figure 7b.Therefore, supply error fluctuations due to variations of output powers of WTG, PV and load demand may increase in this domain.

Results
Using the simulation model, controllers and control approach discussed in the proceeding section, time domain simulated responses of the proposed hybrid power system under various operating conditions are considered in this section for 1200-s intervals.The performance of the proposed control control approach is compared with the recent well-implemented proportional-integral-derivative (PID) control approach optimized with QOHSA presented in [20] under different combinations of power generation and energy storage subsystems to investigate the effectiveness and robustness of the proposed control scheme.Table 1 shows the parameters of the proposed power system [23] used in the simulation.The power rating of all components of the hybrid power system are as follows: 0 < P WTG < 0.83 pu, 0 < P pv < 0.52 pu, 0 < P d < 0.37 pu, 0 < P FC < 0.28 pu, −0.18 < P BESS < 0.18 pu, −0.38 < P FW < 0.38 pu, 0 < P EV < 0.0055 pu, 0.25 < SOC < 0.9 pu.The initial SOC is assumed to be 0.6.

Case 1
This case study is considered as the basic one consisting of all power generation sources and ESS included in the hybrid power system.Figure 8 presents the time domain simulation results of this case whose simulation responses under different operating conditions are discussed respectively in the following section.(1) Base case: Through 0 s < t < 200 s and 1000 s < t < 1200 s, the average solar radiation φ and wind speed V W are about 0.5 pu and 7.5 m/s, respectively.Depending on this, the average output powers of WTG units and the PV array are around 0.3 pu and 0.5 pu, respectively, with P L equal to 1 pu.Since P WTG with P PV are not sufficient yet to meet demand alone, DEG and FC are connected automatically to the system at t = 0 s.Furthermore, BESS and FW begin to release power to the system trying to decrease the power shortage in this time.No surplus power can be used to charge EV in this interval, as shown in Figure 8.The proposed control scheme succeeded to suppress ∆F and also ∆P e fluctuations with a good estimate of the full-order observer as indicated in the figure, which illustrates ∆P e and estimated ∆P e .However, for the conventional approach, DEG and FC produce a larger amount of power compared to that of the proposed scheme with high fluctuations and overshoot, as shown in the figure.Therefore, BESS and FW decrease their released power in this interval, but still with high fluctuations.Furthermore, the charged power of EV fluctuates significantly.In addition, the conventional approach has clear oscillations for ∆F and ∆P e compared to the proposed one.(2) Sudden load drop: At 200 s < t < 300 s, φ and V W have the same values of the base case, but P L drops suddenly from 1 pu to 0.75 pu at t = 200 s.Consequently, P WTG and P PV produce sufficient output power to meet P L in this case.Therefore, DEG and FC are automatically disconnected from the system at t = 200 s.EV starts to charge power from the system.Furthermore, BESS and FW begin to store excess power.In spite of the sudden load decrease, the proposed control approach still ensures adequate damping performance for ∆F and ∆P e deviations with a suitable estimate of the full-order observer.On the other hand, for the conventional control approach, EV charges to its capacity limit during this interval.Large fluctuations of ∆F and ∆P e still appear compared to the proposed approach.(3) Sudden load rise: During 300 s < t < 400 s, φ and V W still have the same values of the base case, but P L rises suddenly from 0.75 pu to 1 pu at t = 300 s.The hybrid power system returns to the same operating point of the base case.Therefore, all responses are the same as those in the interval 0 s < t < 200 s.The control technique still has the ability to damp ∆F and ∆P e fluctuations with robust performance of the observer, while the conventional control approach has a clear amount of fluctuations for ∆F and ∆P e .(4) Sudden wind speed rise: When 400 s < t < 500 s, the values of φ and P L are still 0.5 pu and 1 pu, respectively.However, V W suddenly increases from 7.5 m/s to 10 m/s at t = 400 s.Therefore, P WTG increases to around 0.8 pu, which is enough with P PV to meet load demand in this period.Then, DEG and FC are totally disconnected from the system at t = 400 s.Furthermore, BESS, FW and EV start to store excess power from the system with a larger amount as compared to what happened in the interval 200 s < t < 300 s, as a result of more existing surplus power in the system due to the wind speed increase.Proposed controllers still mitigate ∆F and ∆P e fluctuations with effective performance of the full-order observer.For the conventional approach, oscillations of ∆F and ∆P e can be noticed clearly with a higher amplitude than that of the proposed control scheme.(5) Sudden wind speed drop: At 500 s < t < 700 s, V W suddenly decreases from 10 m/s to 7.5 m/s at t = 500 s.Therefore, P WTG returns to around 0.3 pu with the same values of φ and P L .The proposed power system comes back to the same operating condition of the base case.The simulation results are the same ones in the interval 0 s < t < 200 s.∆F and ∆P e deviations are still in the acceptable limits, as shown in Figure 8, during this interval with the proper estimate of the full-order observer.However, for the conventional approach, the speed of response of DEG and FC was very slow to the sudden decrease of wind speed, while BESS and FW increase their released power to overcome this shortage.There are still large fluctuations of ∆F and ∆P e compared to the proposed approach.(6) Solar radiation sudden drop: During 700 s < t < 800 s, the values of V W and P L are still 7.5 m/s and 1 pu, respectively.However, φ drops to 0.2 pu suddenly.Therefore, P PV also drops to 0.2 pu.Accordingly, DEG and FC are connected automatically to the system producing a larger amount of output power compared with previous intervals to compensate the drop in P PV .Furthermore, BESS and FW begin to release power to the system.There is no surplus power to charge EV in this period.The proposed controllers damp all fluctuations of ∆F and ∆P e with the perfect estimate of the full-order observer.On the other hand, for the conventional approach, DEG and FC produce output more power than that of the proposed approach.Therefore, BESS and FW decrease their released power during this interval.Large fluctuations appear for the charged power of EV,∆F and ∆P e compared to the proposed approach.(7) Solar radiation linear increase: When 800 s < t < 1000 s, φ linearly increases in this interval from 0.2 pu to 0.5 pu, while V W and P L have constant values of 7.5 m/s and 1 pu, respectively.Hence, P PV linearly increases also from 0.2 pu to 0.5 pu.Then, P DEG , P FC and the released powers of BESS and FW linearly and slightly decrease to meet the power balance condition.There is no excess power to feed EV also in this period.The proposed control scheme still suppress ∆F and ∆P e fluctuations, and the full-order observer has a good performance in this interval.
For the conventional approach, no response for BESS and FW happened in this period, and large fluctuations still exist for ∆F, ∆P e and the charged power of EV.

Case 2
Figure 9 illustrates the simulation results of this case study, which consists of all components of the hybrid power system, except WTG units.The simulation response under various operating conditions can be examined, respectively, as follows.(1) Base case: Through 0 s < t < 200 s and 1000 s < t < 1200 s, the values of φ and P L are 0.5 pu and 1 pu, respectively.Therefore, the average output power of the PV array is about 0.5 pu.Due to the WTG units' outage, P PV is not enough to meet load demand.Therefore, DEG and FC are connected automatically to the system at t = 0 s, generating higher output power compared to the same interval in Case 1 to compensate the lack of generated power due to the outage of WTG units.Furthermore, BESS and FW start to release power to the system to face the shortage of generated power in this period.There is no excess power that can be utilized to charge EV, as shown in Figure 9.The proposed control approach damps all of the ∆F and ∆P e fluctuations with a good estimate of the full-order observer.On the other hand, for the conventional control approach, DEG and FC produce a larger amount of power with high fluctuations and overshoot compared to that of the proposed scheme, as shown in the figure.Therefore, BESS and FW decrease their released power in this interval.Furthermore, high fluctuations appear clearly for the charged power of EV, ∆F and ∆P e .(2) Sudden load drop: At 200 s < t < 300 s, the value of φ is still 0.5 pu.Therefore, P PV has the same previous value of 0.5 pu, but P L decreases suddenly from 1 pu to 0.75 pu at t = 200 s.As a result, DEG and FC decrease their generated output power as the mismatch between generation and demand powers is reduced due to the sudden load drop.Furthermore, BESS and FW begin to decrease their released power to the system.No surplus power exists even in this period to charge EV.The proposed control technique still has the ability to mitigate ∆F and ∆P e deviations with a suitable estimate of the full-order observer, as shown in Figure 9.However, for the conventional control approach, large fluctuations still exist for ∆F, ∆P e and the charged power of EV compared to that of the proposed control scheme.(3) Sudden load rise: During 300 s < t < 700 s, P L increases suddenly from 0.75 pu to 1 pu at t = 300 s, and the value of φ is still 0.5 pu.Therefore, the power system comes back to the same operating point of the base case.All simulation responses are the same as those in the interval 0 s < t < 200 s.The control scheme succeeded to damp the ∆F and ∆P e fluctuations with the effective estimate of the observer, while the conventional control technique fails to damp the oscillations of charged power of EV, ∆F and ∆P e , as shown in Figure 9. (4) Solar radiation sudden drop: During 700 s < t < 800 s, P L is still 1 pu, but φ suddenly drops to 0.2 pu at t = 700 s.Therefore, P PV also decreases to 0.2 pu.Consequently, DEG and FC increase their output power to compensate the sudden drop of P PV and with values more than that in the same period of Case 1 due to the WTG units' outage.Furthermore, BESS and FW release more power to meet the load demand in this interval.Still, there is no excess power to charge EV in this period.The proposed control approach still ensures keeping the ∆F and ∆P e deviations within the acceptable limits with the robust estimate of the full-order observer.For the conventional approach, fluctuations of ∆F and ∆P e still appear with a large value compared to that of the proposed technique.( 5) Solar radiation linear increase: When 800 s < t < 1000 s, P L still has a constant value of 1 pu.
However, φ linearly increases from 0.2 pu to 0.5 pu in this period.Therefore, P PV also linearly rises from 0.2 pu to 0.5 pu.Accordingly, P DEG , P FC and the discharged power of BESS and FW linearly decrease to achieve balance between generation and load demand in this interval.
There is no surplus power to charge EV, as shown in Figure 9.The control mechanism is still able to suppress the ∆F and ∆P e fluctuations, and the full-order observer has adequate performance during this interval, as indicated in Figure 9.However, for the conventional approach, there is no response of DEG and FC due to the linear increase in solar radiation, and they still produce their maximum power.Figure 9 shows large fluctuations of ∆F and ∆P e compared to that of the proposed control scheme in this interval.

Case 3
In this case, the hybrid power system contains all of the power generation sources and ESS with complete shading of solar radiation for 300 s.The time domain simulation responses are indicated in Figure 10, which can be clarified in the following section.(1) Base case: When 0 s < t < 200 s and 1000 s < t < 1200 s, the values of V W , φ and P L are 7.5 m/s, 0.5 pu and 1 pu, respectively.Hence, the averages output powers of WTG units and the PV array are around 0.3 pu and 0.5 pu, respectively.DEG and FC are connected automatically to the system, supplying their output power due to the fact that P PV and P WTG cannot meet the load demand alone in this interval.Furthermore, BESS and FW start to release their power to achieve the power balance condition.There is no excess power to charge EV in this period.The control approach is capable of damping ∆F and ∆P e fluctuations with a good estimate of the full-order observer, as shown in Figure 10.However, for the conventional control scheme, DEG and FC produce higher power with large fluctuations compared to that of the proposed one, as shown in the figure.As a result, BESS and FW decrease their released power in this interval.Figure 10 indicates the clear fluctuations of the charged power of EV, ∆F and ∆P e in this period.(2) Sudden load drop: During 200 s < t < 300 s, P L suddenly decreases from 1 pu to 0.75 pu at t = 200 s, while V W and φ have constant values of 7.5 m/s and 0.5 pu, respectively.Therefore, P WTG and P PV still also have the same values of the base case.Now, DEG and FC are automatically disconnected from the system because P WTG and P PV can produce sufficient power to meet the load demand in this period.Furthermore, BESS and FW begin to store the excess power from the system.In the same time, EV starts to charge part of the surplus power.The proposed control scheme still has the ability to keep the ∆F and ∆P e deviations around zero with robust performance of the observer, while the conventional control approach still has high fluctuations for ∆F and ∆P e , as shown in Figure 10.(3) Sudden load rise: At 300 s < t < 400 s, P L suddenly rises from 0.75 pu to 1 pu at t = 300 s, while V W and φ still have their same previous values.Therefore, the hybrid power system came back safely to the operating point of the base case.All time domain simulation results are the same as those in the interval 0 s < t < 200 s.In spite of these sudden changes, the control technique still holds ∆F and ∆P e fluctuations at about zero.Furthermore, the full-order observer guarantees an effective estimate in this interval.On the other hand, Fluctuations of the charged power of EV, ∆F and ∆P e still exist using the conventional control approach.(4) Sudden wind speed rise: During 400 s < t < 500 s, V W suddenly increases from 7.5 m/s to 10 m/s at t = 400 s.Hence, P WTG rises to around 0.8 pu.The values of φ and P L are still 0.5 pu and 1 pu, respectively.Now, P WTG and P PV are sufficient to meet the load demand in this period.Therefore, DEG and FC are totally disconnected from the system at t = 400 s.All of BESS, FW and EV start to store surplus power from the system with higher values compared to those in the period 200 s < t < 300 s due to more excess power existing in the hybrid system in this interval, as a result of wind speed increase.The control mechanism still succeeds to suppress the ∆F and ∆P e deviations with the adequate estimate of the full-order observer.For the conventional approach, EV charges to its capacity limit in this interval with clear fluctuations of ∆F and ∆P e .(5) Sudden wind speed drop with solar radiation complete shading: At 500 s < t < 800 s, V W suddenly decreases from 10 m/s to 7.5 m/s at t = 500 s.Furthermore, complete shading of solar radiation occurs at this interval.Therefore, φ drops to 0 pu.As a result, P WTG and P PV decrease to about 0.3 pu and 0 pu, respectively while P L is still 1 pu.Therefore, DEG and FC are automatically connected to the system to meet the load demand and withstand this severe operating condition.Furthermore, BESS and FW start to release power to the system to achieve the power balance condition.There is no excess power in this period to charge EV.Simulation results indicate that the proposed control approach still damps the ∆F and ∆P e fluctuations with robust performance of the full-order observer even in such a harsh operating condition, while the conventional approach cannot withstand this large change of power balance, and large fluctuations appear in the beginning of this interval in the output power of DEG, released powers of BESS and FC, the charged power of EV, ∆F and ∆P e , as shown in Figure 10.Furthermore, the speed response of FC due to this change of power balance is low compared to that with the proposed control scheme.(6) Solar radiation linear increase: When 800 s < t < 1000 s, the values of V W and P L are still 7.5 m/s and 1 pu, respectively.φ linearly increases from 0 pu to 0.5 pu in this interval.Therefore, P PV also increases linearly from 0 pu to 0.5 pu.Hence, P DEG , P FC and the released powers of BESS and FW linearly decrease due to the P PV increase to meet the power balance between generation and load demand.No surplus power exists for charging EV in this period.The proposed control technique is still able to mitigate ∆F and ∆P e deviations with an efficient estimate of the full-order observer.However, the conventional control approach still fails to damp the fluctuations of the charged power of EV, ∆F and ∆P e .

Case 4
In this case, WTG units and PV arrays are excluded from the hybrid power system.The time domain simulation results are shown in Figure 11, and the dynamic performance of all components of the power system are discussed in the following section.(1) Base case: At 0 s < t < 200 s, the value of P L is 1 pu.Therefore, DEG and FC start to feed the power system by their output powers with more amounts than those in the same period of all of the previous cases to compensate the WTG units and PV array outages.Furthermore, BESS and FW begin to release power to the hybrid power system trying to achieve the power balance condition in this severe operating condition.No excess power exists to charge EV in this interval.Simulation results shown in Figure 11 confirm the ability of the proposed control scheme to damp the ∆F and ∆P e fluctuations with a robust estimate of the full-order observer.On the other hand, for the conventional control approach, DEG and FC produce their maximum limit of power in this period.Therefore, BESS and FW release less power compared to that of the proposed control scheme.Furthermore, there is no excess power to charge EV.However, large fluctuations of ∆F and ∆P e still appear compared to that of the proposed approach.(2) Sudden load drop: During 200 s < t < 300 s, P L suddenly drops from 1 pu to 0.75 pu at t = 200 s.
Hence, DEG and FC begin to decrease their output power as the difference between generation and load demand reduces.Furthermore, BESS and FW start to decrease their released power to the system.Though load reduction occurred, there is no surplus power to charge EV in this interval.
The control technique keeps ∆F and ∆P e fluctuations within the acceptable limits.Furthermore, the full-order observer guarantees an effective estimate for ∆P e , as shown in Figure 11.For the conventional control approach, there is no response for DEG and FC due to the sudden load drop, and they keep producing their maximum limit of power.As a result, BESS and FW continue to decrease their released power in this interval.Furthermore, there is no surplus power to charge EV. Figure 11 shows clearly the fluctuations of ∆F and ∆P e with using the conventional control scheme in this period.(3) Sudden load rise: At 300 s < t < 1000 s, P L suddenly increase from 0.75 pu to 1 pu at t = 300 s.
Therefore, the hybrid power system returns safely to the operating condition of the base case.All simulation responses are the same as those in the period 0 s < t < 200 s.Despite these sudden changes and severe operating conditions, the control approach is still capable of holding ∆F and ∆P e deviations around zero with efficient performance of the observer as indicated in Figure 11.However, the conventional control approach still fails to damp fluctuations of ∆F and ∆P e in this interval, as shown in Figure 11.

Case 5
This case study demonstrates the time domain simulation results of the hybrid power system the same as Case 1, substituting real wind speed data of Ginoza wind farm, Okinawa, Japan.Only 12 h of data from 6 a.m. to 6 p.m. for one day in January are utilized.The simulation responses of all components of the power system can be analyzed in the following section.
(1) At 0 h < t < 5 h, V W starts near 6 m/s and increases till reaching around 8 m/s at t = 5 h.Therefore, P WTG begins to rise also from about 0.2 pu to 0.5 pu in the same period.However, P WTG with P PV is not sufficient yet to meet the load demand.Hence, DEG and FC are automatically connected to feed the power system with decreasing output power as P WTG increases.Furthermore, BESS and FW start to release their power to the system trying to achieve the power balance condition.
No surplus power exists to charge EV in this interval.Simulation results clarify that the proposed control scheme can keep ∆F and ∆P e fluctuations near zero with the effective performance of the proposed full-order observer.For the conventional control approach, DEG and FC produce a larger amount of power with higher overshoot and fluctuations compared to that of the proposed approach.Therefore, BESS and FW decrease their released power in this interval.Clear fluctuations for the charged power of EV, ∆F and ∆P e compared to the that of the proposed scheme can be noticed in Figure 12 in this period.(2) When 5 h < t < 9 h, V W begins to exceed 8 m/s at t = 5 h.At this time, P WTG with P PV begins to be enough to meet the load demand.Therefore, DEG and FC are totally disconnected from the system.Furthermore, BESS and FW start to charge the excess power.At the same time, EV begins to store part of the surplus power from the hybrid power system.The proposed control approach still has the ability to damp the ∆F and ∆P e fluctuations as indicated in Figure 12.Furthermore, the full-order observer ensures robust estimate performance during this period, while higher fluctuations of ∆F and ∆P e still appear in Figure 12 for using the conventional control approach.(3) During 9 h < t < 12 h, V W begins to decrease.Therefore, P WTG also starts to drop, but still is sufficient with P PV to feed the load demand.Hence, DEG and FC are still totally disconnected from the power system.BESS, FW and EV begin to reduce their stored power from the hybrid power system as P WTG decreases.The proposed control technique succeeded to mitigate ∆F and ∆P e deviations within the acceptable limits with adequate performance of the full-order observer, as shown in Figure 12.However, the conventional control scheme is still unable to damp the fluctuations of ∆F and ∆P e compared to the proposed control approach.
It can be clearly seen from Figure 12 that there are two dips on P d , P FC , P BESS and P FW due to sudden load and solar radiation decreases at the first hour.However, these sudden changes do not affect the proposed control scheme, which is still capable of keeping the intended performance for damping all ∆F and ∆P e fluctuations around zero.Furthermore, the proposed full-order observer guarantees robust achievement for ∆P e estimation, even in such severe operating conditions.

Conclusions
new frequency control scheme for a hybrid power system is presented in this paper.The proposed power system consists of a PV array, AE, FC, BESS, FW, DEG, EV, and four units of WTG.The full-order observer is utilized to estimate the power system supply error.Hence, the estimated supply error is considered in the frequency domain.The high frequency component of supply error is suppressed by BESS and FW, which have short time constants, while DEG, FC and EV, which have long time constants, are controlled to mitigate the low frequency component of supply error.Two PI controllers are applied in the hybrid power system.The mission of each one is to control the system frequency and damp one component of supply error.ε-MOGA is used to optimize the controllers' parameters.The performance of the proposed control approach is compared with that of recent well-implemented approaches, such as the PID controller optimized by QOHSA.Five case studies are considered and analyzed to investigate the robustness and effectiveness of the proposed control scheme under various operating conditions.Simulation results confirm the superiority of the proposed control technique to: (1) Control the flow of the output powers of DEG and FC, the stored and released powers of BESS and FW and the charged power of EV at any time according to the variations of wind speed, solar radiation and load demand, so as to suppress the system frequency and supply error fluctuations.(2) Withstand harsh operating conditions, such as complete shading of the PV array, fast wind speed variation and WTG units and/or PV array outage keeping the system frequency and supply error values within acceptable limits.(3) Decrease the fluctuations of the energy supply and storage systems (DEG, FC, FW, BESS and EV).
Therefore, if this control approach is used, a lesser size of these energy sources will be needed, which improves the overall system efficiency and decreases its cost.
Overall, the performance of the hybrid power system is enhanced significantly using the proposed control approach.Extending this study to analyze the performance of the hybrid power system including a detailed non-linear model for all components is the subsequent task in the near future to make the related study results further solid and practical.

Figure 1 .
Figure 1.Single-line diagram of the hybrid power system.PI: proportional-integral controller.

Figure 7 .
Figure 7. (a) Singular value plot of frequency control loop (∆F* → ∆F); and (b) singular value plot of the supply error control loop (∆P * e → ∆P e ).