A parametric study of the cathode catalyst layer of PEM fuel cells using a pseudo-homogeneous model

International Journal of Hydrogen Energy 26 (2001) 991–999

www.elsevier.com/locate/ijhydene

A parametric study of the cathode catalyst layer of PEM fuel cells using a pseudo-homogeneous model Lixin You, Hongtan Liu ∗ Department of Mechanical Engineering, University of Miami, 208 McArthur Engineering Bldg., Coral Gables, FL 33124, USA

Abstract A pseudo-homogeneous model for the cathode catalyst layer performance in PEM fuel cells is derived from a basic mass– current balance by the control volume approach. The model considers kinetics of oxygen reduction at the catalyst=electrolyte interface, proton transport through the polymer electrolyte and oxygen di/usion through porous media. The governing equations, a two-point boundary problem, are solved using a relaxation method. The numerical results compare well with our experimental data. Using the model, in2uences of various parameters such as overpotential, proton conductivity, catalyst layer porosity, and catalyst surface area on the performance of catalyst layer are quantitatively studied. Based on these results, cathode catalyst layer design parameters can be optimized for speci5ed working conditions. ? 2001 International Association for Hydrogen Energy. Published by Elsevier Science Ltd. All rights reserved. Keywords: PEM; Proton exchange membrane; Fuel cells; Catalyst layer; Mathematical model

1. Introduction Proton exchange membrane (PEM) fuel cell systems are considered to be one of the most promising future vehicular power systems because of their inherent advantages, such as relatively low operating temperatures, high power densities, quick response, quick startup capability, and especially pollution-free operations when pure hydrogen is used. Like any electrochemical devices, a PEM fuel cell consists of an anode, a cathode and an electrolyte. The di/erence is that a PEM fuel cell uses a humidi5ed polymer membrane as electrolyte. In the working temperature range of a PEM fuel cell, the hydrogen oxidation and oxygen reduction rates are very low. The use of catalyst, usually platinum, is necessary to increase reaction rates to practical applicable levels. Because the platinum is very costly and its resource is limited, reduction in catalyst loading is one of the major thrust areas for research. Moreover, the overpotential caused by oxygen reduction in the cathode is one of the main voltage ∗ Corresponding author. Tel.: +1-305-284-2019; fax: +1-305284-2580. E-mail address: [email protected] (H. Liu).

losses in the fuel cell operation. Thus, a good understanding of the various parameter e/ects in catalyst layers is essential for the design and operation of fuel cells. Mathematical models are useful for the analysis and optimization of the performance of fuel cells and especially for cathode catalyst layers in the case of PEM fuel cells. Various mathematical models have been developed for PEM fuel cells, with the emphasis on di/erent aspects, such as transport in the membrane [1–3], transport in the electrode [4 –8], the heat and water management [9,10], and some detailed models for PEM cathodes [11–14]. In these models, the active catalyst layer was either incorporated into the whole PEM model or treated as an ultra-thin 5lm. These models did not provide enough emphasis on the e/ects of di/erent complicated catalyst layer parameters on the catalyst layer performance. Springer and Gottesfeld [15] and Eikerling and Kornyshev [16] presented some excellent analytical solutions for the cathode catalyst layer. However, analytical solutions are only available in limiting conditions when some simpli5cations are assumed, because the governing equations are strongly nonlinear. General solutions to the catalyst layer governing equations are obviously essential in the catalyst design and optimization.

0360-3199/01/$ 20.00 ? 2001 International Association for Hydrogen Energy. Published by Elsevier Science Ltd. All rights reserved. PII: S 0 3 6 0 - 3 1 9 9 ( 0 1 ) 0 0 0 3 5 - 0

992

L. You, H. Liu / International Journal of Hydrogen Energy 26 (2001) 991–999

Nomenclature Av AR m C D E F i0 I m˙ n N R T U y

