Determination of Optimal Parameters for Dual-Layer Cathode of Polymer Electrolyte Fuel Cell Using Computational Intelligence-Aided Design

Because of the demands for sustainable and renewable energy, fuel cells have become increasingly popular, particularly the polymer electrolyte fuel cell (PEFC). Among the various components, the cathode plays a key role in the operation of a PEFC. In this study, a quantitative dual-layer cathode model was proposed for determining the optimal parameters that minimize the over-potential difference and improve the efficiency using a newly developed bat swarm algorithm with a variable population embedded in the computational intelligence-aided design. The simulation results were in agreement with previously reported results, suggesting that the proposed technique has potential applications for automating and optimizing the design of PEFCs.


Introduction
As a result of the increasing need for an efficient and clean energy supply, considerable importance has been placed on the advancement and fundamental research of polymer electrolyte fuel cell (PEFC) technology. Among the components of PEFCs, the cathode plays a key role in the operation of PEFCs, in which an oxygen reduction reaction (ORR) occurs and generates heat. Platinum (Pt) loading, ionic conductivity, and the reaction's exchange current density are among the factors that may affect the performance. Numerous studies have been conducted to develop models and approaches that are essential to battery performance and optimization. geometry design problem for optimization of single serpentine and transient characteristics of PEFC with parallel and interdigitated flow fields using a threedimensional, two-phase model. Jung et al. [27] developed an elaborate simulation model of the fuel cell stack system. Askarzadeh and Rezazadeh [28] proposed an innovative global harmony search algorithm for parameter identification of a SR-12 Modular polymer electrode membrane(PEM) Generator. Wang et al. [29] carried out the parameter sensitivity analysis for a three-dimensional, two-phase, non-isothermal model of polymer electrolyte membrane fuel cell. Chen et al. [30] proposed a quantitative approach for predicting the remaining battery life by using an adaptive, bathtub-shaped function. Considering thermoelectric and thermoeconomic objectives, Sayyaadi and Esmaeilzadeh [31] developed a methodology for optimal PEFC control, in which the net power density and energetic efficiency are maximized. Pathak and Basu [32] discussed a mathematical model for the anode and cathode with an anion-exchange membrane for predicting the performance of a fuel cell considering reaction kinetics and ohmic resistance effects. Noorkami et al. [33] investigated the temperature uncertainty as a key parameter in determining the performance and durability of a PEFC. Molaeimanesh and Akbari [34] proposed a threedimensional lattice Boltzmann model of a PEFC cathode, in which the electrochemical reaction on the catalyst layer is able to simulate single-and multispecies reactive flow in a heterogeneous, anisotropic gas diffusion layer. Wang et al. [35] studied a three-dimensional, two-phase, and non-isothermal fuel cell model incorporating the Leverett-Udell correlation and evaluated its performance.
Although there have been a large number of previous studies, the available literature on the analytical modeling of cathode electrodes fails to address two concerns. First, the previous studies do not capture the coupling effects on PEFC performance resulting from the interactions among the design variables. Second, few effective methods have been developed that allow for quantitative analysis, model verification, and parameter optimization. To fill this void, this paper proposes a bat swarm algorithm with a variable population (BAVP) to construct and optimize the quantitative cathode electrode model, which will be embedded into the computational intelligence-aided design (CIAD) [36] framework. This new CIAD framework provides an expanded capability to accommodate a variety of CI algorithms, and it has three advantages: (1) mobilizing computational resources; (2) taking advantage of multiple CI algorithms; and (3) reducing computational costs. This framework has been demonstrated in some of our previous works in diverse areas: applied energy [30], new drug development for public healthcare [37,38], economy and finance [39], sustainable development [40][41][42], aerospace engineering [43], automotive engineering [44], public security [45], and engineering modeling and design [46,47], among others.
Inspired by the echolocation behavior of bats and first proposed by Yang [48] in 2010, the bat swarm algorithm (BA) allocates computational resources by adjusting its population and accelerating the calculation speed. By using echolocation, a swarming bat can quickly respond to changes in the direction and speed of its neighbors during activities such as detecting prey, avoiding obstacles, and locating roosting crevices in dark surroundings. Useful behavioral information is passed among bats and guides them to move from one configuration to another as one unit. By borrowing this intelligence of social behavior, the BAVP is parallel, independent of initial values, and able to achieve a global optimum.
This work has three main contributions. First, our model can effectively assess PEFC performance. Second, a BAVP swarm intelligence method is devised as the search engine to optimize the model parameters embedded in the CIAD framework. Third, two new metrics, the index of moving mean of the average precision (mmAP) and the index of moving mean of variance (mmVAR), are introduced to characterize the dynamic evolutionary behaviors.
The remainder of this paper is organized as follow: Section 2 (Analytical Modeling) discusses the analytical modeling of the cathode electrode; Section 4 (Computational intelligence-aided design) describes the conceptual framework of CIAD and the integrated solver; Section 3 (Bat swarm algorithm with variable population) describes the BAVP algorithm for the optimization; Section 5 (Optimization and Parameter Determination) defines the fitness function for optimizing the analytical model using the model proposed in Section 2; Section 6 (Empirical Results and Discussion) presents the empirical results and further verifies the optimal design; and Section 7 (Conclusions and Future Works) concludes the paper.

