Influence of the Design Parameters in a Proton Exchange Membrane

2258

Energy & Fuels 2007, 21, 2258-2267

Influence of the Design Parameters in a Proton Exchange Membrane (PEM) Fuel Cell on the Mechanical Behavior of the Polymer Membrane Maher A. R. Sadiq Al-Baghdadi* and Haroun A. K. Shahad Al-Janabi Department of Mechanical Engineering, International Technological UniVersity, 115 Dollis Hill Lane, London NW2 6HS, United Kingdom ReceiVed NoVember 23, 2006. ReVised Manuscript ReceiVed March 14, 2007

A full three-dimensional, non-isothermal computational fluid dynamics model of a proton exchange membrane (PEM) fuel cell with straight flow field channels has been developed to simulate the hygro and thermal strains in a polymer membrane, which developed during the cell operation. This comprehensive model accounts for the major transport phenomena in a PEM fuel cell: convective and diffusive heat and mass transfer, electrode kinetics, transport and phase-change mechanism of water, and potential fields. The model is shown to understand the many interacting, complex electrochemical, transport phenomena, and stresses distribution that have limited experimental data. This model is used to study the effect of design parameters on the behavior of the polymer membrane. The results of this study helped us identify critical parameters and shed insight into the physical mechanisms leading to a fuel cell performance under various operating conditions. The results show that the nonuniform distribution of stresses, caused by the temperature gradient in the cell, induces localized bending stresses, which can contribute to delaminating between the membrane and the gas diffusion layers. These stresses in the membrane may explain the occurrence of cracks and pinholes in the membrane under steadystate loading during regular cell operation.

The durability of proton exchange membranes (PEMs) used in fuel cells is a major factor in the operating lifetime of fuel cell systems. Durability is a complicated phenomenon, linked to the chemical and mechanical interactions of the fuel cell components, i.e., electrocatalysts, membranes, gas diffusion layers, and bipolar plates, under severe environmental conditions, such as elevated temperature and low humidity.1 In fuel cell systems, failure may occur in several ways, such as chemical degradation of the ionomer membrane or mechanical failure in the PEM, that results in the gradual reduction of ionic conductivity, the increase in the total cell resistance, and the reduction of voltage and loss of output power.2 Mechanical damage in the PEM can appear as through-the-thickness flaws or pinholes in the membrane or delaminating between the polymer membrane and gas diffusion layers.1,2 Water management is one of the critical operation issues in PEM fuel cells.3 Spatially varying concentrations of water in both vapor and liquid form are expected throughout the cell because of varying rates of production and transport. Devising better water management is therefore a key issue in PEM fuel cell design, and this requires improved understanding of the

parameters affecting water transport in the membrane.4,5 Thermal management is also required to remove the heat produced by the electrochemical reaction to prevent drying out of the membrane, which in turn can result not only in reduced performance but also in eventual rupture of the membrane.6,7 Thermal management is also essential for the control of the water evaporation or condensation rates. As a result of the changes in temperature and moisture, the PEM, gas diffusion layer (GDL), and bipolar plates will all experience expansion and contraction. Because of the different thermal expansion and swelling coefficients between these materials, hygro-thermal stresses are expected to be introduced into the unit cell during operation. In addition, the nonuniform current and reactant flow distributions in the cell can result in nonuniform temperature and moisture content of the cell, which could in turn, potentially cause localized increases in the stress magnitudes. The need for an improved lifetime of PEM fuel cells necessitates that the failure mechanisms be clearly understood and life prediction models be developed, so that new designs can be introduced to improve long-term performance. Weber and Newman1 developed a one-dimensional model to study the stresses development in the fuel cell. They showed that hygro-thermal stresses might be an important reason for

* To whom correspondence should be addressed: Mechanical and Energy Department, Higher Institute of Mechanical Engineering, Yefren, Post Office Box 65943, Libya. E-mail: [email protected]. (1) Webber, A.; Newman, J. A. Theoretical study of membrane constraint in polymer-electrolyte fuel cell. AIChE J. 2004, 50 (12), 3215-3226. (2) Tang, Y.; Santare, M. H.; Karlsson, A. M.; Cleghorn, S.; Johnson, W. B. Stresses in proton exchange membranes due to hygro-thermal loading. J. Fuel Cell Sci. Technol. 2006, 3 (5), 119-124. (3) Sui, P.; Djilali, N. Analysis of water transport in proton exchange membranes using a phenomenological model. J. Fuel Cell Sci. Technol. 2005, 2 (1), 149-155.