catalyst surface area per unit volume, cm2 cm−3 membrane resistance, H cm2 concentration, mol m−3 di/usion coeIcient, m2 s−1 potential, V Faraday constant, 96493 C mol catalyst exchange current density, A cm−2 current density, A cm−2 oxygen consumption rate, mol s−1 number of electrons molar 2ux, mol cm−2 s−1 universal gas constant, J mol−1 K −1 temperature, K voltage, V dimension through the catalyst layer, cm

Greek symbols a anodic charge transfer coeIcient c cathodic charge transfer coeIcient  dimensionless oxygen concentration

     

catalyst layer thickness, cm porosity potential, V overpotential, V ionic conductivity, H−1 cm−1 dimensionless length through catalyst layer

Subscripts and superscripts c cathode catalyst layer e/ e/ective property which accounts for porosity eq equilibrium f ionomer 5lm lim limiting current density O2 oxygen oc open circuit solid electrode surface y dimension through catalyst layer  catalyst layer=membrane interface

In this paper, we present a pseudo-homogeneous model for the cathode catalyst layer. The model is derived from a basic mass–current balance by the control volume approach. A relaxation method is used to solve the two-point boundary problem. The e/ects of di/erent parameters on the catalyst layer are studied in detail. The discussions in this paper may help in the optimization design of cathode catalyst layer in PEM fuel cells. 2. Mathematical model for the cathode catalyst layer Fig. 1 sketches the structure of a typical cathode catalyst layer. Four di/erent media are usually present for the function of catalyst layer: a di/usion path for reactants and products transfer, an ionic conducting medium for proton transfer, an electrical conduction medium for electrons, and a catalytic surface for electrochemical reaction to take place. A pseudo-homogeneous model similar to Weisbrod et al. [11] is used here to study the catalyst layer in this paper. In this model, the catalytic layer can be described as a macro-pseudo-homogeneous 5lm with the following assumptions: (1) The carbon matrix has good electron conductivity, thus the potential drop among solid matrix can be neglected compared to the potential drop in the electrolyte phase; (2) The catalyst Pt particle is uniformly distributed; (3) The water content within the ionomer is constant; and (4) The oxygen di/usion coeIcient is constant. We consider the composite catalyst layer in Fig. 1 as an e/ective, homogeneous medium of uniform thickness.

Fig. 1. Schematic illustration of the cathode catalyst layer.

Consider a small control volume as shown in Fig. 2. The decrease of oxygen molar 2ux should be balanced by the increase of proton current density. From Fig. 2, the oxygen consumption rate in the control volume is m˙ s = Ny − Ny+dy = −

dN dy: dy

(1)

L. You, H. Liu / International Journal of Hydrogen Energy 26 (2001) 991–999

993

Finally, according to Fick’s law, the oxygen 2ux is related to oxygen concentration by NO2 = − DOe/2

d(CO2 )f ; dy

(12)

where DOe/2 is the e/ective oxygen di/usion coeIcient, and can be determined from DOe/2 = 1:5 DO2 :

Thus, four unknowns, Iy ; ; NO2 ; CO2 in the catalyst layer, can be described by Eqs. (3), (9), (11) and (12). The appropriate boundary conditions are listed below. At y = 0 (the gas di/user=catalyst layer interface)

Fig. 2. Oxygen mass 2ux and current balance.

The balance of current density requires dN Iy+dy = Iy + m˙ s 4F = Iy − 4F dy dy or dN dIy = − 4F : dy dy

(2)

(3)

(4)

When electrochemical reaction takes place, the electrode potential deviates from equilibrium potential, and it is de5ned as  = E − Eeq = solid − :

Iy = 0; NO2 =

(14) I ; nF

(15)

: CO2 = COy=0 2

The half-electrochemical reaction in the cathode is O2 + 4H+ + 4e− → 2H2 O:

(13)

(5)

