Multiobjective Genetic Algorithm-Based Optimization of PID Controller Parameters for Fuel Cell Voltage and Fuel Utilization

Nowadays, given the great deal of fossil fuel consumption and associated environmental pollution, solid oxide fuel cells (SOFCs) have shown their great merits in terms of high energy conversion efficiency and low emissions as a stationary power source. To ensure power quality and efficiency, both the output voltage and fuel utilization of an SOFC should be tightly controlled. However, these two control objectives usually conflict with each other, making the controller design of an SOFC quite challenging and sophisticated. To this end, a multi-objective genetic algorithm (MOGA) was employed to tune the proportional–integral–derivative (PID) controller parameters through the following steps: (1) Identifying the SOFC system through a least squares method; (2) designing the control based on a relative gain array (RGA) analysis; and (3) applying the MOGA to a simulation to search for a set of optimal solutions. By comparing the control performance of the Pareto solutions, satisfactory control parameters were determined. The simulation results demonstrated that the proposed method could reduce the impact of disturbances and regulate output voltage and fuel utilization simultaneously (with strong robustness).


Introduction
Nowadays, with a great deal of fossil fuel consumption, many environmental issues have been raised.With the merits of low pollution and high efficiency [1], hydrogen energy-based fuel cells are considered to be ideal alternative energy sources, and they are widely used to supply power to stand-alone power systems and vehicles [2][3][4].Among all fuel cells, solid oxide fuel cells (SOFCs) attract much attention for their multifuel capabilities and potential for cogeneration [5], as well as their ease of coordination with refinery hydrogen networks [6].Numerous studies about SOFCs have been done.Rao et al. [7] analyzed the thermodynamic characteristics and techno-economic performance of a system consisting of an SOFC-gas turbine system, Xia et al. [8] studied low-temperature operations of SOFCs, and Wang et al. [9] modified the anode material to solve the problem of coke formation.Besides these, studies on efficient control strategies to ensure the safety, stability, and efficiency of SOFCs are also very important.In order to maintain the efficiency of an SOFC, fuel utilization in an SOFC must be kept between 70% and 90% [10].In addition, to ensure the power quality of an SOFC, the output voltage must be controlled tightly.However, due to the nonlinearity and couplings of multiple variables and frequent load disturbances, the control of an SOFC is extremely hard, and the controller must be appropriately designed.As a result, designing an efficient control scheme is of great significance in the research on SOFCs.
To solve the problem of maintaining the output voltage and fuel utilization of SOFCs, many advanced control strategies have been proposed, such as model predictive control (MPC) [11][12][13], fuzzy proportional-integral-derivative (PID) control [14], fuzzy logic control [15], and neural network control [16].All of these control strategies have been shown to obtain excellent control performance in numerous simulation studies: However, due to their computational complexity, these control strategies are difficult to implement in practical applications.
With the merits of a simple structure and high reliability, PID controllers are widely used in the process industry: However, to ensure control performance, the control parameters of PID controllers must be tuned well.To this end, numerous artificial intelligence (AI)-based optimization methods have been applied to tune PID controllers.For example, Qin, Sun, and Hua [17] used a fruit fly optimization algorithm (FOA) to optimize the PID controllers of a refrigeration system, Zamani et al. [18] applied a particle swarm optimization (PSO) algorithm to optimize the PID controller of an automatic voltage regulator (AVR), and Elbayomy et al. [19] optimized the PID controller of an electrohydraulic servo actuator system (EHSAS) to control the movable surface of a space vehicle.However, these references all used a single-objective optimization method, where the authors judged a controller's performance through one objective function that weighted the sum of different performance indices, such as the integral absolute error (IAE), overshoot, and settling time.Although in the end this method achieved good results, it ignored the problem that the control objectives were usually conflicting and that the weight factors were decided subjectively by the author.Instead of using only one objective function, a Pareto-based multiobjective optimization algorithm involves two or more objective functions and is thus able to deal with the conflicting objectives of a control system without designating weight factors: Many researchers choose this method to solve optimization problems [20][21][22][23].
This paper aims to (1) identify a dynamic SOFC model as a linearized model, (2) apply two PID controllers to maintain the output voltage and fuel utilization of an SOFC under load disturbance, and (3) enhance the control system performance by optimizing the controller with multiobjective optimization.The remaining part of this paper is organized as follows: In Section 2, a dynamic SOFC model is introduced, while in Section 3, a linearized model is identified and a controller for the model is designed.Section 4 optimizes the controller parameters with genetic algorithm (GA)-based multiobjective optimization, and Section 5 applies the optimal controllers to the nonlinear system and analyzes the optimization results.Conclusions are drawn in Section 6.