(4) Berning, T.; Lu, D. M.; Djilali, N. Three-dimensional computational analysis of transport phenomena in a PEM fuel cell. J. Power Sources 2002, 106 (1-2), 284-294. (5) Berning, T.; Djilali, N. Three-dimensional computational analysis of transport phenomenon in a PEM fuel cellsA parametric study. J. Power Sources 2003, 124 (2), 440-452. (6) Hyunchul, J.; Meng, H.; Wang, C. Y. A single-phase, non-isothermal model for PEM fuel cells. Int. J. Heat Mass Transfer 2005, 48 (7), 13031315. (7) Sivertsen, B. R.; Djilali, N. CFD based modelling of proton exchange membrane fuel cells. J. Power Sources 2005, 141 (1), 65-78.

1. Introduction

10.1021/ef060596x CCC: $37.00 © 2007 American Chemical Society Published on Web 05/08/2007

Design Parameters in a PEM Fuel Cell

Energy & Fuels, Vol. 21, No. 4, 2007 2259

membrane failure and that mechanical stresses might be particularly important in systems that are non-isothermal. However, their model is one-dimensional and does not include the effects of material property mismatch among PEM, GDL, and bipolar plates. Tang et al.2 studied the hygro and thermal stresses in the fuel cell caused by step changes of the temperature and relative humidity. The influence of membrane thickness was also studied, which shows a less significant effect. However, their model is two-dimensional, where the hygro-thermal stresses are absent in the third direction (flow direction). In addition, a simplified temperature and humidity profile was assumed (a constant temperature for each upper and lower surface of the membrane was assumed), with no internal heat generation. An operating fuel cell has varying local conditions of temperature, humidity, and power generation (and thereby heat generation) across the active area of the fuel cell in three dimensions. Nevertheless, no models have yet been published to incorporate the effect of hygro-thermal stresses into actual fuel cell models to study the effect of these conditions on the stresses developed in the membrane. In addition, as a result of the architecture of a cell, the transport phenomena in a fuel cell are inherently three-dimensional, but no models have yet been published to address the hygro-thermal stresses in PEM with the three-dimensional effect. The difficult experimental environment of fuel cell systems has stimulated efforts to develop models that could simulate and predict multidimensional coupled transport of reactants, heat, and charged species using computational fluid dynamic (CFD) methods. A comprehensive computational model should include the equations and other numerical relations needed to fully define fuel cell behavior over the range of interest. 2. Model Description The present work presents a comprehensive three-dimensional, multiphase, non-isothermal model of a PEM fuel cell that incorporates the significant physical processes and the key parameters affecting fuel cell performance. The following assumptions are made: (i) the fuel cell operates under steady-state conditions; (ii) the membrane is fully humidified, so that the ionic conductivity is constant; (iii) the membrane is impermeable to gases, and the crossover of reactant gases is neglected; (iv) the gas diffusion layer is homogeneous and isotropic; (v) the gases entering the cell are fully humidified; (vi) no liquid water enters the cell at the inlets; (vii) the produced water is in the liquid phase; (viii) a two-phase flow is inside the porous media; (ix) the liquid and gas phases share the same pressure field inside the flow channels; and (x) both phases occupy a certain local volume fraction inside the porous media, and their interaction is accounted for through a multifluid approach. The model accounts for both gas and liquid phases in the same computational domain and thus allows for the implementation of a phase change inside the gas diffusion layers. The model includes the transport of gaseous species, liquid water, protons, energy, and water dissolved in the ion-conducting polymer. Water transport inside the porous gas diffusion and catalyst layers is described by two physical mechanisms (viscous drag and capillary pressure forces) and is described by advection within the gas channels. Water transport across the membrane is also described by two physical mechanisms (electro-osmotic drag and diffusion). Water is assumed to be exchanged among three phases (liquid, vapor, and dissolved), and equilibrium among these phases is assumed. A unique feature of the present model is to incorporate the effect of hygro and thermal stresses into an actual three-dimensional fuel cell model. The model features an algorithm that allows for a more realistic representation of the local activation overpotentials, which leads to improved prediction of the local current density distribution. This leads to a high accuracy prediction of the temperature distribution in the cell and therefore thermal stresses. This model

Figure 1. (a) Three-dimensional computational domain and (b) computational mesh of a PEM fuel cell.

also takes into account convection and diffusion of different species in the channels as well as in the porous gas diffusion layer, heat transfer in the solids as well as in the gases, and electrochemical reactions. The model reflects the influence of numerous parameters on the fuel cell performance, including geometry, materials, operating, and others, to investigate the in situ stresses in polymer membranes. The present multiphase model is capable of identifying important parameters for the wetting behavior of the gas diffusion layers and can be used to identify conditions that might lead to the onset of pore plugging, which has a detrimental effect on the fuel cell performance. This model is used to study the effect of design parameters on the behavior of the polymer membrane. 2.1. Computational Domain. A computational model of an entire cell would require very large computing resources and excessively long simulation times. The computational domain in this study is therefore limited to one straight flow channel with the land areas. The full computational domain consists of cathode and anode gas-flow channels and the membrane electrode assembly as shown in Figure 1. 2.2. Model Equations. 2.2.1. Gas-Flow Channels. In the fuel cell channels, the gas-flow field is obtained by solving the steadystate Navier-Stokes equations, i.e., the continuity equation, the mass conservation equation for each phase, yielding the volume fraction (r), and along with the momentum equations for the

