- Email: [email protected]
A novel membrane transport model for polymer electrolyte fuel cell simulations
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y 3 9 ( 2 0 1 4 ) 7 0 7 7 e7 0 8 8
Available online at www.sciencedirect.com
ScienceDirect journal homepage: www.elsevier.com/locate/he
A novel membrane transport model for polymer electrolyte fuel cell simulations L. Karpenko-Jereb a,*, P. Innerwinkler a, A.-M. Kelterer a, C. Sternig a, C. Fink b, P. Prenninger b, R. Tatschl b a b
Institute of Physical and Theoretical Chemistry, Graz University of Technology, Graz, Austria AVL List GmbH, Graz, Austria
article info
abstract
Article history:
This work presents the development of a 1D model describing water and charge transport
Received 29 November 2013
through the polymer electrolyte membrane (PEM) in the fuel cell. The considered driving
Received in revised form
forces are electrical potential, concentration and pressure gradients. The membrane
12 February 2014
properties such as water diffusion and electro-osmotic coefficients, water sorption and
Accepted 17 February 2014
ionic conductivity are treated as temperature dependent functions. The dependencies of
Available online 25 March 2014
diffusion and electro-osmotic coefficients on the membrane water concentration are described by linear functions. The membrane conductivity is computed in the framework
Keywords:
of the percolation theory under consideration that the conducting phase in the PEM is
Polymer electrolyte membrane fuel
formed by a hydrated functional groups and absorbed water. This developed membrane
cell
model was implemented in the CFD code AVL FIRE using 1D/3D coupling. The simulated
Membrane transport model
polarization curves at various humidification of the cathode are found in good agreement
CFD code AVL FIRE
with the experiments thus confirming the validity of the model.
Percolation
Copyright ª 2014, Hydrogen Energy Publications, LLC. Published by Elsevier Ltd. All rights reserved.
Water sorption isotherm Polarization curve
Introduction Currently, hydrogen Polymer Electrolyte Membrane Fuel Cells (PEMFCs) are one of the most promising devices as high efficient sources of electric power for electric vehicles. The main advantages of the fuel cells are their high energy conversion efficiency and low noise, as well as their high power and energy density even under consideration of the H2-tank [1e4]. The main goal of the PEMFC modeling is to calculate the cell potential and its efficiency. For this purpose it is necessary to find out which processes contribute to the potential losses and the extent of their contribution. If Ohmic voltage drops across
GDLs and bipolar plates are neglected, the PEMFC potential is expressed by the equation: Ucell ¼ Uoc hc ha hmem
(1)
Ucell e cell potential; Uoc e open-circuit voltage; hc, ha e potential losses at cathode and anode; hmem e membrane overpotential. The main purpose of a membrane transport model is to calculate the membrane over-potential. Fig. 1 depicts the most important processes in the proton exchange membrane coated with a catalyst layer (PEMCC) taking place in a fuel cell at operation. These are the proton and electro-osmotic transports, the water generation resulted by the electro-
* Corresponding author. Tel.: þ43 316 873 32241; fax: þ43 316 873 32202. E-mail address: [email protected] (L. Karpenko-Jereb). http://dx.doi.org/10.1016/j.ijhydene.2014.02.083 0360-3199/Copyright ª 2014, Hydrogen Energy Publications, LLC. Published by Elsevier Ltd. All rights reserved.
7078
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y 3 9 ( 2 0 1 4 ) 7 0 7 7 e7 0 8 8
chemical reaction, the water diffusion and convection. Water diffusion arises due to a water concentration difference between the cathode side and the anode side. Usually the cathode has a higher water concentration caused by water formation as a result of the chemical reaction between Hþ and O2. The electro-osmosis (often called electro-osmotic flux or water drag) is caused by the conjugated transport of the water molecules with the protons in the electric field: the water molecules mainly move in the protons’ solvation shells. The convective flux of water is caused by a pressure gradient. Theoretical modeling of the transport in an ion-exchange membrane (IEM) is based usually on structural and transport models [2,5,6]. The structural models give a simplified concept of the inner structure of the membranes. Currently, the following structural models of the IEMs are known: clusternetwork [7,8], dusty-fluid [9,10], capillary [11,12], microheterogeneous two-phase [13,14] and percolation [15e17]. The transport models describe the transport of ions, water and gases through the membrane. Table 1 displays a comparison of several membrane transport models published in the scientific papers [18e28]. Based on the governing equations used for calculating membrane overpotential, the membrane transport models can be separated into three main groups: first e StefaneMaxwell equation based models [18,19]; second e NernstePlanck equation based models [20,21]; third - Ohm’s law based models [6,22e28]. Most of the models are known in the third group. From experimental investigations it is known that the temperature and pressure gradients have a weak influence on fuel cell voltage. Therefore most of the membrane transport models only consider two driving forces: the gradients of water concentration and those of electrical potential. Models such as those from of Rowe et al. [21], Senn and Poulikakos [24] additionally take into account the pressure gradient. Table 1 also displays the approaches to the description of the membrane transport properties (water diffusion and electro-osmotic coefficients) and the proton concentration as a function of the membrane water concentration as well as approximations in the calculations of the total water flux through the membrane. Most of the models take into account the dependencies of the water diffusion and electro-osmotic coefficients on the membrane water
concentration. The simplified model of Berg et al. [20] sets the electro-osmotic coefficient (Cdrag) in the perfluorinated membrane equal to unity and considers the dependence of the proton concentration on the membrane hydration level. Usually, the total water flux via the membrane is considered as a sum of diffusion and electro-osmotic flows. Kulikovsky proposed [6,22] a description of the transport phenomena in PEM based on the approximation that the water back diffusion is completely compensated by the electro-osmotic transport. This approach allows calculating the total water flux through the membrane analytically. But this model is not able to predict electrode water flooding effects. The membrane transport models differ also in the formulation of boundary conditions and in the definition of the water concentration at interfaces of PEM/CL, CL/GDL as well as in the calculation of the water concentration profile in the fuel cell. In Berg’s model [20] the water sorption isotherm is used for the calculations of the c water concentrations at the interfaces GDL/CL (Fig. 1), Ca w ,Cw , while in other models [6,22,24,28] it is used for the computation of the water concentrations at the membrane boundaries with CLs, Caw and Ccw . Whereas the influence of temperature on water sorption properties of the membrane is described insufficiently in the models [6,18e33], however the experimental investigations [34e36] showed the temperature effect on the membrane water sorption isotherm, which is especially pronounced for the equilibrium of the membrane with the saturated water vapor and liquid water. The aim of this present work was the development of a membrane transport model for PEMFC simulations, which describes water and charge transports through the polymer electrolyte membrane, and considers the influence of temperature and water concentration on the important membrane properties such as the water sorption isotherm, the water diffusion and electro-osmotic coefficients as well as on ionic conductivity.
Membrane transport model Approximations The developed membrane transport model is based on the approximations summarized below:
Fig. 1 e Transport and chemical phenomena in proton exchange membrane coated by catalyst in an operating PEMFC. These phenomena are considered in the developed membrane transport model.
1. The transport of water and protons is considered only in one dimension. 2. It is a steady-state model. 3. Three driving forces are considered to describe the transport: gradients of electrical potential, concentration and pressure. 4. The model is an isothermal model, i.e. it does not take into account an effect of temperature gradient across the membrane on the water transport. For the calculations in the membrane transport model an average value between the temperatures at the membrane interfaces at the cathode and anode sides is taken. 5. Four types of water motions are considered in the model, as shown in Fig. 1: water formation at the cathode, water diffusion, electro-osmosis and water convection.
7079
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y 3 9 ( 2 0 1 4 ) 7 0 7 7 e7 0 8 8
Table 1 e Comparison of the membrane transport models. Model
Driving forces
1. StefaneMaxwell equation Carnes, Djilali, [18] Baschuk, Li [19] 2. NernstePlanck equation Berg et al. [20] Rowe, Li [21] 3. Ohm’s law Kulikovsky A. [6,22] Bao et al. [23] Senn, Poulikakos [24] Gurau, Mann [25] Wu et al. [26] Sui et al. [27] del Real et al. [28]
Approximations
D4
DCw
Dp
Dw
Cþ
Cdrag
jm w
þ þ
þ þ
f(Cw) f(Cw)
Const Const
f(Cw) f(Cw)
s0 s0
þ þ
þ þ
þ
f(Cw) f(Cw)
f(Cw) Const
¼1 f(Cw)
s0 s0
þ þ þ þ þ þ þ
þ þ þ þ þ þ þ
þ
Const Const Const f(Cw) f(Cw) f(Cw) f(Cw)
Const Const Const Const Const Const Const
f(Cw) f(Cw) f(Cw) f(Cw) f(Cw) f(Cw) f(Cw)
¼0 s0 s0 s0 s0 s0 s0
D4, DCw, Dp e gradients of the electrical potential, water concentration and pressure correspondingly; Dw e water diffusion coefficient; Cþ e proton concentration; Cdrag e electro-osmotic coefficient; jm w e total water flux throughout the membrane; Cw e water concentration in the membrane.
6. The boundary conditions of water transport utilize the modified equations from Berg’s et al. model, 2005 [20]. 7. The membrane over-potential is calculated using Ohm’s Law. 8. Water diffusion and electro-osmotic coefficient, as well as membrane conductivity, are functions of the membrane water concentration Cw and temperature. 9. Water sorption isotherm depends on the temperature. 10. The activation energies of water diffusion, electroosmotic coefficient and membrane conductivity do not depend on membrane water concentration.
At 4 1, if liquid water is formed at both boundaries CL/ PEM, the water concentrations at anode and cathode sites are equal Caw ¼ Ccw . In this case the diffusive flux is zero, but the convective flux is taken into account and consequently the total water flux through the membrane is described by: i pl;a pl;c jm w ¼ Cdrag $ þ Khyd $ F Lmem
The membrane over-potential is calculated based on Ohm’s law and expressed by: hmem ¼
The development of the auxiliary equations was based on analysis of experimental data for perfluorinated membranes of Nafion type.
Governing equations All symbols are explained in the “Nomenclature” at the end of the paper. The water diffusion is expressed by Fick’s Law: Dw
dCw jdif dCw ¼ w ; jdif w ¼ a$Dw $ dz a dz
(2)
(4)
At 4 < 1, water vapor operating condition, the convective flux is negligible and the total water flux through the membrane is found from the sum of diffusive and electro-osmotic fluxes: jm w ¼ a$Dw
dCw i þ Cdrag $ F dz
(7)
Boundary conditions The balance of water fluxes between the polymer electrolyte membrane and catalyst layer is based on the boundary conditions proposed by Berg et al. [20]: at anode catalyst layer a a (8) jm w ¼ a$gH2 O;a Cw Cw at cathode catalyst layer jm w þ
i ¼ a$gH2 O;c Ccw Cc w 2F
(3)
The convective flux is calculated from: pl;a pl;c jconv ¼ Khyd $ w Lmem
i $Lmem smem
(9)
At 4 < 1 the balance of the water fluxes in PEM and catalyst layers are expressed as:
The electro-osmotic flux is expressed by: i jdrag ¼ Cdrag $ w F
(6)
(5)
a$Dw
dCw i þ Cdrag ¼ a$gH2 O;a Caw Ca w F dz
(10)
a$Dw
dCw i i ¼ a$gH2 O;c Ccw Cc þ Cdrag þ w F 2F dz
(11)
g e the mass transfer coefficient of water from the GDL to the catalyst layer is calculated via the equation: gH2 O ¼
DGDL w dn þ 0:5$LCL
(12)
The water activity inside the GDL close to the interface of the membrane Cw , is calculated from the water sorption
7080
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y 3 9 ( 2 0 1 4 ) 7 0 7 7 e7 0 8 8
isotherm, which defines the dependence of membrane water concentration on equilibrium water vapor.
The development of the auxiliary equations The development of the auxiliary equations in order to describe the dependences of membrane properties on the relative humidity or membrane water concentration and temperature was based on the thoroughly analysis of experimental data available in the scientific papers [16,18e29,34,55,62e67]. Most of the experimental investigations on the Cdrag, Dw, smem as a function of Cw, reported in the scientific literature, were performed at 25 C. In order to describe the transport properties of the membrane at other temperatures Arrhenius’ equation was applied with an assumption that the activation energies do not depend on the membrane water concentration.
Water sorption isotherm In order to define the water concentration at the membrane boundaries, knowledge of the equilibrium between the water concentrations in the membrane and catalyst/gas diffusion layers is required. This equilibrium can be described by a water sorption isotherm e the relation between the membrane hydration number (number of water molecules per one functional group) and the equilibrium humidity. The water sorption isotherm by the polymer electrolyte membranes is intensively studied using different experimental methods [34e39]. Special attention is focused on the difference between the equilibriums of the membrane with the saturated water vapor and with liquid water. For the first time Schroeder [40] indicated that gelatin, an amino acid related to natural polyelectrolytes, absorbs significantly less water from saturated water vapor than at contact with liquid water. At present, this phenomenon is well-known as Schroeder’s paradox and currently being intensively discussed in the scientific literature [39,41e46]. The existence of this phenomenon has theoretically been proven [41,42]. However, it is an interesting fact, that several experimental investigations [39,43,44] did not confirm the presence of Schroeder’s paradox. Utilizing the thermodynamic model approaches, Lu [41] as well as Eikerling and Berg [42] showed that the transition of water from gas to liquid phase attains a discontinuous jump in the total water concentration in polymer electrolyte membranes. In the numerical modeling of PEMFCs, this discontinuity of the membrane property can cause difficulties. Therefore, most of the published models do not take into account Schroeder’s paradox and describe the water sorption isotherms by continuous polynomial functions [18e29]. In the present work a mathematical function for the membrane water sorption isotherm was developed based on the analysis of the published experimental investigations [34,35,37,38] and the following assumptions: 1) The maximal membrane water concentration grows with increasing temperature. 2) The appearance of liquid water in the system causes a sharp but continuous rise of the membrane water concentration. The second statement was guided by the following considerations: at convenient operating conditions of the PEMFC,
Fig. 2 e Temperature dependence of the specific water concentration in Nafion membrane equilibrated with liquid water.
water is present in both vapor and condensed (liquid) states. Consequently, both phases contact the membrane simultaneously. Thus, the water at the membrane boundary most likely represents a mixture of vapor with droplets of the liquid water. Based on the theory proposed by Eikerling and Berg [42], we can assume that water droplets will be more easily absorbed by the membrane compared to water vapor. But a few water droplets will not be enough in order to increase the membrane water concentration to the maximal value. Thus, at the occurrence of the liquid water in the real PEMFC most possibly the growth of the water content will transit continuously in the membrane. Fig. 2 displays the dependencies of the membrane water concentration for Nafion membranes, equilibrated with liquid water at various temperatures. These data were reported by Morris [34] and Hinatsu [35]. As can be seen in the figure, the membrane water concentration grows linearly with increasing temperature and the two data sets are found in complete agreement and are well described by the linear function: ¼ 0:138$T 28:31 cmax w
(13)
Fig. 3 presents the experimental data on the water sorption isotherms of Nafion membranes measured in a wide range of water vapor activity [34,35,37,38]. We should note that different methods were applied for the measurements, and different devices for controlling the establishment of the equilibrium between the membrane and water vapor were used. Therefore, the data determined, for example, at 25 C by different research groups do not demonstrate complete agreement. We have tried to find a correlation between the water sorption isotherm shape and temperature from these data. On the one hand, Morris [34] showed that the membrane water sorption is higher at 50 C than at 25 C, but on other hand the data of Hinatsu [35] measured at 80 C lie below the values obtained by Pimenov at 25 C [37] and by Morris at 50 C [34]. Perhaps the functions Cw ¼ f(4) have different trends at different temperatures but currently this trend cannot be
7081
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y 3 9 ( 2 0 1 4 ) 7 0 7 7 e7 0 8 8
accurately formulated from the available experimental data. Therefore, we assume that at 0 < 4 < 0.97 the shape of the water sorption isotherm practically does not depend on the temperature (Fig. 3a), but at 4 0.97 the growth of membrane water concentration increases with the temperature (Fig. 3b). The suggested fitting describing the water sorption isotherm is expressed by the equation: Cw ¼ 1:55 þ 13:71$f 24:37$f2 þ 21:87$f3 f cmax w ;f
(14)
f ðcmax w ; 4Þ
where is a multiplier, which is directly proportional described to the maximal membrane water concentration cmax w by Eq. (13) and the relative humidity 4. As seen from Fig. 3b the maximal water concentrations (at 4 ¼ 1) computed by Eq. (14) are found in quite good agreement with experimental data.
same order of magnitude and lie in a range from 0.39 1010 to 5.80 1010 m2/s. The water diffusion coefficients obtained by Rivin and Zawodzinski grow directly with the rising hydration number of the membrane. The water diffusion coefficient reported by Jayakody agrees with data measured by Zawodzinski, and the value of Volino is compared with the data published by Rivin. In the presented model a linear fit (Eq. (15)) of the experimental relationship of Dw vs. the membrane hydration number is suggested, which starts at zero and represents the average of the presented experimental points. Dw ¼ dTw1 $Cw
(15)
where dTw1 is the slope of the straight line fitting, where dTw1 ¼ 2:65 1011 m2 =s at T ¼ 298.15 K.
Electro-osmotic coefficient vs. membrane water concentration Water diffusion coefficient vs. membrane water concentration Fig. 4 displays the dependencies of the water diffusion coefficient on the hydration number of Nafion membranes measured at 25 C by Volino 1982 [47], Zawodzinski 1991 [48], Rivin 2001 [49], Jayakody 2004 [50], Suresh 2005 [51] and Majsztrik 2008 [52]. As seen in the figure, the values have the
Fig. 5 presents the electro-osmotic coefficients in the perfluorinated membranes as a function of membrane water concentration measured by Luo at 25 C [53], Ise at 27 C [54] and Springer at 30 C [29]. The linear fit was performed on the data reported by Luo [53], because the data evaluation for the development of the model was carried out at 25 C. No electroosmosis is possible in the dry membrane. Therefore the suggested fit passes through the origin of the coordinates. Thus the electro-osmotic coefficient is expressed by the equation: 1 $Cw Cdrag ¼ cTdrag
(16)
1 ¼ 0:11 at T ¼ 298.15 K. for Nafion membrane cTdrag
Dw ⋅ 1010 m2 s 6
4
2
0 0
Fig. 3 e The specific water concentration of Nafion membrane as a function of the equilibrium relative water vapor humidity: a) in the wide range of 4 [ 0L1.0; b) for high 4 [ 0.96L1. The points are experimental data from Refs. [34,35,37,38], lines are the suggested fittings.
4
8
12
16
20
Rivin, 2001
suggested fitting
Suresh, 2005
Volino, 1982
Jayakody, 2004
Majsztrik, 2008
Cw
Zawodzinski, 1991
Fig. 4 e The water diffusion coefficient of Nafion membrane as a function of the specific water concentration, T [ 25 C. The points are experimental data from Refs. [47e52], the straight line is the suggested fitting.
7082
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y 3 9 ( 2 0 1 4 ) 7 0 7 7 e7 0 8 8
Ionic conductivity vs. membrane water concentration As shown in the investigations [15e17] the dependence of membrane conductivity on the water concentration can successfully be described by the equation of the percolation theory. The basic ideas of this theory were formulated by Broadbent and Hammersly in 1957 [55]. According to the percolation theory, after reaching some critical volume fraction of the conducting phase fcr, the conductivity of a composition consisting of dielectric and conductor dramatically begins to grow with the increasing fraction of this phase. This change from dielectric to conductor has been called the percolation transition and the critical content of the conducting phase fcr is a threshold of the percolation. At f > fcr the dependence of the conductivity on the conducting phase fraction is described by the equation: t smem ¼ sT1 $ f fcr
(17)
where s - pre-factor of percolation equation, t e percolation exponent, f e volume fraction of the conducting phase; fcr e critical volume fraction of the conducting phase. The quantities s, fcr and t are the parameters of the percolation equation. Calculations based on the theory of probability as well as experimental investigations have shown that the percolation parameters lie in the intervals f ¼ 0.15 0.03, t ¼ 1.6 0.4 [56e60]. The conducting phase in perfluorinated membranes represents the connected hydrophilic cluster domains consisting of hydrated functional groups and networked absorbed water. The following equation was suggested [16,61] in order to calculate the volume fraction of the conducting phase in ionexchange membranes: f ¼ a$MSO3 þ W
(18)
where a e ion-exchange capacity, i.e. the concentration of acid groups; MSO3 - molecular weight of the acid group and W e water mass fraction in PEM.
Cdrag 2.5 2
Fig. 6 e a) The dependence of the specific conductivity of Nafion membrane on the volume fraction of the conducting phase; b) bi-logarithmic dependence lg(smem) L lg(f L fcr), used for the determination of the exponent s of the percolation theory equation. The points are experimental data [53]; the line on (a) is the suggested fitting.
The derivation of the equation is described in detail in the same references [16,61]. Fig. 6a demonstrates the dependences of Nafion conductivity measured by Luo et al. [53] as a function of the volume fraction of the conducting phase calculated using Eq. (18). From these experimental data the estimated critical volume fraction of the conducting phase is fcr ¼ 0.15 which corresponds well to the theoretical range. The treatment of these data in bi-logarithmic coordinates lg(smem) lg(f fcr) (Fig. 6b) has allowed the determination of the parameters s and t. Taken into account that fcr ¼ 0.15 the obtained values of these parameters are: t ¼ 1.2 and sT1 ¼ 1:16½S=cm at T ¼ 298.15 K (Fig. 6).
Temperature dependencies of the membrane transport characteristics
1.5 1
Since the temperature dependencies of the water diffusion, electro-osmotic coefficients and proton conductivity of the perfluorinated membranes are well described by Arrhenius’ equation as shown in the investigations [29,62e65], the expressions (15e17)" can be transformed in (19e21) # dif E 1 1 (19) DTw2 ¼ dTw1 $Cw $exp a $ T2 T1 R
Luo, 2010 Ise, 1999 Springer, 1991 suggested fitting
0.5 0 0
5
10
15
20
Cw
Fig. 5 e The electro-osmotic coefficient of Nafion membrane as a function of the membrane water concentration. The points are the experimental data from Refs. [29,53,54], the straight line is the suggested fitting.
" 2 1 CTdrag ¼ cTdrag $Cw $exp
# drag Ea 1 1 $ T2 T1 R
t Es 1 1 2 ¼ sT1 $ f fcr $exp a $ sTmem R T2 T1
(20)
(21)
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y 3 9 ( 2 0 1 4 ) 7 0 7 7 e7 0 8 8
7083
Fig. 7 e Temperature dependence of the electro-osmotic coefficient (a) and specific conductivity (b) of Nafion membranes in coordinates: a) ln(Cdrag) L (L1000/T); b) ln(smem) L (L1000/T). The points are experimental data from Refs. [54,66,67], the straight lines are linear trends fitting the experimental data.
Table 2 presents the values of the activation energies of Cdrag and smem, which were determined in this work from the temperature dependencies of the experimental data [54,66,67] using a linear fit of ln(Cdrag,smem) versus (1/T) as shown in Fig. 7. The table also contains the activation energy values of Dw and smem published in Refs. [29,64,65]. As seen from the table the activation energy of the water diffusion coefficient of Nafion membrane reported by Kreuer [64] decreases with the growth of membrane hydration level. Unfortunately, any consideration of the mathematical correlation of this dependence would significantly complicate the numerical solution for the total water flux. For this reason it has been assumed that the activation energies do not depend on water concentration and therefore the average values were taken for the simulations.
calculation of the membrane over-potential is depicted in Fig. 8. In order to determine the total water flux in the membrane and membrane water concentration the following padrag , Esa . Water rameters are needed: dw, cdrag, Edif a , Ea concentrations at the boundaries of GDL/CL are calculated from the water sorption isotherm. The total water flux jm w is
Input T2 , i T1 d wT1 , cdrag , Eadif , Eadrag , Eaσ
C wa* , C wc* Results and discussion Programming of the developed membrane transport model was carried out in Fortran95. The general sequence of the Table 2 e The values of the activation energy of water diffusion coefficient, electro-osmotic coefficient and ionic conductivity of Nafion membrane. Reference Water diffusion coefficient Kreuer K., 1997 [64]
Electro-osmotic coefficient Ren X., Gottesfeld S., 2001 [66] Luo Z. et al., 2010 [54] Ionic membrane conductivity Meier F., Eisenberg G., 2004 [65] Springer T. et al., 1991 [29] Luo Z. et al., 2010 [54] Saito M. et al., 2005 [67] a
Cw
Ea 103, [J/mol]
3 5 10
22.960 19.776 16.690
22 22
6.818a 8.019a
e e 22 21
9.889 10.537 9.980a 8.445a
The values were calculated in the present work from the experimental data (Fig. 7) published in Refs. [53,54,66].
j wm
C wa , C wc , C w ( z ) C wmean
T σ mem 2
η mem Output C w ( z ), j wm ,η mem Fig. 8 e Flow chart of the calculation steps in the developed membrane transport model.
7084
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y 3 9 ( 2 0 1 4 ) 7 0 7 7 e7 0 8 8
Fig. 9 e Computational meshes of a) single channel cell, b) Full coupled cell. Images of the experimental cell: c) cathode bipolar plate and end plate, d) whole assembly.
found iteratively by Newton’s method, and the calculation of the water concentration profile Cw(z) is performed by RungeeKutta’s method. From the average membrane water conthe membrane conductivity smem is centration Cmean w calculated using the equation of the percolation theory and then the membrane over-potential hmem is determined based on Ohm’s Law. The developed PEM Model was implemented in the CFD Code AVL FIRE [68]. The convergence of AVL FIRE with this new PEM Model was carefully tested using a Single Channel Fuel Cell (Fig. 9(a)) in around 300 calculations of FC performance at various operating parameters: p, T, 4, l. The voltageecurrent characteristics simulated for varied p, T and l are presented in Figures S1eS4 (Appendix). The model calculations indicate a visible growth of the fuel cell current density with increasing temperature and outlet pressure. Raising the cell temperatures enhances reaction kinetics and contributes to obtain lower activation losses. Higher pressures mean higher hydrogen and oxygen concentrations, which result in enhanced reaction rates and thus increasing current densities. The change of the anode gas stoichiometric ratio does not lead to a visible effect on the cell performance. The increasing cathode stoichiometric flow ratio results in a light enhancement of the current densities, it is linked to increasing reaction kinetics at the higher amount of oxygen delivered to the cathode catalyst layer. The computational results are found in qualitative agreement with experimental observations [69e71]. Advanced validation of the developed PEM Model was carried out using a Full Coupled Fuel Cell (Fig. 9(b)) at various relative humidities at the cathode inlet. The computational geometry consists of four parts: the first part is the main part of
the co-simulation, containing fuel cell and gas channels with 677,264 computational cells (CCs); the second part with bipolar plates is represented by 994,532 CCs; the third part is the end plates with 36,154 CCs; fourth e the cooling channels divided into 23,408 CCs. In the first part of the simulation cell the mass transport is calculated, in the other three parts mainly the heat transport between the cell elements is computed. More detailed descriptions of the simulation procedure of PEMFC by CFD Code AVL FIRE are found in Refs. [68,72]. The experimental measurements were carried out in the Christian-Doppler (CD) Laboratory for Fuel Cell Systems using a fully automated test rig [73,74]. The test rig is capable of supplying well defined mixtures of hydrogen, nitrogen, oxygen, and air by using mass flow controllers. Bubble humidifiers are utilized for gas humidification, and the humidity levels are checked with downstream humidity sensors. Electrochemical potentials are measured utilizing a USB multiflexer card from National Instruments (NI USB-6218). The fuel cell is composed of carbon bipolar plates with 13 straight parallel flow channels (Fig. 9(c)). The cathode bipolar plate is segmented into 10 individual elements, which are electronically insulated from each other (Fig. 9(d)). The bipolar plates are embedded between two end plates made of stainless steel. Cooling channels are installed in the end plates, which provide a precise temperature control of the fuel cell by applying water as the cooling fluid. The active area of the single fuel cell is 25 cm2. The fuel cell assembly shown in Fig. 9(b) completely corresponds to the fuel cell construction (Fig. 9(d)) used for the experimental measurements. The material characteristics and operating conditions applied in the simulations are presented in Tables 3 and 4 correspondingly.
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y 3 9 ( 2 0 1 4 ) 7 0 7 7 e7 0 8 8
7085
Table 3 e Values of materials parameters applied in the simulations. Parameter GDL Thickness Thermal conductivity Permeability Inplane Through plane Porosity Tortuosity Average contact angle Catalyst layers Thickness Cathode Anode Exchange current density Cathode Anode Transfer coefficient Cathode Anode Membrane Membrane thickness Sulfonic acid group concentration Hydraulic permeability of membrane Activation energy of water diffusion Activation energy of electroosmosis Activation energy of membrane conductivity Constants Faraday constant Universal Gas Constant
Symbol
Value
Unit
lGDL ss
3.00$105 [m] 4.5 [W/mK]
KII Kt ˛ s q
2.33$1012 1.07$1014 0.78 1.5 122
LCL,c LCL,a
1.125$105 [m] 1.125$105 [m]
i0c i0a
3.2.105 1.0.1011
[A/m3] [A/m3]
ke,c ke,a
0.1 0.5
[] []
Lmem a Khyd
35$106 1900.0 1$1011
Edif a drag
[m2] [m2] [] [] [ ]
m mol/m3 mol/ (m Pa s) 19.809$103 J/mol
Ea
7.418$103 J/mol
Esa
9.713$103 J/mol
F R
96485.309 C/mol 8.31451 J/(K mol)
Fig. 10 displays the polarization curves of the fuel cell for varying relative humidity at the cathode inlet, calculated by the CFD code AVL FIRE and measured experimentally. A lower inlet relative humidity leads to a decrease of the membrane water concentration and, hence, to a lower proton conductivity and current density. This performance decrease occurs in both the experiment and simulation. However, at a lower humidity the current density is overestimated by the simulation. For two cell voltage values, which are marked by blue dotted lines on Fig. 10, the distributions of the current density
Table 4 e Operating conditions of PEMFC applied in the simulations. Parameter Inlet temperature, [K] Outlet gas pressure, [bar] At cathode At anode Stoichiometric flow ratio At cathode At anode Relative humidity At cathode At anode Voltage, [V]
Values 343.15 1.0 1.0 2.2 1.5 0.9/0.6/0.3 0.9 0.420e0.855
Fig. 10 e Polarization curves of PEMFC with perfluorinated membrane: simulation results (lines); experimental data (points). The operating conditions la/lc [ 2.2/1.5; T [ 70 C; P [ 101325 Pa. The experimental polarization curves were measured with sulfo-cationic perfluorinated membrane, gas diffusion layer SGL-35BC. The GDL and PEM properties are given in Table 3.
along the gas channel were calculated (Fig. 11). As seen from the figure the dependencies of current densities on the channel length at 4in ¼ 0.9 and 4in ¼ 0.6 are found in quite good agreement with the experimental data. In both the experiment and simulation a lower humidity moves the current density maximum towards the cathode outlet (y ¼ 0 cm). However, in the simulation this happens to a greater extent, especially for 4in ¼ 0.3.
Conclusion In the present work a new 1D membrane transport model was developed. In this model the calculation of the membrane over-potential is based on Ohm’s Law. The boundary conditions describing the equilibrium between the water fluxes through CL and PEM utilized the modified boundary conditions from Berg’s et al. [20] model. The dependencies of water diffusion and electro-osmotic coefficients are considered as a linear function of the membrane water concentration, their temperature functions are described using Arrhenius’ equation. The membrane conductivity is computed from the equation of the percolation theory. The feature of the applied approach that the volume fraction of the conducting phase in the membrane is considered as a sum of the volume fractions of the functional groups and absorbed water, is described in detail in Refs. [16,61]. For the first time, the developed model considers temperature influence on the maximal water concentration in the polymer electrolyte membrane. The jumping growth in the membrane water concentration at the transition from saturated water vapor to liquid water is considered as a continuous process. The developed model was implemented in CFD code AVL FIRE. The validation of the model was carried out using the
7086
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y 3 9 ( 2 0 1 4 ) 7 0 7 7 e7 0 8 8
Acknowledgment The authors would like to thank Prof. Hacker V. e Head of the Christian-Doppler (CD) Laboratory for Fuel Cell Systems at Institute of Chemical Engineering and Environmental Technology, Graz University of Technology, for the experimental data presented here. This work has been financially supported by the Austrian Research Promotion Agency (FFG) and AVL List GmbH in the framework of Bridge1 Program, Project No. 827559 “PEMMODEL” Investigation of PEM (proton exchange membrane) transport model for the 3D CFD AVL-FIRE Simulations) and IV2Splus Program, Project No.835811 “A3 FALCON” (Advanced 3D Fuel Cell Analysis and Condition diagnostics).
Appendix A. Supplementary data Supplementary data related to this article can be found at http://dx.doi.org/10.1016/j.ijhydene.2014.02.083.
Nomenclature
Fig. 11 e The distribution of the current density along the channel from inlet (y [ 0 cm) to outlet (y [ 12 cm). Points are experimental data, lines are simulation results.
experimental data obtained on a segmented single fuel cell. The calculated polarization characteristics of the fuel cell at the different humidifications of the inlet air at the cathode are found in quite good agreement with the experimentally measured curves.
Outlook The presented model has been extended from 1D to 3D. The CFD code AVL FIRE with the implemented model will be commercially available in future releases.
Latin symbols a sulfonic acid group concentration, mol/m3 Cw normalized water concentration in the membrane normalized water concentration in the membrane, cmax w equilibrated with liquid water electro-osmotic coefficient (drag coefficient) in the Cdrag membrane 1 electro-osmotic coefficient in the membrane at cTdrag Cw ¼ 1 and T1 Dw water diffusion coefficient in the membrane, m2/s T1 water diffusion coefficient in the membrane at dw Cw ¼ 1 and T1, m2/s activation energy, J/mol Ea F Faraday constant, (A s)/mol f volume fraction of the conducting phase in the membrane critical volume fraction of the conducting phase fcr i electric current density, A/m2 jw water flux, mol/m2s m total water flux through the membrane, mol/m2s jw hydraulic permeability of the membrane, mol/ Khyd (m Pa s) membrane thickness, m Lmem p pressure, Pa R universal gas constant, J/(mol K) T temperature, K t exponent of the equation of the percolation theory, e U Potential, V y channel length, cm z normal direction to the membrane, m Greek symbols 4 water activity in gas phase, relative humidity water transfer coefficient, m/s gH2 O h potential loss, V
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y 3 9 ( 2 0 1 4 ) 7 0 7 7 e7 0 8 8
l s smem
stoichiometric coefficient pre-factor of the equation of the percolation theory, S/m specific membrane conductivity, S/m
Subscripts and superscripts a anode c cathode m; mem membrane mean related to an average value l related to liquid water CL related to catalyst layer GDL related to gas diffusion layer w; H2O water drag related to electro-osmotic transport of water dif diffusion conv related to convective flux s related to membrane conductivity Abbreviations CFD computational fluid dynamics CL catalyst layer FC fuel cell IEM ion exchange membrane GDL gas diffusion layer PEM polymer electrolyte membrane PEMCC polymer electrolyte membrane coated with catalyst layer PEMFC polymer electrolyte membrane fuel cell RH relative humidity
references
[1] Gou B, Na WK, Diong B. Fuel cells. Modeling, control and applications. CRC Press. Taylor & Francis Group; 2010. [2] Barbir F. PEM fuel cells. Theory and practice. Elsevier Academic Press; 2005. [3] Buechi FN, Inaba M, Schmidt ThJ. Polymer electrolyte fuel cell durability. New York: Springer; 2009. [4] Larminie J, Dicks A. Fuel cell systems explained http://www. amazon.de/Fuel-Systems-Explained-James-Larminie/dp/ 0471490261/ref¼sr_1_1? s¼books&ie¼UTF8&qid¼1374826216&sr¼11&keywords¼PEMþfuelþcell - #Wiley VCH; 2000. [5] Weber AZ, Newman J. Modeling transport in polymerelectrolyte fuel cells. Chem Rev 2004;104:4679e726. [6] Kulikovsky AA. Analytical modelling of fuel cells. Elsevier; 2010. [7] Hsu WY, Gierke TD. Ion transport and clustering in nafion perfluorinated membranes. J Membr Sci 1983;13:307e26. [8] Datye VK, Taylor PL, Hopfinger AJ. Simple model for clustering and ionic transport in ionomer membranes. Macromolecules 1984;17:1704e8. [9] Fimrite J, Struchtrup H, Djilali N. Transport phenomena in polymer electrolyte membranes I. Modeling framework. J Electrochem Soc 2005;152:A1804e14. [10] Thampan T, Malhotra S, Tang H, Datta R. Modeling of conductive transport in proton-exchange membranes for fuel cells. J Electrochem Soc 2000;147:3242e50. [11] Cwirko EH, Carbonell RG. A theoretical analysis of Donnan dialysis across charged porous membranes. J Membr Sci 1990;48:155e79.
7087
[12] Koter S. Transport of simple electrolyte solutions through ion-exchange membranes e the capillary model. J Membr Sci 2002;206:201e15. [13] Gnusin NP, Zabolotzki VI, Nikonenko VV, Meshechkov AI. Development of the principle of the generalize conductivity to describe transport phenomena in disperse systems caused by the different driving forces. Zh Fiz Khim 1980;54:1518e22. [14] Zabolotsky VI, Nikonenko VV. Effect of structural membrane inhomogeneity on transport properties. J Membr Sci 1993;79:181e98. [15] Hsu WY, Berkley JR, Meakin P. Ion percolation and insulatorto-conductor transition in nafion perfluorosulfonic acid membranes. Macromolecules 1980;13:198e200. [16] Berezina NP, Karpenko LV. Percolation effects in ionexchange materials. Colloid J 2000;62:676e84. [17] Sumberova V, Bleha M, Wodzki R. Percolation transition in heterogeneous ion exchange membranes. J App Polym Sci 1992;46:611e7. [18] Carnes B, Djilali N. Analysis of coupled proton and water transport in a PEM fuel cell using the binary friction membrane model. Electrochim Acta 2006;52:1038e52. [19] Baschuk JJ, Li X. Modeling of ion and water transport in the polymer electrolyte membrane of PEM fuel cells. Int J Hydrogen Energy 2010;35:5095e103. [20] Berg P, Promislow K, Pierre J, Stumper J, Wetton B. Water management in PEM fuel cells. J Electrochem Soc 2004;151:A341e53. [21] Rowe A, Li X. Mathematical modeling of proton exchange membrane fuel cells. J Power Sources 2001;102:82e96. [22] Kulikovsky AA. The effect of cathodic water on performance of a polymer electrolyte fuel cell. Electrochim Acta 2004;49:5187e96. [23] Bao C, Ouyang M, Yi B. Analysis of water management in proton exchange membrane fuel cells. Tsinghua Sci Tech 2006;11:54e64. [24] Senn SM, Poulikakos D. Polymer electrolyte fuel cells with porous materials as fluid distributors and comparisons with traditional channeled systems. Trans ASME 2004;126:410e8. [25] Gurau V, Mann Jr JA. A Continuum model for water transport in the ionomer-phase of catalyst coated membranes for PEMFCs. Hindawi Publ Corp Adv Mech Eng; 2010. Article ID 372795. [26] Wu H, Berg P, Li X. Non-isothermal transient modeling of water transport in PEM Fuel Cells. J Power Sources 2007;165:232e43. [27] Sui PC, Kumar S, Djilali N. Advanced computational tools for PEM fuel cell design part 1. Development and base case simulations. J Power Sources 2008;180:410e22. [28] del Real AJ, Arce A, Carlos B. Development and experimental validation of a PEM fuel cell dynamic model. J Power Sources 2007;173:310e24. [29] Springer TE, Zawodzinski TA, Gottesfeld S. Polymer electrolyte fuel cell model. J Electrochem Soc 1991;138:2334e42. [30] Noiying P, Hinaje M, Thounthong P, Rae¨l S, Davat B. Using electrical analogy to describe mass and charge transport in PEM fuel cell membrane transport models. Renew Energy 2012;44:128e40. [31] Hossain M, Islam SZ, Pollard P. Investigation of species transport in a gas diffusion layer of a polymer electrolyte membrane fuel cell through two-phase modelling. Renew Energy 2013;51:404e18. [32] Misran E, Hassan NSM, Daud WRW, Majlan EH, Rosli MI. Water transport characteristics of a PEM fuel cell at various operating pressures and temperatures. Int J Hydrogen Energy 2013;38:9401e8. [33] Cao T-F, Lin H, Chen L, He Y-L, Tao W-Q. Numerical investigation of the coupled water and thermal management in PEM fuel cell. Appl Energy 2013;112:1115e25.
7088
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y 3 9 ( 2 0 1 4 ) 7 0 7 7 e7 0 8 8
[34] Morris DR, Sun X. Water-sorption and transport properties of nafion 117. J App Polym Sci 1993;50:1445e52. [35] Hinatsu JT, Mizuhata M, Takenaka H. Water uptake of perfluorosulfonic acid membranes from liquid water and water vapor. J Electrochem Soc 1994;141:1493e8. [36] Jalani NH, Datta R. The effect of equivalent weight, temperature, cationic forms, sorbates, and nanoinorganic additives on the sorption behavior of Nafion. J Membr Sci 2005;264:167e75. [37] Pimenov AV. Water sorption of sulfo-cationic membranes and prediction of their ion-exchange selectivity; 1992. PhD Thesis. S-Petersburg. [38] Zawodzinski T, Derouin C, Radzinski S, Sherman R, Smith T, Springer T, et al. Water uptake by and transport through Nafion 117 membranes. J Electrochem Soc 1993;140:1041e7. [39] Jeck S, Scharfer P, Kind M. Absence of Schroeder’s paradox: experimental evidence for water-swollen Nafion membranes. J Membr Sci 2011;373:74e9. [40] Schroeder P. Ueber Erstarrungs- und Quellungserscheinungen von Gelatine. Z Phys Chem 1903;45:75e117. [41] Choi P, Datta R. Sorption in proton-exchange membranes. An explanation of Schroeder’s paradox. J Electrochem Soc 2003;150:E601e7. [42] Eikerling M, Berg P. Poroelecroelastic theory of water sorption and swelling in polymer electrolyte membranes. Soft Matter 2011;7:5976e90. [43] Bass M, Freger V. An experimental study of Schroeder’s paradox in Nafion and Dowex polymer electrolytes. Desalination 2006;199:277e9. [44] Onishi LM, Prausnitz JM, Newman J. Water-Nafion equilibria. Absence of Schroeder’s paradox. J Phys Chem B 2007;111:10166e73. [45] Freger V. Hydration of ionomers and Schroeder’s paradox in nafion. J Phys Chem B 2009;113:24e36. [46] Davankov VA, Pastukhov AV. Paradoxes of thermodynamics of swelling equilibria of polymers in liquids and vapors. J Phys Chem B 2011;115:15188e95. [47] Volino F, Pineri M, Dianoux AJ, De Geyer A. Water mobility in a water-soaked Nafion membrane. A high-resolution neutron quasielastic study. J Polym Sci Part A-2 Polym Phys 1982;20:481e96. [48] Zawodzinski T, Neeman M, Sillerud L, Gottesfeld S. Determination of water diffusion coefficients in perfluorosulfonate ionomeric membranes. J Phys Chem 1991;95:6040e4. [49] Rivin D, Kendrick CE, Gibson PW, Schneider NS. Solubility and transport behavior of water and alcohols in Nafion. Polymer 2001;42:623e35. [50] Jayakody JRP, Stallworth PE, Mananga ES, Farrington-Zapata J, Greenbaum SG. High pressure NMR study of water self-diffusion in NAFION-117 membrane. J Phys Chem B 2004;108:4260e2. [51] Suresh G, Scindia YM, Pandey AK, Goswami AA. Selfdiffusion coefficient of water in Nafion-117 membrane with different monovalent counterions: a radiotracer study. J Membr Sci 2005;250:39e45. [52] Majsztrik P, Bocarsly A, Benziger J. Water permeation through nafion membranes: the role of water activity. J Phys Chem B 2008;112:16280e9. [53] Luo Z, Chang Z, Zhang Y, Zh Liu, Li J. Electro-osmotic drag coefficient and proton conductivity in Nafion membrane for PEMFC. Int J Hydrogen Energy 2010;35:3120e4. [54] Ise M, Kreuer KD, Maier J. Electroosmotic drag in polymer electrolyte membranes: an electrophoretic NMR study. Solid State Ionics 1999;125:213e23.
[55] Broadbent SR, Hammersley JM. Percolation processes. I. Crystals and mazes. Math Proc Camb Phil Soc 1957;53(03):629e41. [56] Kirkpatrick S. Percolation and conduction. Rev Mod Phys 1973;45:574e8. [57] Efros AL. Fizika i geometriya besporyadka [Physics and geometry of disorder]. Moscow: Nauka; 1982. [58] Roldughin VI, Vysotskii VV. Percolation properties of metalfilled polymer films, structure and mechanisms of conductivity. Progress in Organic Coatings 2000;39:81e100. [59] Shklovskii BI, Efros AL. Theory of occurrence and conductivity of highly heterogeneous media. Usp Fiz Nauk 1975;117:401e35. [60] Sarychev AK, Vinogradoff AP, Karimov AM. Calculation of percolation conductivity in 3D. J Phys C Solid State Phys 1985;18:L105e8. [61] Karpenko LV. Electro-transport properties of the ionexchange membranes as a function of the polymer composition and the concentration of the equilibrium electrolyte solution. PhD Thesis [in Russian]. Krasnodar, Russia: Kuban State University; 1999. [62] Stumper J, Campbell SA, Wilkinson DP, Johnson MC, Davis M. In-situ methods for the determination of current distributions in PEM fuel cells. Electrochim Acta 1998;43:3773e83. [63] Hsing I-M, Futerko P. Two-dimensional simulation of water transport in polymer electrolyte fuel cells. Chem Eng Sci 2000;55:4209e18. [64] Kreuer KD. On the development of proton conducting materials for technological applications. Solid State Ionics 1997;97:1e15. [65] Meier F, Eigenberger G. Transport parameters for the modelling of water transport in ionomer membranes for PEM-fuel cells. Electrochim Acta 2004;49:1731e42. [66] Ren X, Gottesfeld S. Electro-osmotic drag of water in poly (perfluorosulfonic acid) membranes. J Electrochem Soc 2001;148:A87e93. [67] Saito M, Naoko A, Kikuko H, Okada T. Mechanisms of ion and water transport in perfluorosulfonated ionomer membranes for fuel cells. J Phys Chem B 2004;108:16064e70. [68] Edition 11/2010. AVL List GmbH Manual AVL FIRE version 2010; 2010. p. 91. [69] Wasterlain S, Candusso D, Hissel D, Harel F, Bergman P, Menard P, et al. Study of temperature, air dew point temperature and reactant flow effects on proton exchange membrane fuel cell performances using electrochemical spectroscopy and voltammetry techniques. J Power Sources 2010;195:984e93. [70] Al-Baghdadi MARS, Shanad Al-Janabi HAK. Effect of operating parameters on the hygroethermal stresses in proton exchange membranes of fuel cells. Int J Hydrogen Energy 2007;32:4510e22. [71] Jang J-H, Chiu H-C, Yan W-M, Sun W-L. Effects of operating conditions on the performances of individual cell and stack of PEM fuel cell. J Power Sources 2008;180:476e83. [72] Fink C, Fouquet N. Three-dimensional simulation of polymer electrolyte membrane fuel cells with experimental validation. Electrochim Acta 2011;56:10820e31. [73] Baumgartner WR. Characterization of polymer-electrolytemembrane fuel cells. PhD Thesis. Graz University of Technology; 2007. [74] Baumgartner WR, Parz P, Fraser SD, Wallno¨fer E, Hacker V. Polarization study of PEMFC with four reference electrodes at hydrogen starvation conditions. J Power Sources 2008;182:413e21.
Available online at www.sciencedirect.com
ScienceDirect journal homepage: www.elsevier.com/locate/he
A novel membrane transport model for polymer electrolyte fuel cell simulations L. Karpenko-Jereb a,*, P. Innerwinkler a, A.-M. Kelterer a, C. Sternig a, C. Fink b, P. Prenninger b, R. Tatschl b a b
Institute of Physical and Theoretical Chemistry, Graz University of Technology, Graz, Austria AVL List GmbH, Graz, Austria
article info
abstract
Article history:
This work presents the development of a 1D model describing water and charge transport
Received 29 November 2013
through the polymer electrolyte membrane (PEM) in the fuel cell. The considered driving
Received in revised form
forces are electrical potential, concentration and pressure gradients. The membrane
12 February 2014
properties such as water diffusion and electro-osmotic coefficients, water sorption and
Accepted 17 February 2014
ionic conductivity are treated as temperature dependent functions. The dependencies of
Available online 25 March 2014
diffusion and electro-osmotic coefficients on the membrane water concentration are described by linear functions. The membrane conductivity is computed in the framework
Keywords:
of the percolation theory under consideration that the conducting phase in the PEM is
Polymer electrolyte membrane fuel
formed by a hydrated functional groups and absorbed water. This developed membrane
cell
model was implemented in the CFD code AVL FIRE using 1D/3D coupling. The simulated
Membrane transport model
polarization curves at various humidification of the cathode are found in good agreement
CFD code AVL FIRE
with the experiments thus confirming the validity of the model.
Percolation
Copyright ª 2014, Hydrogen Energy Publications, LLC. Published by Elsevier Ltd. All rights reserved.
Water sorption isotherm Polarization curve
Introduction Currently, hydrogen Polymer Electrolyte Membrane Fuel Cells (PEMFCs) are one of the most promising devices as high efficient sources of electric power for electric vehicles. The main advantages of the fuel cells are their high energy conversion efficiency and low noise, as well as their high power and energy density even under consideration of the H2-tank [1e4]. The main goal of the PEMFC modeling is to calculate the cell potential and its efficiency. For this purpose it is necessary to find out which processes contribute to the potential losses and the extent of their contribution. If Ohmic voltage drops across
GDLs and bipolar plates are neglected, the PEMFC potential is expressed by the equation: Ucell ¼ Uoc hc ha hmem
(1)
Ucell e cell potential; Uoc e open-circuit voltage; hc, ha e potential losses at cathode and anode; hmem e membrane overpotential. The main purpose of a membrane transport model is to calculate the membrane over-potential. Fig. 1 depicts the most important processes in the proton exchange membrane coated with a catalyst layer (PEMCC) taking place in a fuel cell at operation. These are the proton and electro-osmotic transports, the water generation resulted by the electro-
* Corresponding author. Tel.: þ43 316 873 32241; fax: þ43 316 873 32202. E-mail address: [email protected] (L. Karpenko-Jereb). http://dx.doi.org/10.1016/j.ijhydene.2014.02.083 0360-3199/Copyright ª 2014, Hydrogen Energy Publications, LLC. Published by Elsevier Ltd. All rights reserved.
7078
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y 3 9 ( 2 0 1 4 ) 7 0 7 7 e7 0 8 8
chemical reaction, the water diffusion and convection. Water diffusion arises due to a water concentration difference between the cathode side and the anode side. Usually the cathode has a higher water concentration caused by water formation as a result of the chemical reaction between Hþ and O2. The electro-osmosis (often called electro-osmotic flux or water drag) is caused by the conjugated transport of the water molecules with the protons in the electric field: the water molecules mainly move in the protons’ solvation shells. The convective flux of water is caused by a pressure gradient. Theoretical modeling of the transport in an ion-exchange membrane (IEM) is based usually on structural and transport models [2,5,6]. The structural models give a simplified concept of the inner structure of the membranes. Currently, the following structural models of the IEMs are known: clusternetwork [7,8], dusty-fluid [9,10], capillary [11,12], microheterogeneous two-phase [13,14] and percolation [15e17]. The transport models describe the transport of ions, water and gases through the membrane. Table 1 displays a comparison of several membrane transport models published in the scientific papers [18e28]. Based on the governing equations used for calculating membrane overpotential, the membrane transport models can be separated into three main groups: first e StefaneMaxwell equation based models [18,19]; second e NernstePlanck equation based models [20,21]; third - Ohm’s law based models [6,22e28]. Most of the models are known in the third group. From experimental investigations it is known that the temperature and pressure gradients have a weak influence on fuel cell voltage. Therefore most of the membrane transport models only consider two driving forces: the gradients of water concentration and those of electrical potential. Models such as those from of Rowe et al. [21], Senn and Poulikakos [24] additionally take into account the pressure gradient. Table 1 also displays the approaches to the description of the membrane transport properties (water diffusion and electro-osmotic coefficients) and the proton concentration as a function of the membrane water concentration as well as approximations in the calculations of the total water flux through the membrane. Most of the models take into account the dependencies of the water diffusion and electro-osmotic coefficients on the membrane water
concentration. The simplified model of Berg et al. [20] sets the electro-osmotic coefficient (Cdrag) in the perfluorinated membrane equal to unity and considers the dependence of the proton concentration on the membrane hydration level. Usually, the total water flux via the membrane is considered as a sum of diffusion and electro-osmotic flows. Kulikovsky proposed [6,22] a description of the transport phenomena in PEM based on the approximation that the water back diffusion is completely compensated by the electro-osmotic transport. This approach allows calculating the total water flux through the membrane analytically. But this model is not able to predict electrode water flooding effects. The membrane transport models differ also in the formulation of boundary conditions and in the definition of the water concentration at interfaces of PEM/CL, CL/GDL as well as in the calculation of the water concentration profile in the fuel cell. In Berg’s model [20] the water sorption isotherm is used for the calculations of the c water concentrations at the interfaces GDL/CL (Fig. 1), Ca w ,Cw , while in other models [6,22,24,28] it is used for the computation of the water concentrations at the membrane boundaries with CLs, Caw and Ccw . Whereas the influence of temperature on water sorption properties of the membrane is described insufficiently in the models [6,18e33], however the experimental investigations [34e36] showed the temperature effect on the membrane water sorption isotherm, which is especially pronounced for the equilibrium of the membrane with the saturated water vapor and liquid water. The aim of this present work was the development of a membrane transport model for PEMFC simulations, which describes water and charge transports through the polymer electrolyte membrane, and considers the influence of temperature and water concentration on the important membrane properties such as the water sorption isotherm, the water diffusion and electro-osmotic coefficients as well as on ionic conductivity.
Membrane transport model Approximations The developed membrane transport model is based on the approximations summarized below:
Fig. 1 e Transport and chemical phenomena in proton exchange membrane coated by catalyst in an operating PEMFC. These phenomena are considered in the developed membrane transport model.
1. The transport of water and protons is considered only in one dimension. 2. It is a steady-state model. 3. Three driving forces are considered to describe the transport: gradients of electrical potential, concentration and pressure. 4. The model is an isothermal model, i.e. it does not take into account an effect of temperature gradient across the membrane on the water transport. For the calculations in the membrane transport model an average value between the temperatures at the membrane interfaces at the cathode and anode sides is taken. 5. Four types of water motions are considered in the model, as shown in Fig. 1: water formation at the cathode, water diffusion, electro-osmosis and water convection.
7079
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y 3 9 ( 2 0 1 4 ) 7 0 7 7 e7 0 8 8
Table 1 e Comparison of the membrane transport models. Model
Driving forces
1. StefaneMaxwell equation Carnes, Djilali, [18] Baschuk, Li [19] 2. NernstePlanck equation Berg et al. [20] Rowe, Li [21] 3. Ohm’s law Kulikovsky A. [6,22] Bao et al. [23] Senn, Poulikakos [24] Gurau, Mann [25] Wu et al. [26] Sui et al. [27] del Real et al. [28]
Approximations
D4
DCw
Dp
Dw
Cþ
Cdrag
jm w
þ þ
þ þ
f(Cw) f(Cw)
Const Const
f(Cw) f(Cw)
s0 s0
þ þ
þ þ
þ
f(Cw) f(Cw)
f(Cw) Const
¼1 f(Cw)
s0 s0
þ þ þ þ þ þ þ
þ þ þ þ þ þ þ
þ
Const Const Const f(Cw) f(Cw) f(Cw) f(Cw)
Const Const Const Const Const Const Const
f(Cw) f(Cw) f(Cw) f(Cw) f(Cw) f(Cw) f(Cw)
¼0 s0 s0 s0 s0 s0 s0
D4, DCw, Dp e gradients of the electrical potential, water concentration and pressure correspondingly; Dw e water diffusion coefficient; Cþ e proton concentration; Cdrag e electro-osmotic coefficient; jm w e total water flux throughout the membrane; Cw e water concentration in the membrane.
6. The boundary conditions of water transport utilize the modified equations from Berg’s et al. model, 2005 [20]. 7. The membrane over-potential is calculated using Ohm’s Law. 8. Water diffusion and electro-osmotic coefficient, as well as membrane conductivity, are functions of the membrane water concentration Cw and temperature. 9. Water sorption isotherm depends on the temperature. 10. The activation energies of water diffusion, electroosmotic coefficient and membrane conductivity do not depend on membrane water concentration.
At 4 1, if liquid water is formed at both boundaries CL/ PEM, the water concentrations at anode and cathode sites are equal Caw ¼ Ccw . In this case the diffusive flux is zero, but the convective flux is taken into account and consequently the total water flux through the membrane is described by: i pl;a pl;c jm w ¼ Cdrag $ þ Khyd $ F Lmem
The membrane over-potential is calculated based on Ohm’s law and expressed by: hmem ¼
The development of the auxiliary equations was based on analysis of experimental data for perfluorinated membranes of Nafion type.
Governing equations All symbols are explained in the “Nomenclature” at the end of the paper. The water diffusion is expressed by Fick’s Law: Dw
dCw jdif dCw ¼ w ; jdif w ¼ a$Dw $ dz a dz
(2)
(4)
At 4 < 1, water vapor operating condition, the convective flux is negligible and the total water flux through the membrane is found from the sum of diffusive and electro-osmotic fluxes: jm w ¼ a$Dw
dCw i þ Cdrag $ F dz
(7)
Boundary conditions The balance of water fluxes between the polymer electrolyte membrane and catalyst layer is based on the boundary conditions proposed by Berg et al. [20]: at anode catalyst layer a a (8) jm w ¼ a$gH2 O;a Cw Cw at cathode catalyst layer jm w þ
i ¼ a$gH2 O;c Ccw Cc w 2F
(3)
The convective flux is calculated from: pl;a pl;c jconv ¼ Khyd $ w Lmem
i $Lmem smem
(9)
At 4 < 1 the balance of the water fluxes in PEM and catalyst layers are expressed as:
The electro-osmotic flux is expressed by: i jdrag ¼ Cdrag $ w F
(6)
(5)
a$Dw
dCw i þ Cdrag ¼ a$gH2 O;a Caw Ca w F dz
(10)
a$Dw
dCw i i ¼ a$gH2 O;c Ccw Cc þ Cdrag þ w F 2F dz
(11)
g e the mass transfer coefficient of water from the GDL to the catalyst layer is calculated via the equation: gH2 O ¼
DGDL w dn þ 0:5$LCL
(12)
The water activity inside the GDL close to the interface of the membrane Cw , is calculated from the water sorption
7080
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y 3 9 ( 2 0 1 4 ) 7 0 7 7 e7 0 8 8
isotherm, which defines the dependence of membrane water concentration on equilibrium water vapor.
The development of the auxiliary equations The development of the auxiliary equations in order to describe the dependences of membrane properties on the relative humidity or membrane water concentration and temperature was based on the thoroughly analysis of experimental data available in the scientific papers [16,18e29,34,55,62e67]. Most of the experimental investigations on the Cdrag, Dw, smem as a function of Cw, reported in the scientific literature, were performed at 25 C. In order to describe the transport properties of the membrane at other temperatures Arrhenius’ equation was applied with an assumption that the activation energies do not depend on the membrane water concentration.
Water sorption isotherm In order to define the water concentration at the membrane boundaries, knowledge of the equilibrium between the water concentrations in the membrane and catalyst/gas diffusion layers is required. This equilibrium can be described by a water sorption isotherm e the relation between the membrane hydration number (number of water molecules per one functional group) and the equilibrium humidity. The water sorption isotherm by the polymer electrolyte membranes is intensively studied using different experimental methods [34e39]. Special attention is focused on the difference between the equilibriums of the membrane with the saturated water vapor and with liquid water. For the first time Schroeder [40] indicated that gelatin, an amino acid related to natural polyelectrolytes, absorbs significantly less water from saturated water vapor than at contact with liquid water. At present, this phenomenon is well-known as Schroeder’s paradox and currently being intensively discussed in the scientific literature [39,41e46]. The existence of this phenomenon has theoretically been proven [41,42]. However, it is an interesting fact, that several experimental investigations [39,43,44] did not confirm the presence of Schroeder’s paradox. Utilizing the thermodynamic model approaches, Lu [41] as well as Eikerling and Berg [42] showed that the transition of water from gas to liquid phase attains a discontinuous jump in the total water concentration in polymer electrolyte membranes. In the numerical modeling of PEMFCs, this discontinuity of the membrane property can cause difficulties. Therefore, most of the published models do not take into account Schroeder’s paradox and describe the water sorption isotherms by continuous polynomial functions [18e29]. In the present work a mathematical function for the membrane water sorption isotherm was developed based on the analysis of the published experimental investigations [34,35,37,38] and the following assumptions: 1) The maximal membrane water concentration grows with increasing temperature. 2) The appearance of liquid water in the system causes a sharp but continuous rise of the membrane water concentration. The second statement was guided by the following considerations: at convenient operating conditions of the PEMFC,
Fig. 2 e Temperature dependence of the specific water concentration in Nafion membrane equilibrated with liquid water.
water is present in both vapor and condensed (liquid) states. Consequently, both phases contact the membrane simultaneously. Thus, the water at the membrane boundary most likely represents a mixture of vapor with droplets of the liquid water. Based on the theory proposed by Eikerling and Berg [42], we can assume that water droplets will be more easily absorbed by the membrane compared to water vapor. But a few water droplets will not be enough in order to increase the membrane water concentration to the maximal value. Thus, at the occurrence of the liquid water in the real PEMFC most possibly the growth of the water content will transit continuously in the membrane. Fig. 2 displays the dependencies of the membrane water concentration for Nafion membranes, equilibrated with liquid water at various temperatures. These data were reported by Morris [34] and Hinatsu [35]. As can be seen in the figure, the membrane water concentration grows linearly with increasing temperature and the two data sets are found in complete agreement and are well described by the linear function: ¼ 0:138$T 28:31 cmax w
(13)
Fig. 3 presents the experimental data on the water sorption isotherms of Nafion membranes measured in a wide range of water vapor activity [34,35,37,38]. We should note that different methods were applied for the measurements, and different devices for controlling the establishment of the equilibrium between the membrane and water vapor were used. Therefore, the data determined, for example, at 25 C by different research groups do not demonstrate complete agreement. We have tried to find a correlation between the water sorption isotherm shape and temperature from these data. On the one hand, Morris [34] showed that the membrane water sorption is higher at 50 C than at 25 C, but on other hand the data of Hinatsu [35] measured at 80 C lie below the values obtained by Pimenov at 25 C [37] and by Morris at 50 C [34]. Perhaps the functions Cw ¼ f(4) have different trends at different temperatures but currently this trend cannot be
7081
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y 3 9 ( 2 0 1 4 ) 7 0 7 7 e7 0 8 8
accurately formulated from the available experimental data. Therefore, we assume that at 0 < 4 < 0.97 the shape of the water sorption isotherm practically does not depend on the temperature (Fig. 3a), but at 4 0.97 the growth of membrane water concentration increases with the temperature (Fig. 3b). The suggested fitting describing the water sorption isotherm is expressed by the equation: Cw ¼ 1:55 þ 13:71$f 24:37$f2 þ 21:87$f3 f cmax w ;f
(14)
f ðcmax w ; 4Þ
where is a multiplier, which is directly proportional described to the maximal membrane water concentration cmax w by Eq. (13) and the relative humidity 4. As seen from Fig. 3b the maximal water concentrations (at 4 ¼ 1) computed by Eq. (14) are found in quite good agreement with experimental data.
same order of magnitude and lie in a range from 0.39 1010 to 5.80 1010 m2/s. The water diffusion coefficients obtained by Rivin and Zawodzinski grow directly with the rising hydration number of the membrane. The water diffusion coefficient reported by Jayakody agrees with data measured by Zawodzinski, and the value of Volino is compared with the data published by Rivin. In the presented model a linear fit (Eq. (15)) of the experimental relationship of Dw vs. the membrane hydration number is suggested, which starts at zero and represents the average of the presented experimental points. Dw ¼ dTw1 $Cw
(15)
where dTw1 is the slope of the straight line fitting, where dTw1 ¼ 2:65 1011 m2 =s at T ¼ 298.15 K.
Electro-osmotic coefficient vs. membrane water concentration Water diffusion coefficient vs. membrane water concentration Fig. 4 displays the dependencies of the water diffusion coefficient on the hydration number of Nafion membranes measured at 25 C by Volino 1982 [47], Zawodzinski 1991 [48], Rivin 2001 [49], Jayakody 2004 [50], Suresh 2005 [51] and Majsztrik 2008 [52]. As seen in the figure, the values have the
Fig. 5 presents the electro-osmotic coefficients in the perfluorinated membranes as a function of membrane water concentration measured by Luo at 25 C [53], Ise at 27 C [54] and Springer at 30 C [29]. The linear fit was performed on the data reported by Luo [53], because the data evaluation for the development of the model was carried out at 25 C. No electroosmosis is possible in the dry membrane. Therefore the suggested fit passes through the origin of the coordinates. Thus the electro-osmotic coefficient is expressed by the equation: 1 $Cw Cdrag ¼ cTdrag
(16)
1 ¼ 0:11 at T ¼ 298.15 K. for Nafion membrane cTdrag
Dw ⋅ 1010 m2 s 6
4
2
0 0
Fig. 3 e The specific water concentration of Nafion membrane as a function of the equilibrium relative water vapor humidity: a) in the wide range of 4 [ 0L1.0; b) for high 4 [ 0.96L1. The points are experimental data from Refs. [34,35,37,38], lines are the suggested fittings.
4
8
12
16
20
Rivin, 2001
suggested fitting
Suresh, 2005
Volino, 1982
Jayakody, 2004
Majsztrik, 2008
Cw
Zawodzinski, 1991
Fig. 4 e The water diffusion coefficient of Nafion membrane as a function of the specific water concentration, T [ 25 C. The points are experimental data from Refs. [47e52], the straight line is the suggested fitting.
7082
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y 3 9 ( 2 0 1 4 ) 7 0 7 7 e7 0 8 8
Ionic conductivity vs. membrane water concentration As shown in the investigations [15e17] the dependence of membrane conductivity on the water concentration can successfully be described by the equation of the percolation theory. The basic ideas of this theory were formulated by Broadbent and Hammersly in 1957 [55]. According to the percolation theory, after reaching some critical volume fraction of the conducting phase fcr, the conductivity of a composition consisting of dielectric and conductor dramatically begins to grow with the increasing fraction of this phase. This change from dielectric to conductor has been called the percolation transition and the critical content of the conducting phase fcr is a threshold of the percolation. At f > fcr the dependence of the conductivity on the conducting phase fraction is described by the equation: t smem ¼ sT1 $ f fcr
(17)
where s - pre-factor of percolation equation, t e percolation exponent, f e volume fraction of the conducting phase; fcr e critical volume fraction of the conducting phase. The quantities s, fcr and t are the parameters of the percolation equation. Calculations based on the theory of probability as well as experimental investigations have shown that the percolation parameters lie in the intervals f ¼ 0.15 0.03, t ¼ 1.6 0.4 [56e60]. The conducting phase in perfluorinated membranes represents the connected hydrophilic cluster domains consisting of hydrated functional groups and networked absorbed water. The following equation was suggested [16,61] in order to calculate the volume fraction of the conducting phase in ionexchange membranes: f ¼ a$MSO3 þ W
(18)
where a e ion-exchange capacity, i.e. the concentration of acid groups; MSO3 - molecular weight of the acid group and W e water mass fraction in PEM.
Cdrag 2.5 2
Fig. 6 e a) The dependence of the specific conductivity of Nafion membrane on the volume fraction of the conducting phase; b) bi-logarithmic dependence lg(smem) L lg(f L fcr), used for the determination of the exponent s of the percolation theory equation. The points are experimental data [53]; the line on (a) is the suggested fitting.
The derivation of the equation is described in detail in the same references [16,61]. Fig. 6a demonstrates the dependences of Nafion conductivity measured by Luo et al. [53] as a function of the volume fraction of the conducting phase calculated using Eq. (18). From these experimental data the estimated critical volume fraction of the conducting phase is fcr ¼ 0.15 which corresponds well to the theoretical range. The treatment of these data in bi-logarithmic coordinates lg(smem) lg(f fcr) (Fig. 6b) has allowed the determination of the parameters s and t. Taken into account that fcr ¼ 0.15 the obtained values of these parameters are: t ¼ 1.2 and sT1 ¼ 1:16½S=cm at T ¼ 298.15 K (Fig. 6).
Temperature dependencies of the membrane transport characteristics
1.5 1
Since the temperature dependencies of the water diffusion, electro-osmotic coefficients and proton conductivity of the perfluorinated membranes are well described by Arrhenius’ equation as shown in the investigations [29,62e65], the expressions (15e17)" can be transformed in (19e21) # dif E 1 1 (19) DTw2 ¼ dTw1 $Cw $exp a $ T2 T1 R
Luo, 2010 Ise, 1999 Springer, 1991 suggested fitting
0.5 0 0
5
10
15
20
Cw
Fig. 5 e The electro-osmotic coefficient of Nafion membrane as a function of the membrane water concentration. The points are the experimental data from Refs. [29,53,54], the straight line is the suggested fitting.
" 2 1 CTdrag ¼ cTdrag $Cw $exp
# drag Ea 1 1 $ T2 T1 R
t Es 1 1 2 ¼ sT1 $ f fcr $exp a $ sTmem R T2 T1
(20)
(21)
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y 3 9 ( 2 0 1 4 ) 7 0 7 7 e7 0 8 8
7083
Fig. 7 e Temperature dependence of the electro-osmotic coefficient (a) and specific conductivity (b) of Nafion membranes in coordinates: a) ln(Cdrag) L (L1000/T); b) ln(smem) L (L1000/T). The points are experimental data from Refs. [54,66,67], the straight lines are linear trends fitting the experimental data.
Table 2 presents the values of the activation energies of Cdrag and smem, which were determined in this work from the temperature dependencies of the experimental data [54,66,67] using a linear fit of ln(Cdrag,smem) versus (1/T) as shown in Fig. 7. The table also contains the activation energy values of Dw and smem published in Refs. [29,64,65]. As seen from the table the activation energy of the water diffusion coefficient of Nafion membrane reported by Kreuer [64] decreases with the growth of membrane hydration level. Unfortunately, any consideration of the mathematical correlation of this dependence would significantly complicate the numerical solution for the total water flux. For this reason it has been assumed that the activation energies do not depend on water concentration and therefore the average values were taken for the simulations.
calculation of the membrane over-potential is depicted in Fig. 8. In order to determine the total water flux in the membrane and membrane water concentration the following padrag , Esa . Water rameters are needed: dw, cdrag, Edif a , Ea concentrations at the boundaries of GDL/CL are calculated from the water sorption isotherm. The total water flux jm w is
Input T2 , i T1 d wT1 , cdrag , Eadif , Eadrag , Eaσ
C wa* , C wc* Results and discussion Programming of the developed membrane transport model was carried out in Fortran95. The general sequence of the Table 2 e The values of the activation energy of water diffusion coefficient, electro-osmotic coefficient and ionic conductivity of Nafion membrane. Reference Water diffusion coefficient Kreuer K., 1997 [64]
Electro-osmotic coefficient Ren X., Gottesfeld S., 2001 [66] Luo Z. et al., 2010 [54] Ionic membrane conductivity Meier F., Eisenberg G., 2004 [65] Springer T. et al., 1991 [29] Luo Z. et al., 2010 [54] Saito M. et al., 2005 [67] a
Cw
Ea 103, [J/mol]
3 5 10
22.960 19.776 16.690
22 22
6.818a 8.019a
e e 22 21
9.889 10.537 9.980a 8.445a
The values were calculated in the present work from the experimental data (Fig. 7) published in Refs. [53,54,66].
j wm
C wa , C wc , C w ( z ) C wmean
T σ mem 2
η mem Output C w ( z ), j wm ,η mem Fig. 8 e Flow chart of the calculation steps in the developed membrane transport model.
7084
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y 3 9 ( 2 0 1 4 ) 7 0 7 7 e7 0 8 8
Fig. 9 e Computational meshes of a) single channel cell, b) Full coupled cell. Images of the experimental cell: c) cathode bipolar plate and end plate, d) whole assembly.
found iteratively by Newton’s method, and the calculation of the water concentration profile Cw(z) is performed by RungeeKutta’s method. From the average membrane water conthe membrane conductivity smem is centration Cmean w calculated using the equation of the percolation theory and then the membrane over-potential hmem is determined based on Ohm’s Law. The developed PEM Model was implemented in the CFD Code AVL FIRE [68]. The convergence of AVL FIRE with this new PEM Model was carefully tested using a Single Channel Fuel Cell (Fig. 9(a)) in around 300 calculations of FC performance at various operating parameters: p, T, 4, l. The voltageecurrent characteristics simulated for varied p, T and l are presented in Figures S1eS4 (Appendix). The model calculations indicate a visible growth of the fuel cell current density with increasing temperature and outlet pressure. Raising the cell temperatures enhances reaction kinetics and contributes to obtain lower activation losses. Higher pressures mean higher hydrogen and oxygen concentrations, which result in enhanced reaction rates and thus increasing current densities. The change of the anode gas stoichiometric ratio does not lead to a visible effect on the cell performance. The increasing cathode stoichiometric flow ratio results in a light enhancement of the current densities, it is linked to increasing reaction kinetics at the higher amount of oxygen delivered to the cathode catalyst layer. The computational results are found in qualitative agreement with experimental observations [69e71]. Advanced validation of the developed PEM Model was carried out using a Full Coupled Fuel Cell (Fig. 9(b)) at various relative humidities at the cathode inlet. The computational geometry consists of four parts: the first part is the main part of
the co-simulation, containing fuel cell and gas channels with 677,264 computational cells (CCs); the second part with bipolar plates is represented by 994,532 CCs; the third part is the end plates with 36,154 CCs; fourth e the cooling channels divided into 23,408 CCs. In the first part of the simulation cell the mass transport is calculated, in the other three parts mainly the heat transport between the cell elements is computed. More detailed descriptions of the simulation procedure of PEMFC by CFD Code AVL FIRE are found in Refs. [68,72]. The experimental measurements were carried out in the Christian-Doppler (CD) Laboratory for Fuel Cell Systems using a fully automated test rig [73,74]. The test rig is capable of supplying well defined mixtures of hydrogen, nitrogen, oxygen, and air by using mass flow controllers. Bubble humidifiers are utilized for gas humidification, and the humidity levels are checked with downstream humidity sensors. Electrochemical potentials are measured utilizing a USB multiflexer card from National Instruments (NI USB-6218). The fuel cell is composed of carbon bipolar plates with 13 straight parallel flow channels (Fig. 9(c)). The cathode bipolar plate is segmented into 10 individual elements, which are electronically insulated from each other (Fig. 9(d)). The bipolar plates are embedded between two end plates made of stainless steel. Cooling channels are installed in the end plates, which provide a precise temperature control of the fuel cell by applying water as the cooling fluid. The active area of the single fuel cell is 25 cm2. The fuel cell assembly shown in Fig. 9(b) completely corresponds to the fuel cell construction (Fig. 9(d)) used for the experimental measurements. The material characteristics and operating conditions applied in the simulations are presented in Tables 3 and 4 correspondingly.
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y 3 9 ( 2 0 1 4 ) 7 0 7 7 e7 0 8 8
7085
Table 3 e Values of materials parameters applied in the simulations. Parameter GDL Thickness Thermal conductivity Permeability Inplane Through plane Porosity Tortuosity Average contact angle Catalyst layers Thickness Cathode Anode Exchange current density Cathode Anode Transfer coefficient Cathode Anode Membrane Membrane thickness Sulfonic acid group concentration Hydraulic permeability of membrane Activation energy of water diffusion Activation energy of electroosmosis Activation energy of membrane conductivity Constants Faraday constant Universal Gas Constant
Symbol
Value
Unit
lGDL ss
3.00$105 [m] 4.5 [W/mK]
KII Kt ˛ s q
2.33$1012 1.07$1014 0.78 1.5 122
LCL,c LCL,a
1.125$105 [m] 1.125$105 [m]
i0c i0a
3.2.105 1.0.1011
[A/m3] [A/m3]
ke,c ke,a
0.1 0.5
[] []
Lmem a Khyd
35$106 1900.0 1$1011
Edif a drag
[m2] [m2] [] [] [ ]
m mol/m3 mol/ (m Pa s) 19.809$103 J/mol
Ea
7.418$103 J/mol
Esa
9.713$103 J/mol
F R
96485.309 C/mol 8.31451 J/(K mol)
Fig. 10 displays the polarization curves of the fuel cell for varying relative humidity at the cathode inlet, calculated by the CFD code AVL FIRE and measured experimentally. A lower inlet relative humidity leads to a decrease of the membrane water concentration and, hence, to a lower proton conductivity and current density. This performance decrease occurs in both the experiment and simulation. However, at a lower humidity the current density is overestimated by the simulation. For two cell voltage values, which are marked by blue dotted lines on Fig. 10, the distributions of the current density
Table 4 e Operating conditions of PEMFC applied in the simulations. Parameter Inlet temperature, [K] Outlet gas pressure, [bar] At cathode At anode Stoichiometric flow ratio At cathode At anode Relative humidity At cathode At anode Voltage, [V]
Values 343.15 1.0 1.0 2.2 1.5 0.9/0.6/0.3 0.9 0.420e0.855
Fig. 10 e Polarization curves of PEMFC with perfluorinated membrane: simulation results (lines); experimental data (points). The operating conditions la/lc [ 2.2/1.5; T [ 70 C; P [ 101325 Pa. The experimental polarization curves were measured with sulfo-cationic perfluorinated membrane, gas diffusion layer SGL-35BC. The GDL and PEM properties are given in Table 3.
along the gas channel were calculated (Fig. 11). As seen from the figure the dependencies of current densities on the channel length at 4in ¼ 0.9 and 4in ¼ 0.6 are found in quite good agreement with the experimental data. In both the experiment and simulation a lower humidity moves the current density maximum towards the cathode outlet (y ¼ 0 cm). However, in the simulation this happens to a greater extent, especially for 4in ¼ 0.3.
Conclusion In the present work a new 1D membrane transport model was developed. In this model the calculation of the membrane over-potential is based on Ohm’s Law. The boundary conditions describing the equilibrium between the water fluxes through CL and PEM utilized the modified boundary conditions from Berg’s et al. [20] model. The dependencies of water diffusion and electro-osmotic coefficients are considered as a linear function of the membrane water concentration, their temperature functions are described using Arrhenius’ equation. The membrane conductivity is computed from the equation of the percolation theory. The feature of the applied approach that the volume fraction of the conducting phase in the membrane is considered as a sum of the volume fractions of the functional groups and absorbed water, is described in detail in Refs. [16,61]. For the first time, the developed model considers temperature influence on the maximal water concentration in the polymer electrolyte membrane. The jumping growth in the membrane water concentration at the transition from saturated water vapor to liquid water is considered as a continuous process. The developed model was implemented in CFD code AVL FIRE. The validation of the model was carried out using the
7086
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y 3 9 ( 2 0 1 4 ) 7 0 7 7 e7 0 8 8
Acknowledgment The authors would like to thank Prof. Hacker V. e Head of the Christian-Doppler (CD) Laboratory for Fuel Cell Systems at Institute of Chemical Engineering and Environmental Technology, Graz University of Technology, for the experimental data presented here. This work has been financially supported by the Austrian Research Promotion Agency (FFG) and AVL List GmbH in the framework of Bridge1 Program, Project No. 827559 “PEMMODEL” Investigation of PEM (proton exchange membrane) transport model for the 3D CFD AVL-FIRE Simulations) and IV2Splus Program, Project No.835811 “A3 FALCON” (Advanced 3D Fuel Cell Analysis and Condition diagnostics).
Appendix A. Supplementary data Supplementary data related to this article can be found at http://dx.doi.org/10.1016/j.ijhydene.2014.02.083.
Nomenclature
Fig. 11 e The distribution of the current density along the channel from inlet (y [ 0 cm) to outlet (y [ 12 cm). Points are experimental data, lines are simulation results.
experimental data obtained on a segmented single fuel cell. The calculated polarization characteristics of the fuel cell at the different humidifications of the inlet air at the cathode are found in quite good agreement with the experimentally measured curves.
Outlook The presented model has been extended from 1D to 3D. The CFD code AVL FIRE with the implemented model will be commercially available in future releases.
Latin symbols a sulfonic acid group concentration, mol/m3 Cw normalized water concentration in the membrane normalized water concentration in the membrane, cmax w equilibrated with liquid water electro-osmotic coefficient (drag coefficient) in the Cdrag membrane 1 electro-osmotic coefficient in the membrane at cTdrag Cw ¼ 1 and T1 Dw water diffusion coefficient in the membrane, m2/s T1 water diffusion coefficient in the membrane at dw Cw ¼ 1 and T1, m2/s activation energy, J/mol Ea F Faraday constant, (A s)/mol f volume fraction of the conducting phase in the membrane critical volume fraction of the conducting phase fcr i electric current density, A/m2 jw water flux, mol/m2s m total water flux through the membrane, mol/m2s jw hydraulic permeability of the membrane, mol/ Khyd (m Pa s) membrane thickness, m Lmem p pressure, Pa R universal gas constant, J/(mol K) T temperature, K t exponent of the equation of the percolation theory, e U Potential, V y channel length, cm z normal direction to the membrane, m Greek symbols 4 water activity in gas phase, relative humidity water transfer coefficient, m/s gH2 O h potential loss, V
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y 3 9 ( 2 0 1 4 ) 7 0 7 7 e7 0 8 8
l s smem
stoichiometric coefficient pre-factor of the equation of the percolation theory, S/m specific membrane conductivity, S/m
Subscripts and superscripts a anode c cathode m; mem membrane mean related to an average value l related to liquid water CL related to catalyst layer GDL related to gas diffusion layer w; H2O water drag related to electro-osmotic transport of water dif diffusion conv related to convective flux s related to membrane conductivity Abbreviations CFD computational fluid dynamics CL catalyst layer FC fuel cell IEM ion exchange membrane GDL gas diffusion layer PEM polymer electrolyte membrane PEMCC polymer electrolyte membrane coated with catalyst layer PEMFC polymer electrolyte membrane fuel cell RH relative humidity
references
[1] Gou B, Na WK, Diong B. Fuel cells. Modeling, control and applications. CRC Press. Taylor & Francis Group; 2010. [2] Barbir F. PEM fuel cells. Theory and practice. Elsevier Academic Press; 2005. [3] Buechi FN, Inaba M, Schmidt ThJ. Polymer electrolyte fuel cell durability. New York: Springer; 2009. [4] Larminie J, Dicks A. Fuel cell systems explained http://www. amazon.de/Fuel-Systems-Explained-James-Larminie/dp/ 0471490261/ref¼sr_1_1? s¼books&ie¼UTF8&qid¼1374826216&sr¼11&keywords¼PEMþfuelþcell - #Wiley VCH; 2000. [5] Weber AZ, Newman J. Modeling transport in polymerelectrolyte fuel cells. Chem Rev 2004;104:4679e726. [6] Kulikovsky AA. Analytical modelling of fuel cells. Elsevier; 2010. [7] Hsu WY, Gierke TD. Ion transport and clustering in nafion perfluorinated membranes. J Membr Sci 1983;13:307e26. [8] Datye VK, Taylor PL, Hopfinger AJ. Simple model for clustering and ionic transport in ionomer membranes. Macromolecules 1984;17:1704e8. [9] Fimrite J, Struchtrup H, Djilali N. Transport phenomena in polymer electrolyte membranes I. Modeling framework. J Electrochem Soc 2005;152:A1804e14. [10] Thampan T, Malhotra S, Tang H, Datta R. Modeling of conductive transport in proton-exchange membranes for fuel cells. J Electrochem Soc 2000;147:3242e50. [11] Cwirko EH, Carbonell RG. A theoretical analysis of Donnan dialysis across charged porous membranes. J Membr Sci 1990;48:155e79.
7087
[12] Koter S. Transport of simple electrolyte solutions through ion-exchange membranes e the capillary model. J Membr Sci 2002;206:201e15. [13] Gnusin NP, Zabolotzki VI, Nikonenko VV, Meshechkov AI. Development of the principle of the generalize conductivity to describe transport phenomena in disperse systems caused by the different driving forces. Zh Fiz Khim 1980;54:1518e22. [14] Zabolotsky VI, Nikonenko VV. Effect of structural membrane inhomogeneity on transport properties. J Membr Sci 1993;79:181e98. [15] Hsu WY, Berkley JR, Meakin P. Ion percolation and insulatorto-conductor transition in nafion perfluorosulfonic acid membranes. Macromolecules 1980;13:198e200. [16] Berezina NP, Karpenko LV. Percolation effects in ionexchange materials. Colloid J 2000;62:676e84. [17] Sumberova V, Bleha M, Wodzki R. Percolation transition in heterogeneous ion exchange membranes. J App Polym Sci 1992;46:611e7. [18] Carnes B, Djilali N. Analysis of coupled proton and water transport in a PEM fuel cell using the binary friction membrane model. Electrochim Acta 2006;52:1038e52. [19] Baschuk JJ, Li X. Modeling of ion and water transport in the polymer electrolyte membrane of PEM fuel cells. Int J Hydrogen Energy 2010;35:5095e103. [20] Berg P, Promislow K, Pierre J, Stumper J, Wetton B. Water management in PEM fuel cells. J Electrochem Soc 2004;151:A341e53. [21] Rowe A, Li X. Mathematical modeling of proton exchange membrane fuel cells. J Power Sources 2001;102:82e96. [22] Kulikovsky AA. The effect of cathodic water on performance of a polymer electrolyte fuel cell. Electrochim Acta 2004;49:5187e96. [23] Bao C, Ouyang M, Yi B. Analysis of water management in proton exchange membrane fuel cells. Tsinghua Sci Tech 2006;11:54e64. [24] Senn SM, Poulikakos D. Polymer electrolyte fuel cells with porous materials as fluid distributors and comparisons with traditional channeled systems. Trans ASME 2004;126:410e8. [25] Gurau V, Mann Jr JA. A Continuum model for water transport in the ionomer-phase of catalyst coated membranes for PEMFCs. Hindawi Publ Corp Adv Mech Eng; 2010. Article ID 372795. [26] Wu H, Berg P, Li X. Non-isothermal transient modeling of water transport in PEM Fuel Cells. J Power Sources 2007;165:232e43. [27] Sui PC, Kumar S, Djilali N. Advanced computational tools for PEM fuel cell design part 1. Development and base case simulations. J Power Sources 2008;180:410e22. [28] del Real AJ, Arce A, Carlos B. Development and experimental validation of a PEM fuel cell dynamic model. J Power Sources 2007;173:310e24. [29] Springer TE, Zawodzinski TA, Gottesfeld S. Polymer electrolyte fuel cell model. J Electrochem Soc 1991;138:2334e42. [30] Noiying P, Hinaje M, Thounthong P, Rae¨l S, Davat B. Using electrical analogy to describe mass and charge transport in PEM fuel cell membrane transport models. Renew Energy 2012;44:128e40. [31] Hossain M, Islam SZ, Pollard P. Investigation of species transport in a gas diffusion layer of a polymer electrolyte membrane fuel cell through two-phase modelling. Renew Energy 2013;51:404e18. [32] Misran E, Hassan NSM, Daud WRW, Majlan EH, Rosli MI. Water transport characteristics of a PEM fuel cell at various operating pressures and temperatures. Int J Hydrogen Energy 2013;38:9401e8. [33] Cao T-F, Lin H, Chen L, He Y-L, Tao W-Q. Numerical investigation of the coupled water and thermal management in PEM fuel cell. Appl Energy 2013;112:1115e25.
7088
i n t e r n a t i o n a l j o u r n a l o f h y d r o g e n e n e r g y 3 9 ( 2 0 1 4 ) 7 0 7 7 e7 0 8 8
[34] Morris DR, Sun X. Water-sorption and transport properties of nafion 117. J App Polym Sci 1993;50:1445e52. [35] Hinatsu JT, Mizuhata M, Takenaka H. Water uptake of perfluorosulfonic acid membranes from liquid water and water vapor. J Electrochem Soc 1994;141:1493e8. [36] Jalani NH, Datta R. The effect of equivalent weight, temperature, cationic forms, sorbates, and nanoinorganic additives on the sorption behavior of Nafion. J Membr Sci 2005;264:167e75. [37] Pimenov AV. Water sorption of sulfo-cationic membranes and prediction of their ion-exchange selectivity; 1992. PhD Thesis. S-Petersburg. [38] Zawodzinski T, Derouin C, Radzinski S, Sherman R, Smith T, Springer T, et al. Water uptake by and transport through Nafion 117 membranes. J Electrochem Soc 1993;140:1041e7. [39] Jeck S, Scharfer P, Kind M. Absence of Schroeder’s paradox: experimental evidence for water-swollen Nafion membranes. J Membr Sci 2011;373:74e9. [40] Schroeder P. Ueber Erstarrungs- und Quellungserscheinungen von Gelatine. Z Phys Chem 1903;45:75e117. [41] Choi P, Datta R. Sorption in proton-exchange membranes. An explanation of Schroeder’s paradox. J Electrochem Soc 2003;150:E601e7. [42] Eikerling M, Berg P. Poroelecroelastic theory of water sorption and swelling in polymer electrolyte membranes. Soft Matter 2011;7:5976e90. [43] Bass M, Freger V. An experimental study of Schroeder’s paradox in Nafion and Dowex polymer electrolytes. Desalination 2006;199:277e9. [44] Onishi LM, Prausnitz JM, Newman J. Water-Nafion equilibria. Absence of Schroeder’s paradox. J Phys Chem B 2007;111:10166e73. [45] Freger V. Hydration of ionomers and Schroeder’s paradox in nafion. J Phys Chem B 2009;113:24e36. [46] Davankov VA, Pastukhov AV. Paradoxes of thermodynamics of swelling equilibria of polymers in liquids and vapors. J Phys Chem B 2011;115:15188e95. [47] Volino F, Pineri M, Dianoux AJ, De Geyer A. Water mobility in a water-soaked Nafion membrane. A high-resolution neutron quasielastic study. J Polym Sci Part A-2 Polym Phys 1982;20:481e96. [48] Zawodzinski T, Neeman M, Sillerud L, Gottesfeld S. Determination of water diffusion coefficients in perfluorosulfonate ionomeric membranes. J Phys Chem 1991;95:6040e4. [49] Rivin D, Kendrick CE, Gibson PW, Schneider NS. Solubility and transport behavior of water and alcohols in Nafion. Polymer 2001;42:623e35. [50] Jayakody JRP, Stallworth PE, Mananga ES, Farrington-Zapata J, Greenbaum SG. High pressure NMR study of water self-diffusion in NAFION-117 membrane. J Phys Chem B 2004;108:4260e2. [51] Suresh G, Scindia YM, Pandey AK, Goswami AA. Selfdiffusion coefficient of water in Nafion-117 membrane with different monovalent counterions: a radiotracer study. J Membr Sci 2005;250:39e45. [52] Majsztrik P, Bocarsly A, Benziger J. Water permeation through nafion membranes: the role of water activity. J Phys Chem B 2008;112:16280e9. [53] Luo Z, Chang Z, Zhang Y, Zh Liu, Li J. Electro-osmotic drag coefficient and proton conductivity in Nafion membrane for PEMFC. Int J Hydrogen Energy 2010;35:3120e4. [54] Ise M, Kreuer KD, Maier J. Electroosmotic drag in polymer electrolyte membranes: an electrophoretic NMR study. Solid State Ionics 1999;125:213e23.
[55] Broadbent SR, Hammersley JM. Percolation processes. I. Crystals and mazes. Math Proc Camb Phil Soc 1957;53(03):629e41. [56] Kirkpatrick S. Percolation and conduction. Rev Mod Phys 1973;45:574e8. [57] Efros AL. Fizika i geometriya besporyadka [Physics and geometry of disorder]. Moscow: Nauka; 1982. [58] Roldughin VI, Vysotskii VV. Percolation properties of metalfilled polymer films, structure and mechanisms of conductivity. Progress in Organic Coatings 2000;39:81e100. [59] Shklovskii BI, Efros AL. Theory of occurrence and conductivity of highly heterogeneous media. Usp Fiz Nauk 1975;117:401e35. [60] Sarychev AK, Vinogradoff AP, Karimov AM. Calculation of percolation conductivity in 3D. J Phys C Solid State Phys 1985;18:L105e8. [61] Karpenko LV. Electro-transport properties of the ionexchange membranes as a function of the polymer composition and the concentration of the equilibrium electrolyte solution. PhD Thesis [in Russian]. Krasnodar, Russia: Kuban State University; 1999. [62] Stumper J, Campbell SA, Wilkinson DP, Johnson MC, Davis M. In-situ methods for the determination of current distributions in PEM fuel cells. Electrochim Acta 1998;43:3773e83. [63] Hsing I-M, Futerko P. Two-dimensional simulation of water transport in polymer electrolyte fuel cells. Chem Eng Sci 2000;55:4209e18. [64] Kreuer KD. On the development of proton conducting materials for technological applications. Solid State Ionics 1997;97:1e15. [65] Meier F, Eigenberger G. Transport parameters for the modelling of water transport in ionomer membranes for PEM-fuel cells. Electrochim Acta 2004;49:1731e42. [66] Ren X, Gottesfeld S. Electro-osmotic drag of water in poly (perfluorosulfonic acid) membranes. J Electrochem Soc 2001;148:A87e93. [67] Saito M, Naoko A, Kikuko H, Okada T. Mechanisms of ion and water transport in perfluorosulfonated ionomer membranes for fuel cells. J Phys Chem B 2004;108:16064e70. [68] Edition 11/2010. AVL List GmbH Manual AVL FIRE version 2010; 2010. p. 91. [69] Wasterlain S, Candusso D, Hissel D, Harel F, Bergman P, Menard P, et al. Study of temperature, air dew point temperature and reactant flow effects on proton exchange membrane fuel cell performances using electrochemical spectroscopy and voltammetry techniques. J Power Sources 2010;195:984e93. [70] Al-Baghdadi MARS, Shanad Al-Janabi HAK. Effect of operating parameters on the hygroethermal stresses in proton exchange membranes of fuel cells. Int J Hydrogen Energy 2007;32:4510e22. [71] Jang J-H, Chiu H-C, Yan W-M, Sun W-L. Effects of operating conditions on the performances of individual cell and stack of PEM fuel cell. J Power Sources 2008;180:476e83. [72] Fink C, Fouquet N. Three-dimensional simulation of polymer electrolyte membrane fuel cells with experimental validation. Electrochim Acta 2011;56:10820e31. [73] Baumgartner WR. Characterization of polymer-electrolytemembrane fuel cells. PhD Thesis. Graz University of Technology; 2007. [74] Baumgartner WR, Parz P, Fraser SD, Wallno¨fer E, Hacker V. Polarization study of PEMFC with four reference electrodes at hydrogen starvation conditions. J Power Sources 2008;182:413e21.