Dynamic Model for an SOFC
The model presented in this paper was based on three submodels: An electrochemical model, a mass balance model, and an energy balance model.To simplify the model, six assumptions were made as follows: 1.
Both hydrogen and air are preheated to a specific temperature before entering the anode and cathode; 2.
Both hydrogen and air are in full contact with the anode and cathode, and the stoichiometric quantity of oxygen at the cathode is sufficiently large; 3.
A large stoichiometric quantity of oxygen exists at the cathode; 4.
The oxygen concentration in the air is 21%; 5.
The whole stack can be presented as a combination of individual stacks; 6.
Ideal gas laws are employed for both the fuel flow and air flow.
A structural diagram of this model is shown in Figure 1, where the fuel flow, air flow, oxygen ion flow, and electron flow are represented by red, blue, yellow, and green arrows, respectively.Sustainability 2019, 10, x FOR PEER REVIEW 3 of 20  4 2

Electrochemical Model
As we can tell from Figure 1, the SOFC consists of an anode and a cathode, and they are connected by a solid ceramic electrolyte.Through the electrolyte, oxygen ions are transported from the cathode to the anode, and electrons are transported from the anode to the cathode through the external electric circuit.The reaction on the anode and the cathode and the total reaction are shown as follows: 2H O Based on Nernst's equation [14], the open circuit voltage of the fuel cell stack Ecell can be formulized as Here, E0,cell can be calculated as where 0 0,cell E is the reference voltage at 298 K and 1 atm of pressure.
The partial pressures in Equation ( 4), pH2, pO2, and pH2O, are the average of the inlet and outlet partial pressures of the respective components, calculated as follows: Due to ohmic loss, activation loss, and concentration loss, the output voltage of the fuel cell stack Vcell is lower than Ecell, which is calculated as The expression of the three voltage drops are given below:

Electrochemical Model
As we can tell from Figure 1, the SOFC consists of an anode and a cathode, and they are connected by a solid ceramic electrolyte.Through the electrolyte, oxygen ions are transported from the cathode to the anode, and electrons are transported from the anode to the cathode through the external electric circuit.The reaction on the anode and the cathode and the total reaction are shown as follows: Based on Nernst's equation [14], the open circuit voltage of the fuel cell stack E cell can be formulized as Here, E 0,cell can be calculated as where E 0 0,cell is the reference voltage at 298 K and 1 atm of pressure.The partial pressures in Equation (4), p H2 , p O2 , and p H2O , are the average of the inlet and outlet partial pressures of the respective components, calculated as follows: Due to ohmic loss, activation loss, and concentration loss, the output voltage of the fuel cell stack V cell is lower than E cell , which is calculated as The expression of the three voltage drops are given below: where R ohm is the internal resistance of the fuel cell, and it is calculated as follows [24]: According to assumption 5, the output voltage of the fuel cell can be calculated with V cell by where N is the number of fuel cells in a series.

Mass Balance Model
According to assumption 5, the dynamic effective partial pressure of H 2 , O 2 , and H 2 O can be presented as follows [25]: V an RT dp ch while based on assumption 2, the inlet and outlet mass flow rate of H 2 , O 2 , and H 2 O can be described as