2260 Energy & Fuels, Vol. 21, No. 4, 2007

Al-Baghdadi and Al-Janabi

pressure distribution inside the channels. The continuity equation for the gas phase inside the channel is given by ∇(rgFgug) ) 0

(1)

and for the liquid phase inside the channel becomes ∇(rlFlul) ) 0

(2)

Two sets of momentum equations are solved in the channels, and they share the same pressure field. Under these conditions, it can be shown that the momentum equations become 2 ∇(FgugXug - µg∇ug) ) -∇rg P + µg∇ug + ∇[µg(∇ug)T] 3

(

)

(3)

and 2 ∇(FlulXul - µl∇ul) ) -∇rl P + µl∇ul + ∇[µl(∇ul)T] 3

(

)

[

N

M

∑D M ij

j)1

j

[(

∇yj + yj

)

∇M M

]

∇P

+ (xj - yj)

P

]

∇T T

) 0 (5)

where the subscript i denotes oxygen at the cathode side and hydrogen at the anode side and j is the water vapor in both cases. Nitrogen is the third species at the cathode side. The Maxwell-Stefan diffusion coefficients of any two species are dependent upon temperature and pressure. They can be calculated according to the empirical relation based on the kinetic gas theory8 T1.75 × 10-3 P[(

∑V )

1/3

ki

k

+(

∑V )

[

1/3 2

kj

]

1

+

Mi

1 Mj

]

1/2

(6)

Two liquid water transport mechanisms are considered: shear, which drags the liquid phase along with the gas phase in the direction of the pressure gradient, and capillary forces, which drive liquid water from high to low saturation regions.9 Therefore, the momentum equation for the liquid phase inside the gas diffusion layer becomes ul ) -

Kpl Kpl ∂Pc ∇sat ∇P + µl µl ∂sat

(11)

The functional variation of the capillary pressure with saturation is prescribed following Leverett,9 who has shown that

()

Pc ) σ

Kp

1/2

(1.417(1 - sat) - 2.12(1 - sat)2 + 1.263(1 - sat)3) (12)

The liquid phase consists of pure water, while the gas phase has multicomponents. The transport of each species in the gas phase is governed by a general convection-diffusion equation in conjunction with the Stefan-Maxwell equations to account for multispecies diffusion

[

∇ -(1 - sat)Fgyi

N

M

∑D M j)1

[(

∇yj + yj

ij

j

In this equation, the pressure is in atm and the binary diffusion coefficient is in cm2 s-1. The values for (∑Vki) are given by Fuller et al.8 The temperature field is obtained by solving the convective energy equation (7)

The gas and liquid phases are assumed to be in thermodynamic equilibrium; hence, the temperature of the liquid water is the same as the gas-phase temperature. 2.2.2. Gas Diffusion Layers. The physics of multiple phases through a porous medium is further complicated here with the phase change and the sources and sinks associated with the electrochemical reaction. The equations used to describe transport in the gas diffusion layers are given below. Mass transfer in the form of evaporation (m˘ phase > 0) and condensation (m˘ phase < 0) is assumed, so that the mass balance equations for both phases are ∇((1 - sat)Fgug) ) m˘ phase

(8)

∇(satFlul) ) m˘ phase

(9)

and

)

∇M M

]

∇P

+ (xj - yj)

(1 - sat)Fgyiug + DTi

]

∇T T

P

+

) m˘ phase (13)

To account for geometric constraints of the porous media, the diffusivities are corrected using the Bruggemann correction formula10 1.5 Deff ij ) Dij

k

∇(rg(FgCpgugT - kg∇T)) ) 0

(10)

+

rgFgyiug + DTi

Dij )

Kp ug ) -(1 - sat) ∇P µg

(4)

The mass balance is described by the divergence of the mass flux through diffusion and convection. Multiple species are considered in the gas phase only, and the species conservation equation in multicomponent, multiphase flow can be written in the following expression for species i: ∇ -rgFgyi

The momentum equation for the gas phase reduces to Darcy’s law, which is, however, based on the relative permeability for the gas phase (Kp). The relative permeability accounts for the reduction in the pore space available for one phase because of the existence of the second phase.9 The momentum equation for the gas phase inside the gas diffusion layer becomes