The kinetic expression for the oxygen reduction rate per unit volume can be described by using the Butler–Volmer expression with the assumption that reduction current is positive
    a F −c F i = − i0+ exp − exp ; (6) RT RT where i0+ is exchange current density at the equilibrium potential   COf 2 + ref + : (7) i0 = i0 COref2 The change in current is proportional to reaction current density and catalyst surface area per unit volume Av dIy (8) = iAv dy or
    dIy a F −c F + = − Av i0 exp − exp : (9) dy RT RT When water velocity is neglected, the proton current is related to the polymer electrolyte potential di/erence d Iy : (10) = dy e/ From (5) and (10) and assuming the matrix of solid phase to be equipotential (as assumption 1), changes in overpotential can be expressed as d(−) Iy : (11) = dy e/

(16)

At y =  (the catalyst layer=membrane interface) d(−) I : = dy e/

(17)

Because of the resistance of porous electrode, the oxygen transport 2ux across the catalyst layer is limited. The maximum oxygen 2ux can be de5ned as Nmax = DOe/2

COy=0 2 : 

(18)

The corresponding current density is de5ned as the limiting current density Ilim = 4FNmax =

4FDOe/2 COy=0 2 : 

By using the dimensionless variables y = ; c N∗ = =

NO 2 NO2 = ; Nmax DOe/ COy=0 2 2

(19)

(20) (21)

CO 2 ; COy=0 2

(22)

Iy ; Ilim

(23)

I∗ =

the following dimensionless governing equations and boundary conditions can be obtained from Eqs. (3), (9), (11) and (12), dN ∗ dI ∗ ; =− d d
    c F a F dI ∗ Av i0+ exp − = − exp ; d Ilim RT RT

(24) (25)

994

L. You, H. Liu / International Journal of Hydrogen Energy 26 (2001) 991–999

d(−) Ilim  ∗ I ; = d e/ d = − N ∗: d

(26)

Neglecting the overpotential loss in the anode side, the fuel cell output voltage is given as

(27)

U = Uoc − IARm − |c |y=c :

The corresponding dimensionless boundary conditions are given as below. At  = 0, I ∗ = 0;

(28)

N ∗ = I∗ ;

(29)

 = 1:

(30)

At  = 1, d(−) Ilim  ∗ I : = d e/

(31)

3. Results and discussions The governing equations (24) – (27) and boundary conditions (28) – (31) in the cathode active catalyst layer are a two-point boundary problem [17]. A relaxation approach is used to solve the problem. The overpotential of the membrane=active catalyst layer interface is 5rst assumed as a known parameter; the current, overpotential, oxygen 2ux and oxygen concentration distribution in the catalyst layer can then be determined.

(32)

As oxygen reduction reaction rate is low, a certain level of overpotential is necessary to maintain a practical current density. However, a high overpotential decreases the fuel cell output voltage, thus the catalyst layer design should provide high current density at a low overpotential. The catalyst layer is the only place where electrochemical reaction takes place in PEM fuel cells, so its characteristics greatly in2uence fuel cell performances. To decrease the catalyst loading and improve the PEM fuel cell performance, the present catalyst layer is made of a mixture of ionomer, PTFE, carbon and platinum particles. The fraction of each component in the mixture a/ects the e/ective porosity, proton conductivity, electrical conductivity and the catalyst surface area, which in2uence the catalyst layer performance. By modeling the in2uences of these parameters on the catalyst layer performance, it is possible to optimize the fraction of each component. 3.1. Validation of the mathematical model To validate the model proposed in this paper, some modeling results are compared with experiments. Since we are not able to conduct tests on an isolated catalyst layer only, the comparison is made possible by comparing the

◦

Fig. 3. The comparison of PEM fuel cell performance using air and oxygen, respectively. Air or O2 side: Pin = 3 bar, Thyd = 60 C, ◦ ◦ uin = 0:35 m s−1 ; Hydrogen side: Pin = 1 bar, Thyd = 60 C, uin = 0:6 m s−1 ; fuel cell operating temperature Tcell = 80 C.

L. You, H. Liu / International Journal of Hydrogen Energy 26 (2001) 991–999