Energy Balance Model
As Figure 1 shows, hydrogen flows through the anode surface and changes into water, while air is provided by a central air supply tube (AST) and is forced to flow through the cathode.Due to the simplification of the model, the means of heat transfer inside the cell are radiation, convection, and mass flow [26].
Therefore, the energy balance equations for different parts of the cell are shown as follows: 1.
Cell tube: q cell,in = q gen = q chem − q ele , ( q cell,out = q rad + q conv,ann + q f low,air,ann + q conv, f uel + q f low, f uel , ( q conv,ann = h cell A cell,inner (T cell − T air,ann ), q f low,air,ann = M air x air C air (T air,ann,in − T air,ann,out ), Fuel: Air between the cell tube and the AST: q air,ann,net = m air,ann C air dT air,ann dt ; AST: q AST,conv,outer = h AST,outer A AST,outer T air,cell − T AST , q AST,conv,inner = h AST,inner A AST,inner (T AST − T air,AST ), q f low,air,AST = m air C air (T air,AST,in − T air,AST,out ); Air in the AST: Based on the formulation of ( 21)-( 35), the energy balance model for the SOFC was developed.

Control Problems
The main propose of the controller design was to ensure that the output voltage V out and the fuel utilization FU follow their reference values and maintain their stability as strongly as possible.The main problems in the controller design are listed as follows: 1.
High coupling: The SOFC is a coupled system, which makes its controller design challenging and complicated; 2.
Frequent disturbance: Due to the frequent change of the load, a continual disturbance exists in the system, which requires the controller to have strong robustness and accurate control performance; 3.
Conflicting objectives: The accuracy of the control performance of V out and FU are the two conflicting objectives.However, to improve the power quality of the SOFC, both objectives need to be optimized at the same time.

Model Identification
Since the SOFC model is a multi-input-multi-output (MIMO) system, the load resistance R was chosen as the disturbance of the system, the fuel flow m H2 and the cathode pressure P c were chosen as the control variables, and the output voltage V out and fuel utilization FU were chosen as the controlled variables.Thus, the linearized model can be identified as a 3 × 2 MIMO system: where y 1 is the output voltage V out , y 2 is the fuel utilization FU, u 1 is the fuel flow m H2 , u 2 is the cathode pressure P c , d is the load resistance R, G 11 represents the transfer function from m H2 to V out , G 12 represents the transfer function from m H2 to FU, G 21 represents the transfer function from P c to V out , G 22 represents the transfer function from P c to FU, G d1 represents the transfer function from R to V out , and G d2 represents the transfer function from R to FU.
To identify these transfer functions, a step response of this system was analyzed.First, a step signal was set on m H2 , while P c and R were kept as constants, and the step responses of V out and FU were recorded.Second, a step signal was set on P c , while m H2 and R were kept as constants, and then the step responses of V out and FU were recorded.Third, a step signal was set on R, while P c and m H2 were kept as constants, and the step responses of V out and FU were recorded.Finally, based on the step response results, the transfer functions can be readily identified as follows: The step response curves are drawn in Figure 2. As we can tell from the results, the step response curves of the identified model are the same as those of the original model, and thus the identified model is able to be used for control design.
as the control variables, and the output voltage Vout and fuel utilization FU were chosen as the controlled variables.Thus, the linearized model can be identified as a 3 × 2 MIMO system: where y1 is the output voltage Vout, y2 is the fuel utilization FU, u1 is the fuel flow mH2, u2 is the cathode pressure Pc, d is the load resistance R, G11 represents the transfer function from mH2 to Vout, G12 represents the transfer function from mH2 to FU, G21 represents the transfer function from Pc to Vout, G22 represents the transfer function from Pc to FU, Gd1 represents the transfer function from R to Vout, and Gd2 represents the transfer function from R to FU.
To identify these transfer functions, a step response of this system was analyzed.First, a step signal was set on mH2, while Pc and R were kept as constants, and the step responses of Vout and FU were recorded.Second, a step signal was set on Pc, while mH2 and R were kept as constants, and then the step responses of Vout and FU were recorded.Third, a step signal was set on R, while Pc and mH2 were kept as constants, and the step responses of Vout and FU were recorded.Finally, based on the step response results, the transfer functions can be readily identified as follows: The step response curves are drawn in Figure 2. As we can tell from the results, the step response curves of the identified model are the same as those of the original model, and thus the identified model is able to be used for control design.

Relative Gain Array (RGA) Paring
To analyze the relationship between the input and output variables of a system, we applied the relative gain array (RGA) defined by Skogestad and Postlethwaite in [27].First, by setting the time s to 0, we got the steady-state matrix of G11, G12, G21, and G22: Then, the RGA matrix was calculated: ( ) The results show that the output voltage Vout is affected more by the cathode pressure Pc, while the fuel utilization FU is affected more by the fuel flow mH2.Thus, Vout can be controlled by Pc easily, and FU can be controlled by mH2 easily.

Controller Design
In this paper, two PID controllers were applied to control mH2 and Pc independently.The formula for a PID controller is shown as ( ) where kp, ki, and kd are the proportional gain, integral gain, and derivative gain, respectively.In previous research, the main work has been to search for the optimal values of kp, ki, and kd [17,19,28].

Relative Gain Array (RGA) Paring
To analyze the relationship between the input and output variables of a system, we applied the relative gain array (RGA) defined by Skogestad and Postlethwaite in [27].First, by setting the time s to 0, we got the steady-state matrix of G 11 , G 12 , G 21 , and G 22 : Then, the RGA matrix was calculated: The results show that the output voltage V out is affected more by the cathode pressure P c , while the fuel utilization FU is affected more by the fuel flow m H2 .Thus, V out can be controlled by P c easily, and FU can be controlled by m H2 easily.

Controller Design
In this paper, two PID controllers were applied to control m H2 and P c independently.The formula for a PID controller is shown as where k p , k i , and k d are the proportional gain, integral gain, and derivative gain, respectively.In previous research, the main work has been to search for the optimal values of k p , k i , and k d [17,19,28].However, the filter coefficient N is also of great importance, and an appropriate filter coefficient can reduce the output oscillation and improve the control performance.As a result, the filter coefficient was also studied for this paper, and there were eight parameters that needed to be decided on, which were the parameters of PID1 (k p1 , k i1 , k d1 , N 1 ) and the parameters of PID2 (k p2 , k i2 , k d2 , N 2 ).The structure of the control system is shown in Figure 3, and the parameters of the controllers were initialized by the pidTuner Toolbox in Matlab as Sustainability 2019, 10, x FOR PEER REVIEW 8 of 20 The structure of the control system is shown in Figure 3, and the parameters of the controllers were initialized by the pidTuner Toolbox in Matlab as

Introduction of Multiobjective Optimization based on a Genetic Algorithm
Multiobjective optimization usually involves two or more objective functions, and it is widely used in the field, where objectives are conflicting and need to be optimized at the same time.As is shown in Figure 4, the optimization result is a set of Pareto optimal solutions, where none of the objectives can be promoted without deteriorating other objectives.A multiobjective optimization problem (MOP) can be described as where J is the vector of objective functions, n is the number of objective functions, x represents the decision vector, and X is the constraint of the decision vector x.The objective functions are designed by the decision-maker to judge the quality of the solution, and in the field of industrial control, objective functions usually consist of an IAE, overshoot, and settling time.
To solve optimization problems, many methods have been developed, and the GA proposed by Holland in 1975 [29] is one of the most popular ones.Meanwhile, in 1994, Srinivas and Deb [30] developed a new algorithm based on a GA called a nondominated sorting genetic algorithm (NSGA), and then in 2000, an improved version of the NSGA (called the NSGA-II) was developed by Deb et al. [31].With the merits of low computational complexities and the ability of parallel computation, it is widely used for MOPs [21,32].Afterwards, in order to improve the distribution density of the Pareto front and avoid the problem of premature ripening, an improved NSGA-II was developed by introducing the parameter Pareto faction x, which represents the ratio of the optimal solutions to the total population (Equation(43)): where N represents the population size, and NPareto represents the number of Pareto optimal solutions.A flow chart of this algorithm is shown in Figure 5.
In this section, the proposed improved NSGA-II is applied to find the optimal sets of control parameters for the SOFC.

Introduction of Multiobjective Optimization based on a Genetic Algorithm
Multiobjective optimization usually involves two or more objective functions, and it is widely used in the field, where objectives are conflicting and need to be optimized at the same time.As is shown in Figure 4, the optimization result is a set of Pareto optimal solutions, where none of the objectives can be promoted without deteriorating other objectives.A multiobjective optimization problem (MOP) can be described as where J is the vector of objective functions, n is the number of objective functions, x represents the decision vector, and X is the constraint of the decision vector x.The objective functions are designed by the decision-maker to judge the quality of the solution, and in the field of industrial control, objective functions usually consist of an IAE, overshoot, and settling time.
To solve optimization problems, many methods have been developed, and the GA proposed by Holland in 1975 [29] is one of the most popular ones.Meanwhile, in 1994, Srinivas and Deb [30] developed a new algorithm based on a GA called a nondominated sorting genetic algorithm (NSGA), and then in 2000, an improved version of the NSGA (called the NSGA-II) was developed by Deb et al. [31].With the merits of low computational complexities and the ability of parallel computation, it is widely used for MOPs [21,32].Afterwards, in order to improve the distribution density of the Pareto front and avoid the problem of premature ripening, an improved NSGA-II was developed by introducing the parameter Pareto faction x, which represents the ratio of the optimal solutions to the total population (Equation (43)): where N represents the population size, and N Pareto represents the number of Pareto optimal solutions.A flow chart of this algorithm is shown in Figure 5.In this section, the proposed improved NSGA-II is applied to find the optimal sets of control parameters for the SOFC.

Optimization of the Control Parameters of the SOFC
As was mentioned in Section 3, there were eight parameters in the two PID controllers that needed to be optimized: k p1 , k i1 , k d1 , N 1 , k p2 , k i2 , k d2 , and N 2 .Thus, the decision vector could be described as x = (k p1 , k i1 , k d1 , N 1 , k p2 , k i2 , k d2 , N 2 ).The integral absolute error (IAE) of the output voltage and fuel utilization were employed as two objective functions, and the MOP could be formulated as where and the constraint for each parameter is shown as The population size N was set as 50, the maximum generation Max was set as 50, the Pareto fraction x was set as 0.35, the mutation function was a constraint-dependent type, and the crossover probability was 0.8.The stopping criteria were designed as follows: 1.
The number of generations is higher than the maximum number of iterations; 2.
The average relative change in the spread of the Pareto solutions over 100 generations is less than 10 −4 , and the spread is smaller than the average spread over the last 100 generations.
For each generation, four load disturbances occurred in each simulation process, and two objective functions, IAE 1 and IAE 2 , of 50 individuals were calculated based on the simulation results.
The distribution of the final population in the objective function space is shown in Figure 6, which clearly illustrates the conflicting relationship between IAE 1 and IAE 2 .Any decrease in IAE 2 led to an increase in IAE 1 .In other words, any improvement in the control performance of fuel utilization led to a deterioration in the control performance of the output voltage.Five typical points were chosen from the Pareto front to get the optimal parameters of the controllers, and the value of each parameter and the objective function are displayed in Table 1.In Table 1, point A possessed the lowest IAE 1 and the highest IAE 2 , while point E possessed the highest IAE 1 and the lowest IAE 2 .As we can tell, there were trade-offs between the two objectives, and to decide the final optimal solution, a decision-making process should be executed based on the decision-maker's insight.
Figure 7 denotes the variation of the objective functions with the generation.As we can tell from the result, IAE 1 decreased from 22.75 to 1.20 and converged at about 1.197 after 30 generations, while IAE 2 decreased from 0.8 to 0.011 and converged at about 0.0108 after 25 generations.Both IAE 1 and IAE 2 were greatly minimized after about 30 generations, and thus the optimization was successful.

Results and Discussion
In this section, to analyze the control performance of the five sets of optimal control parameters and the initial control parameters, a simulation was performed based on the linearized model implemented, and the results were compared and discussed.

Comparison of the Optimal Points
In this part, a step signal was added to the load resistance, and the disturbance responses of the controllers on the five optimal points and the original controller were compared.The best optimal points were chosen based on the comparison results.

Simulation Results
During the simulation process, a step signal was added to the load resistance to make it change from 5 Ω to 6 Ω at 20 s, while the reference values of V out and FU were 168 V and 80%, respectively.The disturbance responses of the output power P out , the output voltage, the fuel utilization, the fuel flow, and the cathode pressure are shown in Figure 8, and detailed data from the simulation are shown in Table 2.

Results and Discussion
In this section, to analyze the control performance of the five sets of optimal control parameters and the initial control parameters, a simulation was performed based on the linearized model implemented, and the results were compared and discussed.

Comparison of the Optimal Points
In this part, a step signal was added to the load resistance, and the disturbance responses of the controllers on the five optimal points and the original controller were compared.The best optimal points were chosen based on the comparison results.

Simulation Results
During the simulation process, a step signal was added to the load resistance to make it change from 5  to 6  at 20 s, while the reference values of Vout and FU were 168 V and 80%, respectively.The disturbance responses of the output power Pout, the output voltage, the fuel utilization, the fuel flow, and the cathode pressure are shown in Figure 8, and detailed data from the simulation are shown in Table 2.

Discussion
From the simulation results, we can observe that compared to the initial controller before optimization, the controllers tuned to the five optimal points all had better control performance.
In Figure 8b, the controllers on the five optimal points all had faster dynamic responses on Pout compared to the initial controller.At optimal point A, the speed of the dynamic response was the fastest, while optimal point E possessed the slowest dynamic response.In Figure 8c,d, the controllers on the five optimal points all had faster control effects on Pc and mH2.
In Figure 8e, it is clear that the control performance of Vout could be dramatically improved by applying the controllers on the five optimal points, and the settling time and IAE of the output voltage gradually increased from optimal point A to E, while the optimal point C possessed the biggest overshoot among all the optimal points (0.0211%).In Figure 8f, both the settling time and IAE of the response of FU on the five optimal points were better than those before optimization.Among all of the optimal points, the controller on point A possessed the longest settling time (1.225 s), while the settling time for the controller before optimization was 4.832 s.In addition, the controller on optimal point E possessed the best control performance with the shortest settling time, 0.062 3 s, and the lowest IAE, 0.0017.
In conclusion, the controllers on the five optimal points all had a smaller overshoot, lower settling times, and a lower IAE than the initial controller.Based on the simulation results, the control parameters on optimal point D were chosen as the best solution in order to reduce the control action

Discussion
From the simulation results, we can observe that compared to the initial controller before optimization, the controllers tuned to the five optimal points all had better control performance.
In Figure 8b, the controllers on the five optimal points all had faster dynamic responses on P out compared to the initial controller.At optimal point A, the speed of the dynamic response was the fastest, while optimal point E possessed the slowest dynamic response.In Figure 8c,d, the controllers on the five optimal points all had faster control effects on P c and m H2 .
In Figure 8e, it is clear that the control performance of V out could be dramatically improved by applying the controllers on the five optimal points, and the settling time and IAE of the output voltage gradually increased from optimal point A to E, while the optimal point C possessed the biggest overshoot among all the optimal points (0.0211%).In Figure 8f, both the settling time and IAE of the response of FU on the five optimal points were better than those before optimization.Among all of the optimal points, the controller on point A possessed the longest settling time (1.225 s), while the settling time for the controller before optimization was 4.832 s.In addition, the controller on optimal point E possessed the best control performance with the shortest settling time, 0.062 3 s, and the lowest IAE, 0.0017.
In conclusion, the controllers on the five optimal points all had a smaller overshoot, lower settling times, and a lower IAE than the initial controller.Based on the simulation results, the control parameters on optimal point D were chosen as the best solution in order to reduce the control action and improve the control performance of both the output voltage and fuel utilization.
Additionally, there were some defects in the optimization process.Compared to the improved NSGA-II we applied in this paper, the NSGA-III proposed by Deb and Jain [33] in 2014 is more proper for an MOP with a Pareto front, whose objective values may have a different scale.Moreover, the population size and execution times applied were the default values in the Matlab Optimization Toolbox, which could be optimized to save on optimization time.

Continuous Disturbance Response Simulation
In this part, continuous step signals were added to the load resistance and the current to simulate the realistic working conditions of an SOFC.The disturbance responses of the optimal controller based on optimal point D and the original controller were compared, and conclusions were drawn based on the simulation results.

Resistance Disturbance Response Simulation
During the simulation process, a continuous step signal was added to the load resistance to make it change from 5 Ω to 5.5 Ω, 6 Ω, 6.5 Ω, and 5 Ω at 200 s, 400 s, 600 s, and 800 s, respectively.The reference values of V out and FU were kept at 168 V and 80%, respectively.The disturbance responses of the output power P out , the output voltage, the fuel utilization, the fuel flow, and the cathode pressure are shown in Figure 9, and detailed data from the simulation are shown in Table 3.
Sustainability 2019, 10, x FOR PEER REVIEW 14 of 20 Additionally, there were some defects in the optimization process.Compared to the improved NSGA-II we applied in this paper, the NSGA-III proposed by Deb and Jain [33] in 2014 is more proper for an MOP with a Pareto front, whose objective values may have a different scale.Moreover, the population size and execution times applied were the default values in the Matlab Optimization Toolbox, which could be optimized to save on optimization time.

Continuous Disturbance Response Simulation
In this part, continuous step signals were added to the load resistance and the current to simulate the realistic working conditions of an SOFC.The disturbance responses of the optimal controller based on optimal point D and the original controller were compared, and conclusions were drawn based on the simulation results.

Resistance Disturbance Response Simulation
During the simulation process, a continuous step signal was added to the load resistance to make it change from 5  to 5.5 , 6 , 6.5 , and 5  at 200 s, 400 s, 600 s, and 800 s, respectively.The reference values of Vout and FU were kept at 168 V and 80%, respectively.The disturbance responses of the output power Pout, the output voltage, the fuel utilization, the fuel flow, and the cathode pressure are shown in Figure 9, and detailed data from the simulation are shown in Table 3.During the simulation process, a continuous step signal was added to the current to make it change from 33 A to 27 A, 29 A, 31 A, and 33 A at 200 s, 400 s, 600 s, and 800 s, respectively.The reference values of Vout and FU were kept at 168 V and 80%, respectively.The disturbance responses of the output power Pout, the output voltage, the fuel utilization, the fuel flow, and the cathode pressure are shown in Figure 10, and detailed data from the simulation are shown in Table 4.
(a)  During the simulation process, a continuous step signal was added to the current to make it change from 33 A to 27 A, 29 A, 31 A, and 33 A at 200 s, 400 s, 600 s, and 800 s, respectively.The reference values of V out and FU were kept at 168 V and 80%, respectively.The disturbance responses of the output power P out , the output voltage, the fuel utilization, the fuel flow, and the cathode pressure are shown in Figure 10, and detailed data from the simulation are shown in Table 4.During the simulation process, a continuous step signal was added to the current to make it change from 33 A to 27 A, 29 A, 31 A, and 33 A at 200 s, 400 s, 600 s, and 800 s, respectively.The reference values of Vout and FU were kept at 168 V and 80%, respectively.The disturbance responses of the output power Pout, the output voltage, the fuel utilization, the fuel flow, and the cathode pressure are shown in Figure 10, and detailed data from the simulation are shown in Table 4.

Discussion
From the simulation results, it is evident that compared to the original controller, the optimal controller based on the controller parameters on optimal point D had better control performance.
In Figures 9b and 10b, the step responses of the optimal controller possessed faster dynamic responses than the original one, which led to fewer fluctuations in the output power.Figure 9c,d and Figure 10c,d, also illustrate that the optimal controller had more rapid dynamic responses to cathode pressure and fuel flow, which permitted faster control action in the controller.
From Figures 9e and 10e, it is clear that the control performance of the optimal controller was much better than the original controller from the perspective of the overshoot and settling time.
The overshoot of all of the disturbance responses of the optimal controller was 0, whereas the original controller had a slight overshoot.As for the settling time, the settling time of the optimal controller under resistance disturbance varied from 8.5 s to 16.4 s, while that of the original controller under resistance disturbance varied from 18.1 s to 46.2 s.In addition, the settling time of the optimal controller under current disturbance varied from 12.1 s to 17.3 s, while that of the original controller under resistance disturbance varied from 37.6 s to 50.3 s.
In Figures 9f and 10f, although the overshoot of both the original controller and optimal controller was 0, the settling time of the optimal controller was much shorter due to the rapid control action.The settling time of the optimal controller under resistance disturbance varied from 1.1 s to 2.4 s, and as for the original controller, it varied from 4.8 s to 6.8 s.The settling time of the optimal controller under current disturbance varied from 1.8 s to 2.2 s, and the original controller varied from 5.2 s to 6.7 s.
In conclusion, from the simulation results of the continuous disturbance response, it is clear that the optimal controller on optimal point D had the merits of stability, a smaller overshoot, and reduced settling time, and consequently, it is the ideal controller for the SOFC model under consideration.

Conclusions
An SOFC is an efficient energy conversion device that has been widely used.To deal with the characteristics of couplings of multivariables, nonlinearity of the model, and frequent load disturbances, this paper employed a GA-based optimization algorithm to find the optimal parameters for the PID controllers of an SOFC.In this paper, a 3 × 2 MIMO system was identified based on the SOFC model, and two PID controllers were designed to separately control V out and FU.Then, the optimization objectives were designed, and a GA was applied for multiobjective optimizations of the controller parameters for the two controllers.Finally, a set of optimal points on the Pareto front were found, and the corresponding control parameters of the five points were applied to control the linearized SOFC model.The results showed that the controller parameters at the five optimal points resulted in better control performance compared to before optimization in terms of the IAE, the overshoot, and the settling time.The parameters on optimal point D were considered to be the best for their smooth control action.Therefore, PID controllers with the controller parameters on optimal point D are an ideal choice for SOFC control.

Figure 1 .
Figure 1.A structural diagram of a solid oxide fuel cell (SOFC).

Figure 1 .
Figure 1.A structural diagram of a solid oxide fuel cell (SOFC).

Figure 2 .
Figure 2. Step response curves of the identified model and the original model.(a) Comparison of step response curves of G 11 ; (b) Comparison of step response curves of G 21 ; (c) Comparison of step response curves of G 12 ; (e) Comparison of step response curves of G 22 ; (e) Comparison of step response curves of d 1 ; (f) Comparison of step response curves of d 2 .

Figure 3 .
Figure 3. Structure of the control system.

Figure 3 .
Figure 3. Structure of the control system.

Figure 6 .
Figure 6.The distribution of the final generation in the objective function space.

Figure 7 .
Figure 7.The trend lines of the objective functions: (a) The trend line of IAE1; (b) the trend line of IAE2.

Figure 6 .
Figure 6.The distribution of the final generation in the objective function space.

Figure 6 .
Figure 6.The distribution of the final generation in the objective function space.

Figure 7 .
Figure 7.The trend lines of the objective functions: (a) The trend line of IAE1; (b) the trend line of IAE2.

Figure 7 .
Figure 7.The trend lines of the objective functions: (a) The trend line of IAE 1 ; (b) the trend line of IAE 2 .

Figure 8 .
Figure 8. Simulation results of a step response: (a) Step input of R; (b) disturbance response of P out ; (c) disturbance response of P c ; (d) disturbance response of m H2 ; (e) disturbance response of V out ; (f) disturbance response of FU.

Figure 9 .
Figure 9. Simulation results of a continuous step response: (a) Step input of R; (b) disturbance response of Pout; (c) disturbance response of Pc; (d) disturbance response of mH2; (e) disturbance response of Vout; (f) disturbance response of FU.

Figure 9 .
Figure 9. Simulation results of a continuous step response: (a) Step input of R; (b) disturbance response of P out ; (c) disturbance response of P c ; (d) disturbance response of m H2 ; (e) disturbance response of V out ; (f) disturbance response of FU.

Figure 9 .
Figure 9. Simulation results of a continuous step response: (a) Step input of R; (b) disturbance response of Pout; (c) disturbance response of Pc; (d) disturbance response of mH2; (e) disturbance response of Vout; (f) disturbance response of FU.

Figure 10 .
Figure 10.Simulation results of a continuous step response: (a) Step input of the current; (b) disturbance response of P out ; (c) disturbance response of P c ; (d) disturbance response of m H2 ; (e) disturbance response of V out ; (f) disturbance response of FU.

Table 1 .
The value of the control parameters and objective functions at the five optimal points.

Table 1 .
The value of the control parameters and objective functions at the five optimal points.

Table 1 .
The value of the control parameters and objective functions at the five optimal points.

Table 2 .
Detailed data from the simulation results.

Table 2 .
Detailed data from the simulation results.

Table 3 .
Detailed data from the simulation results.

Table 3 .
Detailed data from the simulation results.

Table 3 .
Detailed data from the simulation results.

Table 4 .
Detailed data from the simulation results.

Table 4 .
Detailed data from the simulation results.