(14)

The heat transfer in the gas diffusion layers is governed by the energy equation as follows: ∇((1 - sat)(FgCpgugT - keff,g∇T)) ) β(Tsolid - T) - m˘ phase∆Hevap (15) where the term (β(Tsolid - T)), on the right-hand side, accounts for the heat exchange to and from the solid matrix of the GDL. The gas and liquid phases are assumed to be in thermodynamic equilibrium; i.e., the liquid water and the gas phase are at the same temperature. The potential distribution in the gas diffusion layers is governed by ∇(λe∇φ) ) 0

(16)

To account for the magnitude of the phase change inside the GDL, expressions are required to relate the level of over- and undersatu(8) Fuller, E. N.; Schettler, P. D.; Giddings, J. C. A new methode for prediction of binary gas-phase diffusion coefficients. Ind. Eng. Chem. 1966, 58 (5), 18-27. (9) Berning, T.; Djilali, N. A 3D, multi-phase, multicomponent model of the cathode and anode of a PEM fuel cell. J. Electrochem. Soc. 2003, 150 (12), A1589-A1598.



Table 1. Electrode and Membrane Parameters for Base-Case Operating Conditions parameter

symbol

value

electrode porosity electrode electronic conductivity membrane ionic conductivity (humidified Nafion 117) transfer coefficient, anode side transfer coefficient, cathode side cathode reference exchange current density anode reference exchange current density electrode thermal conductivity membrane thermal conductivity electrode hydraulic permeability entropy change of the cathode side reaction heat-transfer coefficient between the solid and gas phases protonic diffusion coefficient fixed-charge concentration fixed-site charge electro-osmotic drag coefficient droplet diameter condensation constant scaling parameter for evaporation electrode Poisson’s ratio membrane Poisson’s ratio electrode thermal expansion membrane thermal expansion electrode Young’s modulus membrane Young’s modulus electrode density membrane density membrane humidity swelling-expansion tensor

λe λm Ra Rc iref o,c iref o,a keff kmem kp ∆S β DH+ cf zf nd Ddrop C ω j νGDL νmem PGDL Pmem ΨGDL Ψmem FGDL Fmem λhmem

0.4 100 17.1223 0.5 1 1.8081 × 10-3 2465.598 1.3 0.455 1.76 × 10-11 -326.36 4 × 106 4.5 × 10-9 1200 -1 2.5 1.0 × 10-8 1.0 × 10-5 0.01 0.25 0.25 -0.8 × 10-6 123 × 10-6 1 × 1010 249 × 106 400 2000 23 × 10-4

ration as well as the amount of liquid water present to the amount of water undergoing a phase change. In the present work, the procedure of Berning and Djilali9 was used to account for the magnitude of the phase change inside the GDL. 2.2.3. Catalyst Layers. The catalyst layer is treated as a thin interface, where sink and source terms for the reactants are implemented. Because of the infinitesimal thickness, the source terms are actually implemented in the last grid cell of the porous medium. At the cathode side, the sink term for oxygen is given by MO2 i SO2 ) 4F c

(17)

Whereas the sink term for hydrogen is specified as MH2 i 2F a

SH2 ) -

(18)

The production of water is modeled as a source term and hence can be written as SH2O )

MH2O i 2F c

[

]

T(-∆s) + ηact i neF

(20)

The local current density distribution in the catalyst layers is modeled by the Butler-Volmer equation12,13

( )[ ( )[

ic ) iref o,c

ia ) iref o,a

CO2

exp

COref2

CH 2

CHref2

(

1/2

exp

)

(

)

(

RcF RaF ηact,c + exp η RT RT act,c

(

reference 4 7 10 4 13 5 5 10 10 12 11 4 4 4 4 10 9 9 9 2 2 2 17 2 17 2 2 2

S m-1 S m-1 A m-2 A m-2 W m-1 K-1 W m-1 K-1 m2 J mol-1 K-1 W m-3 m2 s-1 mol m-3 m

K-1 K-1 Pa Pa kg m-3 kg m-3 %-1

Table 2. Geometrical and Operational Parameters for Base-Case Conditions parameter

symbol

channel length channel width channel height land area width gas diffusion layer thickness wet membrane thickness (Nafion 117) catalyst layer thickness hydrogen reference mole fraction oxygen reference mole fraction anode pressure cathode pressure inlet fuel and air temperature relative humidity of inlet fuel and air (fully humidified conditions) air stoichiometric flow ratio fuel stoichiometric flow ratio

L W H Wland δGDL δmem

0.05 1 × 10-3 1 × 10-3 1 × 10-3 0.26 × 10-3 0.23 × 10-3

value

unit m m m m m m

δCL

m

xOref2 Pa Pc Tcell ψ

0.0287 × 10-3 xHref20.846 39 0.177 74 3 3 353.15 100

ξc ξa

2 2

atm atm K %

2.2.4. Membrane. The balance between the electro-osmotic drag of water from anode to cathode and back diffusion from cathode to anode yields the net water flux through the membrane

(19)

The generation of heat in the cell is due to entropy changes as well as the activation overpotential of irreversibility11 q˘ )

unit

)]