995

Table 1 Cathode catalyst layer parameters for base case Parameter

Value

Catalyst layer porosity,  Catalyst layer thickness,  Catalyst surface area per unit volume, Av Exchange current density, i0+ Cathodic transfer coeIcient, c Anodic transfer coeIcient, a Oxygen di/usion coeIcient in plain medium, D Proton conductivity,  Temperature, T Overpotential at catalyst= membrane interface Oxygen concentration at inlet

0.25 12:9 × 10−4 cm 1:4 × 105 cm2 cm−3 4:84 × 10−8 A cm−2 0.52 0.54 0:03 cm2 s−1 0:05 (H cm)−1 350 K 0:3 V 10 mol m−3

modeling results with experimental data on the fuel cell performance. The detailed mathematical models that include the 2ow channels, the gas di/users, the catalyst layers and the membrane are presented elsewhere [18]. Fig. 3 compares the PEM fuel cell performance when air and oxygen are used as the cathode gas, respectively. Good agreements between the experimental data and the modeling results are obtained. 3.2. The e9ects of overpotential Fig. 4(a) – (c) show the variation of proton current, oxygen concentration and overpotential distribution across the catalyst layer with di/erent overpotential at the membrane=catalyst layer interface. Other parameters are given in Table 1. The electrochemical reaction rate depends on the local oxygen concentration and overpotential. From the gas di/user=catalyst layer interface ( = 0) to the catalyst layer=membrane interface ( = 1), the dimensionless current density increases; dimensionless oxygen concentration decreases due to its consumption; the absolute value of overpotential increases due to ohmic losses across the catalyst layer. Fig. 4(a) shows the variation of dimensionless proton current density with overpotential. When the overpotential is small, the oxygen reduction rate is low, thus only a small current density can be produced. With the increase of overpotential, current density increases, and more and more currents are produced in a region close to the catalyst=gas di/user interface. While ||1 is close to 0.35, the current density reaches the “limiting current density”, which is due to oxygen mass transport limitation. Any further increase of ||1 does not increase the overall current density. At a very high overpotential, only a thin layer of the catalyst close to the gas di/user is active. Fig. 4(b) shows the variation of dimensionless oxygen concentration across the catalyst layer with overpotential. When the overpotential is small (||1 = 0:1), the oxygen reduction rate is very low, little oxygen is consumed, and the

Fig. 4. (a) The variation of current density with overpotential, (b) the variation of oxygen concentration with overpotential, (c) the variation of overpotential across catalyst layer.

996

L. You, H. Liu / International Journal of Hydrogen Energy 26 (2001) 991–999

oxygen concentration remains almost unchanged. With the increase of overpotential, more current is produced and more oxygen is consumed. While ||1 is close to 0.35, oxygen concentration becomes zero at the catalyst=membrane interface, and current density reaches the limiting current density. With further increase of overpotential, the oxygen reduction rate increases, and most oxygen is consumed in a thin region of the catalyst layer, the active region. The region where no reaction occurs is called the inactive region. When an inactive region exists, the part of the catalyst is not used, which increases the non-uniform degradation of the catalyst in the active region. In addition, high overpotential decreases the e/ective fuel cell output, but cannot increase current density because of the oxygen mass-transport limitation. Fig. 4(c) shows the variation of local overpotential across the catalyst layer. The gradient of overpotential is inversely proportional to the current density when proton conductivity is assumed as constant. For low overpotential, the corresponding current density is small, the potential loss in the electrolyte is small, and the overpotential across the width of the catalyst layer remains almost unchanged. For high overpotential, the corresponding current density is high, the electrolyte phase potential loss across the catalyst layer is high, and so the absolute overpotential increases across the catalyst layer. Even in the inactive region, where no current is generated, the potential loss still exists due to the fact that all the protons must transfer through this region. 3.3. The e9ects of proton conductivity Fig. 5(a) – (c) show the in2uences of proton conductivity on the catalyst layer performance. Fig. 5(c) shows the in2uences of proton conductivity on the distribution of overpotential. The comparison is based on the same overpotential at the catalyst layer=membrane interface, which decides the total activation and concentration overpotential. The highest possible current is desired in this case to improve fuel cell performance. Because of the ohmic resistance, there is a potential loss in the electrolyte phase from the catalyst layer=membrane interface to the catalyst layer=di/user interface. With the decrease of proton conductivity, electrolyte phase potential loss increases, but the solid carbon matrix potential remains almost unchanged because of its much higher electron conductivity; therefore, the potential di/erence between liquid electrolyte and solid matrix (i.e., overpotential) decreases. As shown in Fig. 5(c), the overpotential at the catalyst layer=di/user interface is only about 0:18 V when  = 0:005 (H cm)−1 . The limiting current density does not change with proton conductivity, so dimensionless current density as shown in Fig. 5(a) is able to re2ect the magnitude of current density. As can be seen from Fig. 5(a) and (b), when proton conductivity is high, oxygen reduction is di/usion controlled, and most of the reaction occurs in a thin region close to the catalyst layer=gas di/user interface; while when proton conductivity is low, the reaction is proton conduction