Analytical Modeling
A schematic diagram of a dual-layer configuration of a cathode electrode is shown in Fig. 1, in which five specific areas are numbered and are explained below. The left side of the electrode attaches to the PEM, and the right side connects to the diffusion media [19,20].
N Includes the assumption that the oxygen concentrations, temperatures, electronic phase potentials, and equilibrium potential are the same between the two layers and are uniform within each layer. The electrodes are thin layers (ƒ10mm) coated on the PEM surface containing a catalyst (typically Pt), carbon(C), an ionomer electrolyte and void space. In general, there are three phases in the electrode: (i) void space for the transport of gaseous reactants, (ii) ionomer content for the transfer of protons, and (iii) carbon support for conducting electrical current. In addition to the electrochemical catalyst, which is essential for all functions, equation (1) is given.
N The two sub-layers are denoted as 'Layer 1' and 'Layer 2'. Five parameters are considered in this model that include the ionic conductivity s m , the catalyst specific area a, the exchange current density i, the ionic resistance R d , the current density I d , the thickness of sub-layer d, the interface location of the two sub-layers l, and so on.
The ionic conductivity factors s m of 'Layer 1' and 'Layer 2' are s 1 and s 2 , respectively, and they are determined from the electrolyte water content l, the ionomer tortuosity t m , the Nafion content E m and the temperature T, as given in equation (2). The ratio of the ionic conductivity factors of 'Layer 1' and 'Layer 2' is given in equation (3). The catalyst specific area a describes the active catalyst surface area per unit volume. The exchange current density i depends on factors such as temperature and the electrochemical characteristics of the catalyst. As shown in equation (4), the factor of the catalyst specific area and the exchange current density multiplication ai is determined by factors such as the structural feature of the electrode, including the reaction interface roughness and the mean radius of the catalyst particles, and is the most important factor for catalyst cost reduction, where E a is the activation energy for the ORR; R g is the universal gas constant, 8.314 J/mol K; and s is the liquid water saturation. The ratio of ai of 'Layer 1' and 'Layer 2' is given in equation (5).
N A lumped variable DU is defined in equation (6), in which R d~d =s m is the overall ionic resistance across the cathode electrode, and I d~{ j d d is the current density based on the transfer current density j d at the interface between the two electrodes.
N The relative location of the interface between the two sub-layers is defined in equation (7), in which is the total thickness of the dual-layer electrodes.
N The thickness ratio of the two sub-layers r d is defined in equation (8) r d~d Considering the cathode electrode in one dimension (x direction), the two indices (g 1 and g 2 ) of the over-potential difference of 'Layer 1' and 'Layer 2' are given in equations (9) and (10), respectively, where P, Y and V are defined in equations (11) to (13) [21].
Bat swarm algorithm with variable population Because the BAVP is inspired by the echolocation characteristics of bat swarms, it can be idealized to include the four following assumptions: N 1 As shown in Fig. 2, all artificial bats (ABs) utilize the same echolocation mechanism to measure distance, and each AB individual B i is able to detect the difference between prey (food) and obstacles. N 2 Each individual B i can generate ultrasounds to echolocate the prey and obstacles with a velocity of u i,j and a position of x i,j at time j, which are stated in Equations (15) and (14), respectively, where x Ã is the current global best position.
N 4 As shown in Equation (17), the population P j of ABs varies from time j to another, which accelerates the optimization process, in which P N is the nonreplaceable population and P R j is the replaceable population at time j.
N Step (4), generate new local solutions x 0 using Equation (18), where E [ [21,1] is a random-walk factor. As defined in Equation (19), A i,j is the loudness of the bat B i at time j, in which a [ [0,1] is a reduction factor.  in Equation (20).
N Step (6), continue running the calculation until the terminal conditions have been satisfied.

Computational intelligence-aided design
Computational intelligence (CI) is a set of nature-inspired approaches that provides numerous capabilities for solving complex problems. Compared to the traditional optimization methods, CI does not need to reformulate the problem to search a non-linear or non-differentiable space. Another advantage of CI is its flexibility in formulating the fitness function, which can be expressed as a function of the system output. This feature is particularly appealing if an explicit objective function is difficult to obtain. Fig. 4 illustrates the CIAD framework, and the entire optimization process can be summarized in the following three main steps: N Step 1, pre-process. In this step, quantitative models under specific conditions are obtained for engineering applications.
N Step 2, optimal design. This step defines the fitness functions according to the design objectives.
N Step 3, post-process. This step produces the final results and completes the post-processing tasks. Specifically, this step reports the optimal solution, analyses and visualizes the results, and presents the recommendations to policy makers. The 'CI integrated solver'(CIS) is employed to optimize the parameters for the fitness function, and the details of the CIS are given in Fig. 5.
As shown in Fig. 5, the conceptual framework of CIS consists of three parts: data input, the CI integrated solver and result output, as follows. Part 1: Data Input (@point A). This part prepares the data input for the CI integrated solver. It collects, filters, stores, and pre-processes data originating from various sources, such as statistical yearbooks, research analyses, and government reports. Part 2: CI Integrated Solver. In this part, a set of nature-inspired computational approaches are integrated into one solver to optimize complex real-world problems, which primarily involves one or more of the following methods: artificial neural networks, a genetic algorithm, fuzzy logic, simulated annealing, artificial immune algorithms, and swarm intelligence algorithms. In this paper, a BAVP algorithm (@point C) is  Dual-Layer Fuel Cell via CIAD embedded in this solver, and the details of this algorithm are discussed in Section 3. Part 3: Result Output (@point B). This part reports the final results from Part 2. As shown in Fig. 4, the data-flow from Steps 4 to 3 is the input of the 'CI integrated solver' interconnected with point A, and the data-flow from Steps 3 to 2 is the output of the 'CI integrated solver' interconnected with point B.

Optimization and Parameter Determination
To determine the optimal parameters for the over-potential difference g, this section introduces two trend indices mmAP and mmVAR for evolutionary optimization, which are given in equations (21) and (22), respectively. As stated in Equation (21), the index of mmAP is a moving average score of the mean value of vector f j , where i~1,2, Á Á Á ,p, p is the population of the data set, and MEAN : ð Þ is the average function. The index of mmVAR is a moving average score of the VAR value of vector f j , as given in Equation (22), where VAR : ð Þ is the variance function. The two indices are employed to assess the short-term fluctuations by capturing the longer-term trend across the evolutionary process. Dual-Layer Fuel Cell via CIAD In Fig. 6, the solid line represents the mmAP scores for each vector f j as given in Equation (21). The dashed lines are the mmAP + mmVAR for each vector f j as given in equation (22), which defines the limits of evolutionary paths of the optimization process (generation versus fitness f ) as the upper and lower boundaries.
The fitness function is in a reciprocal form of the over-potential difference function, as given in Equation (23). The fitness function is defined as the mmAP reciprocal function of the over-potential difference function g, in which maximizing is a way to minimize g, and the goal of this function is to determine the optimal combination of five parameters, DU, r s , r ai , r d and l, that simultaneously minimizes the objective of g. eps is the floating-point relative accuracy, which prevents singularity in the case where g is approaching 0 and is approaching ?.
The over-potential difference function g is given in Equation (24), where DU is the lumped variable, given in Equation (6); r s is the ratio of the ionic conductivity of the two sub-layers, given in Equation (3); r ai is the ratio of ai, given in Equation Dual-Layer Fuel Cell via CIAD (5); r d is the ratio of the thickness, given in Equation (8); and l is the location factor, given in Equation (7).

Empirical Results and Discussion
Maximizing the fitness function yields the minimum of g, which is performed using the specially designed toolboxes SwarmBat [49] and SECFLAB [50]. The computer specifications for the simulations are a 2.1 GHz Intel dual-core processor, Windows XP Professional v5.01 Build 2600 service pack 3, a 2.0 GB 800 MHz dual channel DDR2 SDRAM, and MATLAB R2008a. The initial parameters are listed in Table 1, in which the max-generation number is 100, and it serves as the termination condition in each test. The test number is also 100. The frequency range is set to [20000,500000]Hz. The reduction factor a is 0.9. The population is 50, in which the non-replaceable P N and replaceable population P R are 40 and 10, respectively. The random step is 0.01, and the ranges of DU, r s , r ai , r d , and l are [0,10], [0,10], [0,10], [0,2] and [0,3], respectively. Table 2 presents the optimal combinations (MEAN+VAR) of DU, r s , r ai , r d and l, which indicates that the over-potential is non-uniform within the cathode Table 3. Impacts of five coupled variables on g in 3D figures. Dual-Layer Fuel Cell via CIAD and at particularly high values of the lumped parameter DU and is sensitive to the spatial variation l. Fig. 7 shows the mmAP curves with the upper and lower mmVAR boundaries, in which the mmVAR boundaries stick to the mmAP fitness curves and the fitness increases very quickly; it reaches a plateau from generations 1 to 60 (or so), and it remains steady from generation 60 to 100. Note that all lines converge in generation 100. Fig. 8 shows the fitness mmVAR over the entire simulation. The curves decline quickly within approximately 60 generations and finally reach 0 in generation 100. Figs. 7 and 8 indicate that the proposed optimization algorithm is efficient and accurate.
As also listed in Table 3, to demonstrate the impacts of the five coupled variables on g, Figs. 9 and 15 provide seven '3D' figures to evaluate these impacts.
Specifically, Figs. 9 and 10 show that r s and r ai have similar positive effects on g, indicating that when r s and r ai increase, g increases, and vice versa. Furthermore, when r s and r ai remain constant, DU has limited effects on g. Fig. 11 shows that g is sensitive to r d [ [3,4] and DU v 0.05. Fig. 12 indicates that g increases faster with l w 0.5 and DU v 0.05 and that better g values are obtained with larger l and smaller DU.  Dual-Layer Fuel Cell via CIAD Fig. 13 shows that g is sensitive to smaller r d or DU. In Fig. 14, there is a plateau within r s v 1 and l v 0.5; furthermore, g is sensitive to larger l when r s v 1. When r s ± w 1 and increases, g decreases within the full l range. Fig. 15 shows that g increases when both r ai and l become smaller, which implies that g is unstable with small values of r ai and l.

Conclusions and Future Works
In this study, an analytical model that incorporates five parameters is proposed to explore the transport and electrochemical phenomena in dual-layered cathode electrodes of polymer electrolyte fuel cells. These parameters include the lumped variable DU, the ratio of the ionic conductivity of two sub-layers r s , the ai ratio of the two sub-layers r ai , the ratio of the thickness r d and the relative location factor l. Moreover, a theoretical study on the spatial distribution of reaction rates across the electrode is presented.
The proposed model is utilized to define a design objective: determining the optimal combinations of the five parameters to minimize the over-potential difference g. Based on the trend indices mmAP and mmVAR, a fitness function  Dual-Layer Fuel Cell via CIAD was constructed with the five variables as discussed above, which are optimized by the bat swarm algorithm with a variable population.
The numerical solutions obtained in this study were applied to optimize the electrode performance through a set of optimal dual-layer configurations, and the research findings are summarized in the following three points: 1. The proposed dual-layered cathode electrode model for the determination of the optimal parameters provides a strong argument for implementing the solutions to explore the impacts of each layer's properties on their performance. 2. Based on the developed dual-layered cathode electrode model, a bat swarm algorithm with a variable population is developed, which directly affects the determination of the optimal parameters due to its high efficiency and accuracy. 3. The proposed two trend indices mmAP and mmVAR were utilized to smooth out short-term fluctuations and highlight longer-term trends until the maximum generation fitness point was reached, which helps to measure the computational performance of the BAVP or to deploy other algorithms.

Dual-Layer Fuel Cell via CIAD
Our future research will focus on developing new types of CI algorithms, such as the swarm dolphin algorithm [51], the swarm wolf algorithm and their hybrid derivatives, to optimize further geometrical parameters and optimal combinations for improving the efficiency of polymer electrolyte fuel cells with multiple-layer configurations. To achieve a 'state-of-practice' design framework for the fuel cell, further experimental research is needed to establish an advanced model for chemical-dynamical coupled behavior and the potentials of fuel cell commercialization.