RcF RaF ηact,a + exp η RT RT act,a

)]

i NW ) ndMH2O - ∇(FDW∇cW) F

(23)

The water diffusivity in the polymer can be calculated as follow:14

[ (3031 - T1)]

DW ) 1.3 × 10-10 exp 2416

(24)

The variable cW represents the number of water molecules per sulfonic acid group (i.e., mol of H2O/equiv of SO3-1). The water content in the electrolyte phase is related to the water vapor activity via15,16

(21) cW ) 0.043 + 17.81a - 39.85a2 + 36.0a3 (0 < a e 1) (22)

cW ) 14.0 + 1.4(a - 1) cW ) 16.8

(1 < a e 3)

(a g 3)

(25)


Al-Baghdadi and Al-Janabi πmem ) πM + πT + πS

(29)

where πM is the contribution from the mechanical forces and πT and πS are the thermal- and swelling-induced strains, respectively. The thermal strains resulting from a change in temperature of an unconstrained isotropic volume are given by πT ) Pmem(T - TRe f)

(30)

Similarly, the swelling strains caused by moisture uptake are given by πS ) λhmem(R - RRe f)

Figure 2. Comparison of the model and the experimental polarization curves.

The water vapor activity given by a)

xWP Psat

(26)

(27)

The potential loss in the membrane is due to the resistance to proton transport across the membrane and is governed by ∇(λm∇φ) ) 0

The mechanical boundary conditions are noted in Figure 1a. The initial conditions corresponding to a zero stress state are defined; all components of the cell stack are set to the reference temperature of 20 °C and relative humidity of 35% (corresponding to the assembly conditions).2,17 In addition, a constant pressure of 1 MPa is applied on the surface of the lower graphite plate, corresponding to a case where the fuel cell stack is equipped with springs to control the clamping force.

3. Results and Discussion

Heat transfer in the membrane is governed by ∇(kmem∇T) ) 0

(31)

(28)

2.2.5. Strain Tensor in a Proton Exchange Membrane. Using the hygro-thermoelasticity theory, the effects of the temperature and moisture as well as the mechanical forces on the behavior of elastic bodies have been addressed. In the present work, the total strain tensor is determined using the same expression by Tang et al.2

Figure 3. Strain distribution inside the membrane at base-case conditions.

The governing equations were discretized using a finite volume method and solved using a general-purpose computational fluid dynamic code. Stringent numerical tests were performed to ensure that the solutions were independent of the grid size. A computational mesh consisting of about 150 000 elements was found to provide sufficient spatial resolution. The values of the electrochemical transport parameters for the basecase operating conditions are listed in Table 1. The geometric and base-case operating conditions are listed in Table 2. It is important to note that because this model accounts for all major transport processes and the modeling domain comprises all of the elements of a complete cell, no parameters needed to be adjusted to obtain physical results. The multiphase model is validated by comparing model results to experimental data provided by Wang et al.,12 as seen in Figure 2. The importance of the phase change to the accurate modeling of the fuel cell performance is illustrated. Performance curves with and without a phase change as well as experimental data are shown in Figure



Figure 4. Total displacement [contour (µ, in m)] and deformed shape plot (scale enlarged 300 times) for the membrane at base-case conditions on the y-z plane at x ) 10 mm.