Fig. 5. (a) The variation of current density with proton conductivity, (b) the variation of oxygen concentration with proton conductivity, (c) the variation of overpotential with proton conductivity.

L. You, H. Liu / International Journal of Hydrogen Energy 26 (2001) 991–999

997

controlled, and so most of the reaction occurs in a thin region closest to the catalyst layer=membrane interface. When proton conductivity is high, ( = 0:5 (H cm)−1 ), 90% of the current is generated in the region  ¡ 0:5. When proton conductivity is low ( = 0:005 (H cm)−1 ), oxygen consumption and current generation rate are low in the region  ¡ 0:5 though the oxygen concentration is relatively high in this region. In this case, current generation rate is much higher in the region 0:5 ¡  ¡ 1, where almost 90% of the current is generated. It should be pointed out that the overall current density is lower when the proton conductivity is low. High proton conductivity can be achieved by increasing water content and the ionomer fraction in the catalyst layer. While water content is mainly related to operating conditions, ionomer fraction in the catalyst layer is a design parameter. It is also shown that the overall current density approaches the limiting current density at about  = 0:75 when  = 0:5 (H cm)−1 , which means an optimum ionomer fraction exists in the catalyst layer, because enough carbon particle, catalyst particle and porosity for oxygen transport should also be provided in the catalyst layer. 3.4. The e9ects of e9ective porosity of the catalyst layer Fig. 6(a) and (b) show the in2uences of e/ective porosity on the performance of the catalyst layer. Oxygen can reach the catalyst surface through the hydrophobic pores of solid matrix. The catalyst layer must have enough carbon matrix for good electron conductivity, enough Pt particle surface for electrochemical reaction, and enough polymer electrolyte for proton transport, thus pores in the catalyst layer for oxygen transport are limited. The oxygen transport resistance depends on the e/ective porosity. If the e/ective porosity of the catalyst layer is high, oxygen transfer resistance is low and the limiting current density is high. On the contrary, if the porosity of the catalyst layer is low, oxygen transfer resistance is high and the limiting current density is low. As seen in Fig. 6(a), when  = 0:15, the overall current density is the corresponding limiting current density, and it is much lower than the corresponding current densities when  = 0:25 and 0:35. For  = 0:35, the overall current density is far from the limiting current density. The corresponding oxygen concentration in Fig. 6(b) shows the same trend. When  = 0:15, oxygen is quickly consumed, thus only a small portion of the catalyst layer is active. When  = 0:35, the whole catalyst layer is active. 3.5. The e9ects of the active catalyst surface area Fig. 7(a) – (c) show the in2uences of the active catalyst surface area per unit volume on the catalyst layer performance. A high active catalyst surface area can be realized by using 5ner Pt particles and=or by increasing catalyst loading, which usually comes with higher capital and production costs. A high active catalyst surface area may provide more sites for the electrochemical reaction. Thus, if