2 for the base-case operating conditions. A comparison of the two curves demonstrates that the effects of liquid water accumulation become apparent even at relatively low values of current density. Furthermore, when liquid water effects are not included in the model, the cell voltage does not exhibit an increasingly steep drop as the cell approaches its limiting current density. This drop off in performance is clearly demonstrated by experimental data but cannot be accurately modeled without the incorporation of a phase change. When the effects of a phase change are included, the current model is able to more closely simulate performance, especially in the region where masstransport effects begin to dominate. The performance characteristics of the fuel cell based on a certain parameter can be obtained by varying that parameter while keeping all other parameters constant at base-case conditions. Results obtained from these parametric studies will allow us to identify the critical parameters for fuel cell performance. The behavior of the polymer membrane at various design conditions is compared using the total displacement and deformation shape of the membrane for each case. Results with different design conditions for the cell operate at a nominal current density of 1.2 A cm-2 are discussed in the following subsections. 3.1. Base-Case Operating Conditions. Strain tensor distribution inside the membrane during the cell operating at base-case conditions can be shown in Figure 3. It can be seen that the maximum strain occurs where the temperature is highest, which is near the cathode-side inlet area (from x ) 0 to 20 mm). The maximum strain appears in the lower surface of the membrane (10) Nguyen, P. T.; Berning, T.; Djilali, N. Computational model of a PEM fuel cell with serpentine gas flow channels. J. Power Sources 2004, 130 (1-2), 149-157. (11) Lampinen, M. J.; Fomino, M. Analysis of free energy and entropy changes for half-cell reactions. J. Electrochem. Soc. 1993, 140 (12), 35373546. (12) Wang, L.; Husar, A.; Zhou, T.; Liu, H. A parametric study of PEM fuel cell performances. Int. J. Hydrogen Energy 2003, 28 (11), 1263-1272. (13) Parthasarathy, A.; Srinivasan, S.; Appleby, J. A.; Martin, C. R. Pressure dependence of the oxygen reduction reaction at the platinum microelectrode/nafion interface: Electrode kinetics and mass transport. J. Electrochem. Soc. 1992, 139 (10), 2856-2862. (14) Siegel, N. P.; Ellis, M. W.; Nelson, D. J.; von Spakovsky M. R. A two-dimensional computational model of a PEMFC with liquid water transport. J. Power Sources 2004, 128 (2), 173-184.

Figure 5. Total displacement [contour (µ, in m)] and deformed shape plot (scale enlarged 300 times) for the membrane on the y-z plane at x ) 10 mm for two different GDL porosities: (a) 0.3 and (b) 0.5.

(cathode side), implying that major heat generation takes place near this region. Figure 4 shows the total displacement and deformation shape for the membrane. The figure illustrates the effect of strain on the membrane. Because of the different thermal expansion and swelling coefficients between gas diffusion layers and membrane materials with nonuniform temperature distributions in the cell during operation, hygro and thermal strains and deformation are introduced. The nonuniform distribution of strain, caused by the temperature gradient in the membrane electrode assembly (MEA), induces localized bending stresses, which can contribute to delamination between the membrane and GDLs. It can be (15) Hu, M.; Gu, A.; Wang, M.; Zhu, X.; Yu, L. Three dimensional, two phase flow mathematical model for PEM fuel cell. Part I. Model development. Energy ConVers. Manage. 2004, 45 (11-12), 1861-1882. (16) Hu, M.; Gu, A.; Wang, M.; Zhu, X.; Yu, L. Three dimensional, two phase flow mathematical model for PEM fuel cell. Part II. Analysis and discussion of the internal transport mechanisms. Energy ConVers. Manage. 2004, 45 (11-12), 1883-1916. (17) Product information. DuPont Nafion PFSA membranes N-112, NE1135, N-115, N-117, and NE-1110 perfluorosulfonic acid polymer. NAE101, 2005.


Figure 6. Total displacement [contour (µ, in m)] and deformed shape plot (scale enlarged 300 times) for the membrane on the y-z plane at x ) 10 mm for two different values of the GDL thermal conductivity: (a) 0.5 W m-1 K-1 and (b) 2.9 W m-1 K-1.

seen that the total displacement and degree of the deformation in the membrane are directly related to the temperature. In addition, the deformation that occurs in the membrane under the land areas is much smaller than under the channel areas because of the clamping force effect. 3.2. Effect of GDL Porosity. The porosity of the gas diffusion layer has two comparing effects on the fuel cell performance. Because the porous region provides the space for the reactants to diffuse toward the catalyst region, an increase in the porosity means that the onset of mass-transport limitations occurs at higher current densities; i.e., it leads to higher limiting currents. The adverse effect of a high porosity is an expected increase in the contact resistance. A higher gas diffusion layer porosity improves the mass transport within the cell, and this leads to the reduction of the mass-transport loss. The molar oxygen fraction at the catalyst layer increases with a more even distribution with an increase in the porosity. This is because a higher value of the porosity provides less resistance for oxygen to reach the catalyst layer. A higher porosity evens out the local


Figure 7. Total displacement [contour (µ, in m)] and deformed shape plot (scale enlarged 300 times) for the membrane on the y-z plane at x ) 10 mm for two different values of the channel width: (a) 0.8 mm and (b) 1.2 mm.

current density distribution. For a lower value of the porosity, a much higher fraction of the total current is generated under the channel area. This can lead to local hot spots inside the membrane electrode assembly. These hot spots can lead to a further drying out of the membrane, thus increasing the electric resistance, which in turn leads to more heat generation and therefore an increase in the total displacement and degree of the deformation inside the membrane, as seen in Figure 5. Thus, it is important to keep the current density relatively even throughout the cell. As mentioned above, another loss mechanism that is important when considering different gas diffusion layer porosities is the contact resistance. Contact resistance occurs at all interfaces inside the fuel cell. The magnitude of the contact resistance depends upon various parameters, such as the surface material and treatment and the applied stack pressure. The electrode porosity has a negative effect on electron conduction. Because the solid matrix of the gas diffusion layer provides the pathways for electron transport, the higher volume porosity increases resistance to electron transport in the gas diffusion layers.



Figure 8. Total displacement [contour (µ, in m)] and deformed shape plot (scale enlarged 300 times) for the membrane on the y-z plane at x ) 10 mm for two different GDL thicknesses: (a) 0.2 mm and (b) 0.3 mm.

Figure 9. Total displacement [contour (µ, in m)] and deformed shape plot (scale enlarged 300 times) for the membrane on the y-z plane at x ) 10 mm for two different membrane thicknesses: (a) 0.2 mm and (b) 0.26 mm.

3.3. Effect of GDL Thermal Conductivity. Thermal management is required to remove the heat produced by the electrochemical reaction to prevent drying out of the membrane and excessive thermal stresses that may result in the rupture of the membrane. The small temperature differential between the fuel cell stack and the operating environment make thermal management a challenging problem in PEM fuel cells. The total displacement and degree of the deformation inside the membrane can be shown in Figure 6 for two cases of the gas diffusion layer thermal conductivity values. The maximum temperature with a higher gradient appears in the cathode side catalyst layer of the lower thermal conductivity. Heat generated in the catalyst layer is primarily removed through the gas diffusion layer to the current collector rib by lateral conduction. This process is controlled by the gas diffusion layer thermal conductivity. Therefore, the membrane temperature is strongly influenced by the gas diffusion layer thermal conductivity, indicating a significant role played by lateral heat conduction through the gas diffusion layer in the removal of waste heat to the ambient environment. Therefore, a gas diffusion layer

material having higher thermal conductivity is strongly recommended for fuel cells designed to operate with high power. 3.4. Effect of the Width of Gas Channels. A reduction in the land area width by increasing the width of the gas-flow channel enhances the mass transport of the reactions to the catalyst layer that lies under the land area. The result is an increase in the molar oxygen fraction at the catalyst layer with a more even distribution. It is expected that this will mainly affect the limiting current density and, to a lesser degree, the voltage drop because of mass-transport limitations. The channel width has a large impact on the local current density distribution. For the narrow channel, the local current density can exceed more than 40% of the nominal current density with a sharp dropoff under the land area, where the local current density is about 40% lower than the nominal current density. The wider channel makes for a much more evenly distributed current throughout the cell. However, the temperature peak appears in the cathodeside catalyst layer of the wider channel case, implying that major heat generation takes place where the electrochemical activity is highest, and this leads to an increase in the total displacement



Figure 10. Cell power density with corresponding efficiency and maximum displacement in the membrane for different operating conditions.

and degree of the deformation inside the membrane, as seen in Figure 7. This is because the increase in the width of the gasflow channel means that the velocity of the incoming gas has to be decreased, with all remaining parameters remaining constant, and this will decrease the velocity of the gases in the gas diffusion layer and, hence, reduce the convection heat transfer in this region. Finally, a reduced width of the land area increases the contact resistance between the bipolar plates and the gas diffusion electrodes. Because this is an ohmic loss, it is directly correlated to the land area width. 3.5. Effect of GDL Thickness. The effect of the gas diffusion layer thickness on the fuel cell performance is again mostly on the mass transport, because the ohmic losses of the electrons inside the gas diffusion layer are relatively small as a result of the high conductivity of the carbon-fiber paper. A thinner gas diffusion layer increases the mass transport through it, and this leads to a reduction in the mass-transport loss. The molar oxygen fraction at the catalyst layer increases with a decrease of the gas diffusion layer thickness because of the reduced resistance to oxygen diffusion by the thinner layer. The distribution of the local current density of the cathode side directly depends upon the oxygen concentration. The thicker gas diffusion layer results in a more even distribution of the local current density because of the more even distribution of the molar oxygen fraction at the catalyst layer. Therefore, the maximum temperature gradient appears in the cathode-side catalyst layer of the thinner gas diffusion layer case, and this leads to an increase in the total displacement and degree of the deformation inside the membrane, as seen in Figure 8. 3.6. Effect of the Membrane Thickness. The effect of the membrane thickness on the fuel cell performance is mostly on the resistance of the proton transport across the membrane. The potential loss in the membrane is due to the resistance to proton transport across the membrane from the anode catalyst layer to the cathode catalyst layer. Therefore, a reduction in the membrane thickness means that the path traveled by the protons will be decreased, thereby reducing the membrane resistance, and this leads to a reduction in the potential loss in the