Fig. 6. (a) The variation of current density with porosity, (b) the variation of oxygen concentration with porosity.

oxygen and proton 2uxes are suIciently high, current density increases with the catalyst surface area. However, as discussed above, when the oxygen 2ux is limited and the electrolyte phase potential drop increases with current density, the excessive surface area may not contribute to the current generation. In contrast, a catalyst surface area that is too low may not provide enough sites for oxygen reduction. As shown in Fig. 7(a), with a higher surface area, the current density is higher at the same cross section of the catalyst layer, though the overall current density is still controlled by the limiting current density. Fig. 7(b) shows the variations of oxygen concentrations with surface areas. When Av = 1:4 × 106 , the active region is limited to  ¡ 0:6, which means that a too-high surface area will cause an “inactive region”; whereas when Av = 1:4 × 104 , the electrochemical reaction will not be strong enough, and the overall current density will be lower than Av = 1:4 × 105 . Thus, there also

998

L. You, H. Liu / International Journal of Hydrogen Energy 26 (2001) 991–999

exists an optimum catalyst surface area. Fig. 7(c) provides the variations of overpotential with catalyst surface areas. It should be pointed out that, although the water content is assumed as constant in this work, the change of liquid saturation in the gas di/user and in the catalyst layer is considered in our two-phase 2ow model [17]. The governing equations can still be solved by the approach presented earlier, the e/ects of water can be easily incorporated into the model by considering the changes in the e/ective oxygen di/usion coeIcient and proton conductivity. There are other parameters that in2uence the performance of the catalyst layer, including exchange current density, the cathodic transfer coeIcient and the anodic transfer coeIcient. Increasing the exchange current density and the cathodic transfer coeIcient and decreasing the anodic transfer coeIcient can increase the current density when the overpotential is given. These parameters are related to the electrochemical reaction, and they can be improved by controlling electrochemical reaction conditions, such as temperature, concentration, and favored electrochemical reaction paths. 4. Conclusions

Fig. 7. (a) The variation of current density with catalyst surface area, (b) the variation of oxygen concentration with catalyst surface area, (c) the variation of overpotential with catalyst surface area.

From the basic mass–current balance, a pseudohomogeneous model for the cathode catalyst layer of a PEM fuel cell is developed. Although the model is relatively straightforward, fuel cell performances predicted by this model compared well with our experimental results. Parametric e/ects on the performance of the cathode catalyst layer in PEM fuel cells are systematically studied. Based on the modeling results, the following can be concluded: (1) When overpotential is low, the overall current density increases with overpotential; when overpotential reaches a critical value, the overall current density reaches the limiting current density. Any further increase of overpotential decreases fuel cell output voltage but cannot increase current density. (2) Generally, high proton conductivity is favorable to electrochemical reaction. When proton conductivity is high, the overpotential across the catalyst layer varies little and most current is generated in the portion of the catalyst layer adjacent to gas di/usion layer, especially at high overpotential. When proton conductivity is low, the absolute value of overpotential decreases from the membrane to the gas di/user layer, which causes the reaction rate drop signi5cantly toward the gas di/user layer. (3) The limiting current density and overall current density increase with e/ective porosity. However, the increase of the catalyst layer e/ective porosity is limited because suIcient carbon, catalyst, and ionomer must be also included in the catalyst layer to ensure good electrical and proton conductivity and enough reaction sites. (4) A high active catalyst surface area is essential to ensure a high enough reaction rate. However, if the limiting

L. You, H. Liu / International Journal of Hydrogen Energy 26 (2001) 991–999

current density is low, too high a catalyst surface area may not produce a higher power output. At high current densities, part of the catalyst layer becomes inactive. Thus, it is a waste to increase the catalyst surface area beyond the optimal value for speci5c design parameters and operating conditions. Generally, we can conclude from this parametric study that there are optimized parameter combinations for the catalyst layer design. Since, the paths for transfer of electrons, protons, and reactants all compete for a larger fraction of the catalyst layer, an increase in any one will inevitably be at the scari5cations of the other two. It seems that the best solution is to minimize the thickness of the catalyst layer while optimizing the combinations of the parameters. We hope this work will elicit more experimental research in the design of catalyst layers. Acknowledgements The 5nancial support of the US Department of Energy’s CARAT program under Contract DE- FC02-98EE50531 is gratefully acknowledged. References [1] Springer TE, Zawodzinski TA, Gottesfeld S. Polymer electrolyte fuel cell model. J Electrochem Soc 1991;138: 2334–42. [2] Eikerling M, Kornyshev AA, Stimming U. Electrophysical properties of electrolyte membranes: a random network model. J Phys Chem B 1997;101:10807–20. [3] Verbrugge MW, Hill RF. Ion and solvent transport in ion-exchange membranes. J Electrochem Soc 1990a; 137: 886 –99. [4] Verbrugge MW, Hill RF. Analysis of promising per2uorosulfonic acid membranes for fuel-cell electrolytes. J Electrochem Soc 1990b; 137:3770 –7. [5] Bernardi DM, Verbrugge MW. Mathematical model of a gas di/usion electrode bonded to a polymer electrolyte. AIChE J 1991;37:1151–63.

999

[6] Bernardi DM, Verbrugge MW. A mathematical model for the solid–polymer–electrode fuel cell. J Electrochem Soc 1992;139:2477–91. [7] Singh D, Lu DM, Djilali N. Numerical analysis of transport processes in proton exchange membrane fuel cells. Proceedings of the First International Energy and Environment Symposium, Karadeniz Technical University, Trabzon, Turkey, 1996. p. 1–34. [8] Amphlett JC, Baumert RM, Mann RF, Peppley BA, Roberge PR. Performance modeling of the Ballard Mark IV solid polymer electrolyte fuel cell. J Electrochem Soc 1995;142: 1–15. [9] Fuller TF, Newman J. Water and thermal management in solid–polymer–electrolyte fuel cells. J Electrochem Soc 1993;140:1218–25. [10] Nguyen TV, White RE. A water and heat management model for proton-exchange-membrane fuel cells. J Electrochem Soc 1993;140:2178–86. [11] Weisbrod KR, Grot SA, Vanderborgh N. Through-theelectrode model of a proton exchange membrane fuel cell. In: Gottesfeld S, editor. Proceedings of the First International Symposium on Proton Conducting Membrane Fuel Cells, vol. 95-23. New Jersey: The Electrochemical Society, 1995. p. 152– 66. [12] Gurau V, Liu HT, Kakac S. Two-dimensional model for proton exchange membrane fuel cells. AIChE J 1998;44:2410–22. [13] Broka K, Ekdunge P. Modelling the PEM fuel cell cathode. J Appl Electrochem 1997;27:281–9. [14] Devers D, Wohr M, Yasuda K et al. Simulation of a polymer electrolyte fuel cell electrode. J Appl Electrochem 1997;27:1254–64. [15] Springer TE, Gottesfeld S. Pseudo-homogeneous catalyst layer model for polymer electrolyte fuel cell. In: White WE et al., editors. Modeling of batteries and fuel cells, vol. 91-10. New Jersey: The Electrochemical Society, 1991. p. 197–208. [16] Eikerling M, Kornyshev AA. Modeling the performance of the cathode catalyst layer of polymer electrolyte fuel cells. J Electroanal Chem 1998;453:89–106. [17] Press WH, Teukolsky SA, Vetterling WT et al. Numerical recipes in FORTRAN, 2nd ed. Cambridge: Cambridge University Press, 1992. p. 745 –78. [18] You LX. The two-phase 2ow, transport mechanism and performance studies for PEM fuel cells. Dissertation, University of Miami, 2001.