membrane, which in turn leads to less heat generation in the membrane and, as a result, a reduction in the total displacement and degree of the deformation inside the membrane, as seen in Figure 9. These results suggested that reducing the membrane thickness played a significant role in promoting cell performance. However, there is a limitation to this reduction, because of the effect of increased gas crossover with very thin membranes. 3.7. Cell Performance. The effects of the design parameters on the cell operation have been discussed and studied in sections 3.1-3.6 and can be summarized in this section. Figure 10 shows the fuel cell power density with the corresponding efficiency and maximum displacement in the membrane for different operating conditions for the cell to operate at a nominal current density of 1.2 A cm-2. Operating at lower or higher gas diffusion layer porosities from the basecase value (cases 2 and 3) will decrease the cell power and efficiency, but with a higher porosity value (case 3), the total displacement in membrane will decrease. Operating at a higher gas diffusion layer thermal conductivity (case 5) will decrease the maximum temperature gradient inside the cell, and this leads to a reduction in the total displacement in the membrane. The cell power and efficiency remain constant for both cases of higher and lower values of thermal conductivity. Operating at a narrow gas-flow channel (case 6) will increase the cell power and efficiency and lower the cost. In addition, the total displacement in the membrane will decrease. However, the pressure drop inside the cell will increase with a narrow gasflow channel. Operating with a thinner gas diffusion layer (case 8) will slightly increase the cell power and efficiency and lower the cost. However, the total displacement in the membrane will increase. Operating with a thinner membrane (case 10) will increase the cell power and efficiency and lower the cost. In addition, the total displacement in the membrane will decrease. Conclusions A full three-dimensional, multiphase computational fluid dynamics model of a PEM fuel cell with straight flow channels


has been developed to simulate the hygro and thermal strains in the polymer membrane, which developed during the cell operation. This comprehensive model accounts for the major transport phenomena in a PEM fuel cell: convective and diffusive heat and mass transfer, electrode kinetics, transport and phase-change mechanism of water, and potential fields. A parametric study using this model has been performed. The study quantifies the impact of design parameters on the behavior of the polymer membrane. The model is shown to be able to (1) understand the many interacting, complex electrochemical, transport phenomena, and stresses, which cannot be studied experimentally, (2) identify limiting steps and components, and (3) provide a computer-aided tool for the design and optimization of future fuel cells to improve their lifetime with a much higher power density and lower cost. In addition, the results show that the model is capable of identifying important parameters for the wetting behavior of the gas diffusion layers and can be used to identify conditions that might lead to the onset of pore plugging, which has a detrimental effect on the fuel cell performance, especially in the mass-transport-limited region. Acknowledgment. This project was supported by International Technological University (ITU), London, U.K.

Nomenclature a ) water activity CH2 ) local hydrogen concentration (mol m-3) CH2ref ) reference hydrogen concentration (mol m-3) CO2 ) local oxygen concentration (mol m-3) CO2ref ) reference oxygen concentration (mol m-3) Cp ) specific heat capacity (J kg-1 K-1) D ) diffusion coefficient (m2 s-1) F ) Faraday’s constant ) 96487 (C mol-1)

Energy & Fuels, Vol. 21, No. 4, 2007 2267 i ) local current density (A m-2) iref o,a ) anode reference exchange current density iref o,c ) cathode reference exchange current density keff ) effective electrode thermal conductivity (W m-1 K-1) kmem ) membrane thermal conductivity (W m-1 K-1) kp ) hydraulic permeability (m2) m˘ phase ) mass transfer: evaporation (m˘ phase ) m˘ evap) and condensation (m˘ phase ) m˘ cond) (kg s-1) M ) molecular weight (kg mol-1) NW ) net water flux across the membrane (kg m-2 s-1) nd ) electro-osmotic drag coefficient ne ) number of electrons transfer P ) pressure (Pa) Pc ) capillary pressure (Pa) q˘ ) heat generation (W m-2) R ) universal gas constant ) 8.314 (J mol-1 K-1) s ) specific entropy (J mol-1 K-1) sat ) saturation T ) temperature (K) u ) velocity vector (m s-1) xi ) mole fraction yi ) mass fraction Ra ) charge transfer coefficient, anode side Rc ) charge transfer coefficient, cathode side ) porosity λe ) electrode electronic conductivity (S m-1) λm ) membrane ionic conductivity (S m-1) µ ) viscosity (kg m-1 s-1) F ) density (kg m-3) σ ) surface tension (N m-1) π ) strain Pmem ) membrane thermal expansion (K-1) λhmem ) membrane humidity swelling-expansion tensor (%-1) R ) relative humidity (%) EF060596X

Influence of the Design Parameters in a Proton Exchange Membrane

Recommend Documents