- Email: [email protected]
Effects of gas diffusion layer deformation on the transport phenomena and performance of PEM fuel cells with interdigitated flow fields
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 x x x ( 2 0 1 8 ) 1 e1 4
Available online at www.sciencedirect.com
ScienceDirect journal homepage: www.elsevier.com/locate/he
Effects of gas diffusion layer deformation on the transport phenomena and performance of PEM fuel cells with interdigitated flow fields Shian Li a,b, Bengt Sunden b,* a b
Marine Engineering College, Dalian Maritime University, Dalian, China Department of Energy Sciences, Lund University, Lund, Sweden
article info
abstract
Article history:
In this study, a three-dimensional, non-isothermal, two-phase flow mathematical model is
Received 7 December 2017
developed and applied to investigate the effect of the GDL deformation on transport phe-
Received in revised form
nomena and performance of proton exchange membrane (PEM) fuel cells with interdigi-
19 June 2018
tated flow fields. The thickness and porosity of the GDL is decreased after compression,
Accepted 8 July 2018
and the corresponding transport parameters (permeability, mass diffusivity, thermal
Available online xxx
conductivity and electrical conductivity) are affected significantly. The alterations in geometry and transport parameters of the GDL are considered in the mathematical model.
Keywords:
The oxygen concentration, temperature, liquid water saturation and volumetric current
PEM fuel cells
density distributions of PEM fuel cells without compression are investigated and then
Interdigitated flow fields
compared to the PEM fuel cells with various assembly forces. The numerical results show
GDL deformation
that the cell performance is considerably improved with increasing assembly forces.
Numerical modeling
However, the pressure drops in the gas flow channels are also substantially increased. It is
Cell performance
concluded that the assembly force should be as small as possible to decrease the parasitic losses with consideration of gas sealing concern. © 2018 Hydrogen Energy Publications LLC. Published by Elsevier Ltd. All rights reserved.
Introduction Proton exchange membrane (PEM) fuel cells are considered as very promising alternative energy conversion devices with a wide range of applications due to their attractive advantages [1]. Although the commercialization of PEM fuel cells is still limited by the cost, durability and stability issues, great efforts have been taken by researchers to improve the cell performance, accelerate and promote the process of commercialization. The gas flow field design of PEM fuel cells has significant impact on the cell performance. The gas flow channels
fabricated on the current collectors provide and guide the reactant gases to the reaction sites in the catalyst layers (CLs). Many configurations of the flow fields in PEM fuel cells have been extensively investigated by researchers, e.g., parallel flow fields, serpentine flow fields, interdigitated flow fields, and pin-type flow fields [2,3]. Recently, the PEM fuel cells with novel flow fields were proposed and investigated, such as the flow field with in-line and staggered blockages [4], and convergent-divergent serpentine flow fields [5]. For the interdigitated flow fields, the interdigitated channel forces the reactants to flow from one channel to adjacent channels through the gas diffusion layers (GDLs) under the ribs. The GDL is one of the fuel cell components, which
* Corresponding author. n). E-mail address: [email protected] (B. Sunde https://doi.org/10.1016/j.ijhydene.2018.07.064 0360-3199/© 2018 Hydrogen Energy Publications LLC. Published by Elsevier Ltd. All rights reserved. n B, Effects of gas diffusion layer deformation on the transport phenomena and perforPlease cite this article in press as: Li S, Sunde mance of PEM fuel cells with interdigitated flow fields, International Journal of Hydrogen Energy (2018), https://doi.org/10.1016/ j.ijhydene.2018.07.064
2
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 x x x ( 2 0 1 8 ) 1 e1 4
provides pathways for reactants/products, conducts electrical current, and provides mechanical support. Numerical and experimental studies related to PEM fuel cells with interdigitated flow fields can be found in the open literature [6e11]. He et al. [6] numerically investigated the two-phase transport characteristics in the cathode GDL of a PEM fuel cell by using a two-dimensional, isothermal and multicomponent transport model. A three-dimensional, isothermal and single-phase model was developed to analyze the transport and electrochemical reactions in PEM fuel cells with conventional and interdigitated flow fields [7]. It has been concluded that the limiting current density of PEM fuel cells with interdigitated flow fields is substantially improved at high current densities as compared to that of the conventional flow fields. Wang et al. [8] experimentally studied the effects of operating conditions including the operating pressure, cell temperature, gas flow rate and humidification on the performance of PEM fuel cells with interdigitated flow fields. A similar work was carried out by Yan et al. [9] to investigate the influences of operating conditions on performance of PEM fuel cells with conventional and interdigitated flow fields. Huang et al. [10] reported influence of a porosity gradient in GDLs on the performance of PEM fuel cells with interdigitated flow fields by using a threedimensional, non-isothermal and two-phase model. It was concluded that the cell performance improvement was not significant, but the pressure drop was remarkably decreased. Recently, Cooper et al. [11] performed an experimental study on the performance improvement of PEM fuel cells with interdigitated flow fields as the ratio of channel length to width was decreased. In order to avoid gas leakage and minimize the contact resistance between different layers, an assembly force is applied on fuel cell stacks. With compression, the deformation of GDL is observed due to its porous feature, and the geometry and transport parameters are changed. Accordingly, the transport processes and cell performance are affected. Hottinen et al. [12] compared effect of homogeneous and inhomogeneous compressions of GDLs on the local transport characteristics and overall performance of PEM fuel cells. In that study, experimentally determined transport properties of GDLs were adopted in a two-dimensional, nonisothermal and single-phase model. Su et al. [13] performed numerical simulations on PEM fuel cells with compressed GDLs by using a three-dimensional, isothermal and singlephase model. The variations of porosity and permeability of GDLs caused by the assembly force were taken into account. A three-dimensional, isothermal and single-phase model was applied to investigate the performance of PEM fuel cells as deformation of GDLs was presented [14]. The changes of contact resistance between the GDLs and current collectors for different assembly forces were considered in the mathematical model, as well as the porosity of GDLs. Zhou et al. [15] reported the effect of assembly force on liquid water saturation in deformed GDLs and the average current density of PEM fuel cells for different types of GDLs and membranes. In that study, the finite element method was used to obtain the geometry of GDLs for different compressions, and a threedimensional, non-isothermal and two-phase model was applied to describe the complicated transport processes within the PEM fuel cells. The water transport and
performance of PEM fuel cells were experimentally investigated for different clamping forces by Cha et al. [16]. It was concluded that the ohmic resistance was decreased and the mass transport resistance was increased as the clamping forces increased. Wang et al. [17] numerically investigated the performance of PEM fuel cells with non-homogeneous deformation of the GDLs by using a three-dimensional, nonisothermal and single-phase model. The cell performance was significantly affected by the decrease in porosity caused by the GDL deformation. Toghyani et al. [18] studied the impact of clamping pressure on the thermal and electrochemical performance of PEM fuel cells with three different widths of channel. Shi et al. [19] studied the influence of compression on the water management of PEM fuel cells with different GDLs by using a two-dimensional, isothermal and two-phase model. Recently, the performance of PEM fuel cells with interdigitated flow fields and inhomogeneous compression of GDL has been numerically analyzed by using a twodimensional, isothermal and two-phase model [20]. The effect of GDL deformation on transport characteristics and cell performance in PEM fuel cells with parallel channels has already been widely and extensively investigated during the past decade [12e18]. However, studies on the effect of GDL deformation on the performance of PEM fuel cells with interdigitated flow fields are still very few [19,20]. A comprehensive study is needed to improve the understanding of transport characteristics of PEM fuel cells with interdigitated flow fields by using a more complete mathematical model. In the present study, a three-dimensional, non-isothermal, twophase flow model based on the finite volume method has been developed and applied to investigate the effect of GDL deformation on transport phenomena of PEM fuel cells with interdigitated flow fields. In addition, the influences of deformation on geometry and transport parameters are taken into account. The oxygen concentration, temperature, liquid water saturation, volumetric current density, and cell performance of PEM fuels with/without compression are presented and compared.
Model description Schematic illustrations of PEM fuel cells using interdigitated flow fields with/without compression are provided in Fig. 1. It is clear that the deformed GDLs intrude into the flow channel, and the cross-sectional area of the flow channel is decreased. Fig. 2 illustrates the three-dimensional computational domain considered in this study. Due to the symmetry of configurations in the y-coordinate direction only, half inlet and outlet channels are included in the computational domain. The detailed information of the fuel cell geometry and operating condition parameters are summarized in Table 1 and Table 2. The fuel cell operating temperature and pressure are 353 K and 1 atm, respectively. In addition, 100% relative humidity is applied for both the anode and cathode reactant gases. In the PEM fuel cell mathematical model, the fluid flow is laminar; ideal gas law is applied for the reactant gases; the CLs are homogeneous and isotropic; the reactant gases cannot diffuse across the membrane; the generated water in the cathode CL is in dissolved phase [23].
n B, Effects of gas diffusion layer deformation on the transport phenomena and perforPlease cite this article in press as: Li S, Sunde mance of PEM fuel cells with interdigitated flow fields, International Journal of Hydrogen Energy (2018), https://doi.org/10.1016/ j.ijhydene.2018.07.064
3
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 x x x ( 2 0 1 8 ) 1 e1 4
Fig. 1 e Schematic illustrations of PEM fuel cells with interdigitated flow fields: (a) without compression; (b) with compression.
Fig. 2 e Schematic illustrations of the computational domain of a PEM fuel cell with interdigitated flow fields: (a) threedimensional computational domain; (b) representative plane and line.
Table 1 e Fuel cell geometric dimensions [15,21,22]. Parameter Fuel cell length Fuel cell width Flow channel width Flow channel height Rib width Current collector height GDL thickness CL thickness Membrane thickness
Table 2 e Fuel cell operating conditions.
Value
Units
50 2 1 1 1 1.5 0.2 0.01 0.05
mm mm mm mm mm mm mm mm mm
Parameter Reactant gas Operating pressure, Pa/Pc Stoichiometric ratio, xa/xc Relative humidity, RHa/RHc Operating temperature, Ta/Tc
Anode
Cathode
Hydrogen/Water 1 atm 1.5 100% 353 K
Air/Water 1 atm 1.5 100% 353 K
n B, Effects of gas diffusion layer deformation on the transport phenomena and perforPlease cite this article in press as: Li S, Sunde mance of PEM fuel cells with interdigitated flow fields, International Journal of Hydrogen Energy (2018), https://doi.org/10.1016/ j.ijhydene.2018.07.064
4
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 x x x ( 2 0 1 8 ) 1 e1 4
With the above mentioned assumptions, the governing equations of mass, momentum, species, energy, charge, liquid water and dissolved water are presented as follows [15,23]: Mass conservation equation:
where Yi and Deff,i are the mass fraction and effective diffusivity for the i th species, respectively. Si represents the amount of consumption or generation of the i th species. The effective diffusivity (Deff,i) is obtained by the Bruggeman correlation (see Table 3) which accounts for the tortuous path in the GDLs and CLs. Energy conservation equation:
V$ðr! u Þ ¼ Smass
u T ¼ V$ keff VT þ ST V$ rCp !
Governing equations
(1)
where r and ! u are the mixture fluid density and superficial velocity, respectively. Smass represents the source term of the mass conservation equation. The mass changes of hydrogen, oxygen and water vapor due to the processes of electrochemical reactions and phase change are completely considered in the corresponding layers. Momentum conservation equation: V$ðr! u! u Þ ¼ V$ðmV! u Þ VP þ Smom
(4)
where Cp and keff are the specific heat and effective thermal conductivity, respectively. The irreversible, reversible, ohmic heat generation and phase change terms are all included in the energy equation source term, ST. Charge conservation equation: V$ seff ;s Vfs þ Ss ¼ 0
(5)
V$ seff ;m Vfm þ Sm ¼ 0
(6)
(2)
where P and m are the mixture pressure and dynamic viscosity, respectively. Smom denotes the source term of the momentum equation. Species conservation equation: V$ðr! u Yi Þ ¼ V$ rDeff ;i VYi þ Si
(3)
where seff,s is the effective electrical conductivity, seff,m the effective protonic conductivity, fs the electrical potential, fm the protonic potential. The Butler-Volmer equation and spherical agglomerate model are adopted to describe the
Table 3 e Parameters used in the mathematical model. Parameter
Value
mg cm2 kg m3 mg cm2 kg m3 kg m3 kg mol1 mm e e A m2
10(0.03741*T16.96)
A m2
ref
0.5 1 56.4
e e mol m3
[27]
3.39
mol m3
4.56 103 0.101325e(666/Tþ14.1) 100/1.7/0.3 0.25 20000/5000/2000 0.104 326.36 2.36 106 3.517 104 0.0625 110 /95 100 100 1 1.0e-13 9.15 105 2.82 105 2.2 105 2.56 105
Pa m3 mol1 Pa m3 mol1 W m1 K1 W m1 K1 S m1 J mol1 K1 J mol1 K1 J kg1 Pa s N m1 e s1 s1 s1 m2 m2 s1 m2 s1 m2 s1 m2 s1
ref
Anode reference exchange current density, ia [22] ref
Cathode reference exchange current density, ic [25] Anode transfer coefficient, aa [26] Cathode transfer coefficient, ac [26] Reference hydrogen concentration, cH2 [27] Reference oxygen concentration,
ref c O2
Units
0.4 2.145 104 0.6 1.8 103 1.98 103 1.1 1 0.4 0.5 100
Platinum loading, mpt [24] Platinum density, rpt [22] Carbon loading, mc [22] Carbon density, rc [22] Dry membrane density, rm [22] Membrane equivalent weight, Mm [22] Radius of agglomerate, ragg [22] Volume fraction of ionomer, Li [22] Volume fraction of ionomer in agglomerate, Li,agg [22]
Hydrogen Henry's constant, HH2 [23] Oxygen Henry's constant, HO2 [23] Thermal conductivity of CC/GDL/CL, kCC/GDL/CL [22] Thermal conductivity of membrane, km [22] Electrical conductivity of CC/GDL/CL, ss,CC/GDL/CL [22] Entropy of hydrogen oxidation, DSa [28] Entropy of oxygen reduction, DSc [28] Latent heat of condensation/evaporation, Dhlg [29] Liquid water viscosity, ml [29] Surface tension, s [29] Contact angle of GDL/CL, qGDL/CL [29] Condensation rate, gcon [23] Evaporation rate, gevap [23] Dissolved water phase change rate, g [23] Permeability of CL, KCL [22] Binary diffusivity, DH2 H2 O [30] Binary diffusivity, DO2 H2 O [30] Binary diffusivity, DO2 N2 [30] Binary diffusivity, DH2 ON2 [30]
n B, Effects of gas diffusion layer deformation on the transport phenomena and perforPlease cite this article in press as: Li S, Sunde mance of PEM fuel cells with interdigitated flow fields, International Journal of Hydrogen Energy (2018), https://doi.org/10.1016/ j.ijhydene.2018.07.064
5
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 x x x ( 2 0 1 8 ) 1 e1 4
hydrogen oxidation reaction and oxygen reduction reaction in the anode and cathode CLs, respectively. Liquid water transport equation: Krl mg ! u ¼ V$ðrl Ds VsÞ þ Sl V$ rl Krg ml
(7)
where Krl and Krg are the relative permeability of liquid and gas phase, respectively. Ds denotes the capillary diffusion coefficient. The source term, Sl, represents the phase change process between water vapor and liquid water, which is determined by the water vapor partial pressure and saturation pressure. Dissolved water transport equation:
V$
n r d sm Vfm ¼ V$ m Dl Vl þ Sd F Mm
(8)
where nd and Dl are the electro-osmotic drag coefficient and water diffusivity in the membrane, respectively. The generated water in the cathode CL is in dissolved phase, and the absorption/desorption processes of ionomer in the anode and cathode CLs take place between water vapor and dissolved water. The parameters and complementary equations used in the mathematical model are summarized in Tables 3 and 4. In addition, the source terms of the governing equations are given in Table 5.
The catalyst layer composition and volume fraction The catalyst layer consists of void space, Pt/C particles, and ionomer phase. The volume fraction of each component is determined by the following expressions [26,36,37].
Table 4 e Complementary equations and definitions. Description
Units 1:5 1:5
m s1 2
Effective mass diffusivity [31]
Deff ;i ¼ ð1 sÞ
Mass diffusivity [17]
1 Xi j¼1;jsi Xj =Di;j 1:5 P0 T Di;j ¼ Di;j ðT0 ; P0 Þ T0 P keff ¼ ð1 εÞks þ εkf !0:5 PH ref ja ¼ ð1 sÞia aeff ref 2 ½eaa Fha =RT eac Fha =RT c H 2 HH 2
Binary mass diffusivity [30] Effective thermal conductivity [21] Electrochemical kinetics [31]
Over-potential [31] Open circuit voltage [22] Active surface area [26] Proton conductivity [32] Effective conductivity [22]
ε
Di;m
m2 s1
Di;m ¼ Pn
m2 s1 W m1 K1 A m3
ha ¼ fs fm ; hc ¼ fs fm Voc
V V
Voc ¼ 1:229 8:456 104 ðT 298:15Þ þ 4:31 105 TLnðPH2 P0:5 O2 Þ mpt 3 2 3 ð227:79f 158:57f 201:53f þ 159:5Þ 10 aeff ¼ tCL 1 1 Þ 1268ð 303 T sm ¼ ð0:514l 0:326Þe
m1 S m1 S m1
seff ;s ¼ ð1 εGDL Þ1:5 ss ; seff ;s ¼ ð1 εCL Li Þ1:5 ss seff ;m ¼
Li1:5 sm
Relative permeability [23]
Krl ¼ Ks3 ; Krg ¼ Kð1 sÞ3
m2
Capillary diffusivity [23]
Ks3 dPc Ds ¼ ml ds ε 0:5 ð1:417s 2:12s2 þ 1:263s3 Þ Pc ¼ scosðqÞ K log10 Psat ¼ 2:1794 þ 0:02953ðT 273:15Þ 9:1837 105 ðT 273:15Þ2 þ
m2 s1
Capillary pressure [23] Saturation pressure [23] Electro-osmotic drag coefficient [32] Dissolved water diffusivity [32]
Equilibrium water content [33] Water activity [23] Oxygen diffusivity in liquid water [34]
1:4454 107 ðT 273:15Þ3 l nd ¼ 2:5 22 8 2:05l 3:25 1 1 < 10 ½2416ð303 TÞ 6:65 1:25l Dl ¼ 10 e : 2:563 0:33l þ 0:0264l2 0:000671l3 1:41 þ 11:3a 18:8a2 þ 16:2a3 ða < 1Þ leq ¼ 10:1 þ 2:94ða 1Þ ða 1Þ Xwv P þ 2s a ¼ Psat DO2 ;w ¼ 7:4 1012
TðjMH2 O Þ0:5 mH2 O VO0:62
Oxygen diffusivity in ionomer [23] DO2 ;i ¼ 2:88 10 Permeability of the GDL [35] Inertial coefficient [35]
K ¼
10
1 T1 ½2933ð313 e
Pa Pa e
ð2≪l < 3Þ ð3≪l < 4Þ ð4 < lÞ
m2 s1
e e m2 s1 m2 s1 m2
d2f ε3 2
16kCK ð1 εÞ
b ¼ 2:88 106
1 ε1:5 K
m1
n B, Effects of gas diffusion layer deformation on the transport phenomena and perforPlease cite this article in press as: Li S, Sunde mance of PEM fuel cells with interdigitated flow fields, International Journal of Hydrogen Energy (2018), https://doi.org/10.1016/ j.ijhydene.2018.07.064
6
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 x x x ( 2 0 1 8 ) 1 e1 4
Table 5 e Source terms in the governing equations [23]. Description
Units kg m3 s1
Smass ¼ SH2 þ Swv Anode CL Smass ¼ SO2 þ Swv Cathode CL Smass ¼ Swv Anode and cathode GDLs m ! Anode and cathode GDLs and CLs u þ brj! u j! u Smom ¼ K ja SH2 ¼ MH2 Anode CL 2F jc SO2 ¼ MO2 Cathode CL 4F Swv ¼ Sl Svd MH2 O Anode and cathode CLs Swv ¼ Sl Anode and cathode GDLs TDSa ST ¼ j a ha ja þ seff ;m kVfm k2 þ seff ;s kVfs k2 þ Sl Dhlg Anode CL 2F TDSc jc þ seff ;m kVfm k2 þ seff ;s kVfs k2 þ Sl Dhlg Cathode CL ST ¼ j c hc 4F ST ¼ seff ;m kVfm k2 membrane ST ¼ seff ;s kVfs k2 þ Sl Dhlg ST ¼ seff ;s kVfs k2
kg m2 s2 kg m3 s1 kg m3 s1 kg m3 s1 W m3
Anode and cathode GDLs
Anode and cathode CCs A m3
Ss ¼ ja Anode CL Ss ¼ þ jc Cathode CL Sm ¼ þ ja Anode CL Sm ¼ jc Cathode CL Sl ¼ Sphase Anode and cathode GDLs and CLs 8 εð1 sÞ > > MH2 O ðPwv Psat Þ Pwv Psat < gcond RT Sphase ¼ > εs > : gevap MH2 O ðPwv Psat Þ Pwv < Psat RT Sd ¼ Svd Anode CL Sd ¼ Svd þ Sl Cathode CL r Svd ¼ g m ðleq lÞ Anode and cathode CLs Mm jc Cathode CL Sl ¼ 2F
A m3 kg m3 s1 kg m3 s1
mol m3 s1
The volume fraction of the void space (εCL) is given by: εCL ¼ 1 LPt=C Li
(9)
where LPt/C and Li are the volume fractions of Pt/C particles and ionomer phase, respectively. The volume fraction of the Pt/C particles (LPt/C) is determined by: LPt=C ¼
mpt tCL
1 1f 1 þ rpt f rc
!
mpt mpt þ mc
(11)
where mc is the carbon loading. The volume fraction of the ionomer phase (Li) is expressed by: Li ¼
r3agg
h
LPt=C r3agg Li;agg þ 1 Li;agg
di ¼ ragg
ffi # "sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 3 Li 1 Li;agg Li;agg þ 1 1 LPt=C
(13)
The total volume of the liquid water is: (10)
where mpt is the platinum loading, tcl the thickness of the CL, rpt the density of platinum, rc the density of carbon. The parameter f is defined as the platinum loading divided by the sum of platinum loading and carbon loading. f¼
thickness di can be determined. It is assumed that the agglomerate is only occupied by Pt/C particles and ionomer phase. The thickness of ionomer film (di) is obtained by:
raggþdi
3
r3agg
i
Vw ¼ sεCL VCL
(14)
where s is the liquid water saturation, VCL is the total volume of the CL. The volume of the liquid water covering the individual agglomerate can be obtained: Vw;i ¼
sεCL N
(15)
where N is the number of agglomerate particles per catalyst layer volume. N¼
3LPt=C 4pr3agg 1 Li;agg
(16)
(12)
The ionomer phase consists of two parts: one part is inside the agglomerate, the other part is the ionomer film. When LPt/C, Li, ragg, and Li,agg are given, a unique value of the ionomer
The thickness of the liquid water film (dw) is obtained by: rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 3 3sεCL 3 ragg þ di þ dw ¼ ragg þ di 4pN
(17)
n B, Effects of gas diffusion layer deformation on the transport phenomena and perforPlease cite this article in press as: Li S, Sunde mance of PEM fuel cells with interdigitated flow fields, International Journal of Hydrogen Energy (2018), https://doi.org/10.1016/ j.ijhydene.2018.07.064
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 x x x ( 2 0 1 8 ) 1 e1 4
The cathode agglomerate model In the present analyses, the spherical agglomerate model [37,38] is used to calculate the volumetric current density. Fig. 3 depicts the structure of the agglomerates within the cathode CL. The identical agglomerates are uniformly distributed in the cathode CL, and the individual agglomerate is evenly covered by ionomer and liquid water films.
jc ¼ 4F
ragg þ di þ dw PO2 1 di dw þ þ HO2 Er kc ð1 εCL Þ ragg aagg;i DO2 ;i aagg;w DO2 ;w
1 (18)
where PO2 is the partial pressure of oxygen, HO2 the Henry's constant, Er the effective factor, kc the reaction rate constant, aagg,i the ionomer effective agglomerate surface area, and aagg,w the liquid water effective agglomerate surface area. The effectiveness factor of the spherical agglomerate is defined as: Er ¼
1 1 1 FL tanhð3FL Þ 3FL
(19)
where the Thiele's modulus FL for chemical reactions is given as: FL ¼
ragg 3
sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi kc L1:5 i;agg DO2
(20)
The reaction rate constant is computed as:
ref
kc ¼
ic aeff 4Fð1
ref εCL ÞcO2
eaa Fhc =RT þ eac Fhc =RT
(21)
The effective agglomerate surface area is obtained from: aagg;i ¼
2 3LPt=C εCL ragg þ di r3agg 1 Li;agg
aagg;w ¼
2 3LPt=C εCL ragg þ di þ dw r3agg 1 Li;agg
(22)
(23)
7
and non-uniformly distributed due to the deformation resulting from the assembly force. It is assumed that the thickness and porosity are decreased only caused by the change in pore size. Therefore, the GDL porosity is determined by the following expression [12,15]: ε ¼ 1 ð1 ε0 Þ
d0 d
(24)
where ε0 is the initial porosity, d0 the initial thickness, d the thickness after compression, ε the porosity after compression. The local porosity is varying with the thickness of the GDL, and the transport parameters associated with porosity are also changed.
Numerical implementation The developed mathematical model for PEM fuel cells is implemented in the commercial software ANSYS FLUENT. The governing equations of charge, liquid water transport, and dissolved water transport are implemented by using user defined scalar (UDS) equations with an under relaxation technique developed by the authors.
Boundary conditions At the inlet of the flow channels, the mass flow rate, temperature, and species mass fractions are prescribed. The mass flow rates of reactants at both anode and cathode sides are calculated at a reference current density of 1.0 A/cm2 and are given by Eqs. (25) and (26) The anode and cathode stoichiometric ratios (xa and xc) are provided in Table 2. In addition, the inlet liquid water saturation is assigned as zero. At the outlet of the flow channels, a pressure-outlet boundary condition is applied. The operating temperature and a constant electric potential, fs ¼ 0, are specified at the anode terminal. At the cathode terminal, the operating temperature and a constant electric potential, fs ¼ Vcell, are applied. A symmetry boundary condition is adopted at the outer surface of the x-z plane.
Determination of porosity of the GDLs after compression
Qa ¼
xa MH2 Iref Am 2FYH2
(25)
For simplicity, a constant and uniform porosity is adopted in the numerical investigation. However, the porosity is varied
Qc ¼
xc MO2 Iref Am 4FYO2
(26)
Fig. 3 e Schematic illustration of the ionomer and liquid water films covering the agglomerate.
n B, Effects of gas diffusion layer deformation on the transport phenomena and perforPlease cite this article in press as: Li S, Sunde mance of PEM fuel cells with interdigitated flow fields, International Journal of Hydrogen Energy (2018), https://doi.org/10.1016/ j.ijhydene.2018.07.064
8
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 x x x ( 2 0 1 8 ) 1 e1 4
The algorithm Semi-Implicit Method for Pressure Linked Equations (SIMPLE) is used to handle the pressure-velocity coupling. The second order upwind discretization scheme is adopted for momentum, species, energy and UDS equations. A grid independence study was carefully performed to balance the accuracy and computational resources. A mesh system with the total number of control volumes of 208,000 was selected in the following simulations. To assess the accuracy of the numerical model, the polarization curve obtained from the model is widely adopted to compare with the experimental data in the open literature. The numerical performance of the developed model are evaluated by comparison with the experimental data reported by Yan et al. [39], as shown in Fig. 4. It is found that the numerical method can reproduce the experimental data well. The information of geometry and operating parameters given in Ref. [39] is relatively complete, and therefore it is selected for model validation. For the experimental study, a 25 cm2 fuel cell with a platinum loading of 0.4 mg/cm2 for both the anode and cathode CLs was tested under the operating temperature of 353 K, operating pressure of 1 atm, and a relatively humidity of 100%. In addition, the present model has been applied to investigate the effects of agglomerate model parameters on local transport characteristics and overall cell performance [40].
because the Young's modulus of the GDLs is much smaller than that of the remaining components and the CLs are very thin. The thickness and porosity of the GDLs under different assembly forces are shown in Fig. 5. As shown in Fig. 5a, the thickness of the GDL under the rib is substantially decreased with increasing assembly forces, while the thickness of the GDL under the channel middle region is not affected by the assembly forces. In addition, the local porosity of the deformed GDLs which is determined by Eq. (24), is given in Fig. 5b. When the assembly force is increased from 0 MPa to 2.0 MPa, the thickness is decreased from 0.28 mm to 0.13 mm and the porosity is decreased from 0.78 to 0.53 at the rib region. Note that the variations of thickness and porosity of GDLs are considered in the y-coordinate direction and the porosity is a function of the y-coordinate. The cell performance is commonly presented and compared in terms of the polarization and power density curves. Fig. 6 shows the effect of assembly forces on the performance of PEM fuel cells with interdigitated flow fields.
Results and discussions In the preceding section, the PEM fuel cell mathematical model and GDL deformation were described in detail. In this section, effect of compression on the thickness and porosity of the GDLs is presented. Then the current densities and power densities of PEM fuel cells with/without compression are investigated and the local transport characteristics at an operating voltage 0.7 V are also analyzed. For simplicity, the configurations of the deformed GDLs under different assembly forces reported by Zhou et al. [15] are employed in the present study. It is assumed that the deformation after compression only takes place in the GDLs,
Fig. 4 e Comparison of the numerical results with the experimental data [39].
Fig. 5 e Distribution of the GDL characteristics under different assembly forces: (a) thickness; (b) porosity.
Fig. 6 e Polarization and power density curves of PEM fuel cells under different assembly forces.
n B, Effects of gas diffusion layer deformation on the transport phenomena and perforPlease cite this article in press as: Li S, Sunde mance of PEM fuel cells with interdigitated flow fields, International Journal of Hydrogen Energy (2018), https://doi.org/10.1016/ j.ijhydene.2018.07.064
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 x x x ( 2 0 1 8 ) 1 e1 4
It can be seen that the current density is increased with increasing assembly force. For an operating voltage 0.7 V, the current densities for four different cases are 0.2383 A/cm2, 0.2422 A/cm2, 0.2565 A/cm2 and 0.2575 A/cm2, respectively. In addition, the maximum power densities are 0.4227 W/ cm2, 0.4263 W/cm2, 0.4464 W/cm2 and 0.4503 W/cm2, respectively. In a similar work carried out by Mahmoudi et al. [20], the cell performance was decreased as the compression ratio was increased. Note that a constant pressure drop between the inlet and outlet channels were adopted in that study. He et al. [6] numerically investigated the influence of pressure drop on the cell performance. It was concluded that the cell performance was substantially improved with increasing pressure drop. The pressure drops between the cathode inlet and outlet flow channels for the four cases at the operating voltage 0.7 V are 146.2 Pa, 428.4 Pa, 1111.3Pa and 3104.1 Pa, respectively. The pressure drop of the PEM fuel cell with 2.0 MPa assembly force is approximately 21.2 times higher than that of the PEM fuel
9
cell without an assembly force. This is because the mass flow rate is constant for all cases in this study, and the thickness, porosity and permeability of the GDLs under the rib regions and the cross sectional area of the flow channels are decreased due to the deformation of the GDLs caused by the assembly force. The distributions of oxygen concentration, temperature and liquid water saturation at the cathode GDL-CL interface plane and the volumetric current density at the middle plane of cathode CL of the PEM fuel cell without compression are given in Fig. 7. It is clear that the oxygen concentration is gradually decreased from the inlet region to the outlet region in Fig. 7a. As shown in Fig. 7b, the temperature is slightly increased along the flow direction, and the minimum temperature is observed under the rib region. In Fig. 7c, the liquid water saturation is gradually increased along the flow direction due to the consumption of oxygen and generation of water. In addition, it is observed that the volumetric current density is gradually decreased along the flow direction and the
Fig. 7 e Transport characteristics of PEM fuel cells without compression: (a) oxygen concentration distribution at the cathode GDL-CL interface; (b) temperature distribution at the cathode GDL-CL interface; (c) liquid water saturation distribution at the cathode GDL-CL interface (c) volumetric current density distribution at the middle plane of cathode CL. n B, Effects of gas diffusion layer deformation on the transport phenomena and perforPlease cite this article in press as: Li S, Sunde mance of PEM fuel cells with interdigitated flow fields, International Journal of Hydrogen Energy (2018), https://doi.org/10.1016/ j.ijhydene.2018.07.064
10
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 x x x ( 2 0 1 8 ) 1 e1 4
Fig. 8 e Oxygen concentration distribution at the middle plane of cathode GDL of PEM fuel cells under different assembly forces: (a) without compression; (b) 1.0 MPa; (c) 1.5 MPa; (d) 2.0 MPa.
Fig. 9 e Temperature distribution at the middle plane of cathode GDL of PEM fuel cells under different assembly forces: (a) without compression; (b) 1.0 MPa; (c) 1.5 MPa; (d) 2.0 MPa. n B, Effects of gas diffusion layer deformation on the transport phenomena and perforPlease cite this article in press as: Li S, Sunde mance of PEM fuel cells with interdigitated flow fields, International Journal of Hydrogen Energy (2018), https://doi.org/10.1016/ j.ijhydene.2018.07.064
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 x x x ( 2 0 1 8 ) 1 e1 4
maximum volumetric current density appears under the rib region and close to the inlet channel side. The oxygen concentration distributions in the middle plane of the GDLs of PEM fuel cells under different assembly
11
forces are shown in Fig. 8. It is observed that the oxygen concentration gradually decreases from the inlet region to the outlet region for all cases due to the consumption of oxygen by the oxygen reduction reaction. The maximum oxygen
Fig. 10 e Liquid water saturation distribution at the middle plane of cathode GDL of PEM fuel cells under different assembly forces: (a) without compression; (b) 1.0 MPa; (c) 1.5 MPa; (d) 2.0 MPa.
Fig. 11 e Volumetric current density of PEM fuel cells under different assembly forces. n B, Effects of gas diffusion layer deformation on the transport phenomena and perforPlease cite this article in press as: Li S, Sunde mance of PEM fuel cells with interdigitated flow fields, International Journal of Hydrogen Energy (2018), https://doi.org/10.1016/ j.ijhydene.2018.07.064
12
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 x x x ( 2 0 1 8 ) 1 e1 4
concentration appearing at the region close to the inlet flow channel is increased when the assembly force is gradually increased. This can be explained by the increase of velocity in the flow channel, which is caused by the decrease of cross sectional area of the flow channel with increasing assembly forces. It is also observed that the region of minimum oxygen concentration in the outlet flow channel is enlarged. The temperature distributions in the middle plane of the GDLs of PEM fuel cells under different assembly forces are presented in Fig. 9. It is clear that the minimum temperature and maximum temperature appear under the rib region and flow channel region, respectively. Due to the presence of rib, the generated heat in CLs and GDLs can be easily and effectively transferred through it. When the assembly force is applied, the thickness and porosity of the GDLs are both decreased and the effective thermal conductivity is consequently increased. Accordingly, the temperature under the rib region is gradually decreased with increasing assembly forces. The liquid water saturation distributions in the middle plane of the GDLs of PEM fuel cells under different assembly forces are shown in Fig. 10. It can be seen that the maximum liquid water saturation appears under the rib regions for all cases, and it is gradually increased with increasing assembly forces. This is because the porosity and permeability of the GDLs are both decreased under a higher assembly force. In addition, a lower local temperature is also observed at this region as shown in Fig. 9. The volumetric current density distributions at the middle plane of CLs (Line-1) are shown in Fig. 11. The minimum and maximum volumetric current densities appear at the outlet channel region and under the rib region for all cases, respectively. The volumetric current density is gradually increased when the assembly force is increased from 0 Mpa to 1.5 Mpa. When the assembly force is increased from 1.5 MPa to 2.0 MPa, the maximum volumetric current density is increased and the minimum volumetric current density is decreased. This indicates that a higher variation of the volumetric current density is obtained under the assembly force of 2.0 MPa.
Conclusions In this study, the effect of GDL deformation caused by the assembly force on transport characteristics and performance of PEM fuel cells with interdigitated flow fields were numerically investigated by using a three-dimensional, nonisothermal, two-phase flow mathematical model based on the finite volume method. The following conclusions were drawn: The thickness and porosity of GDLs are decreased due to the deformation. The corresponding transport parameters are altered and consequently the cell performance is also affected. The oxygen concentration is increased and decreased at the inlet channel and outlet channel regions, respectively. The temperature under the rib region is decreased due to the increased effective thermal conductivity. The liquid water saturation under the rib region is increased due to the decreased local porosity and temperature. In addition, the cell performance improvement is obtained and accompanied with a relatively high pressure drop penalty. It is
concluded that the assembly force should be as small as possible to decrease the parasitic losses with consideration of the gas sealing problem. The predicted results provided detailed transport characteristics of PEM fuel cells with interdigitated flow fields under compression, and may serve as guide lines to improve the cell performance and enable proper design. The present numerical model can be further developed with consideration of the anisotropic properties in GDLs and the contact resistance between different layers. In addition, different types of GDLs can be further studied and compared, as well as the fuel cells with/without micro-porous layer.
Acknowledgments The work was carried out at the Department of Energy Sciences, Lund University. The first author gratefully acknowledges the financial support from China Scholarship Council (CSC).
references
[1] Wang Y, Chen KS, Mishler J, Cho SC, Adroher XC. A review of polymer electrolyte membrane fuel cells: technology, applications, and needs on fundamental research. Appl Energy 2011;88:981e1007. [2] Li XG, Sabir I. Review of bipolar plates in PEM fuel cells: flow field design. Int J Hydrogen Energy 2005;30:359e71. [3] Manso AP, Marzo FF, Barranco J, Garikano X, Mujika MG. Influence of geometric parameters of the flow fields on the performance of a PEM fuel cell. Int J Hydrogen Energy 2012;37:15256e87. [4] Heidary H, Kermani MJ, Advani SG, Prasad AK. Experimental investigation of in-line and staggered blockages in parallel flow field channels of PEM fuel cells. Int J Hydrogen Energy 2016;41:6885e93. [5] Chowdhury MZ, Akansu YE. Novel convergent-divergent serpentine flow fields effect on PEM fuel cell performance. Int J Hydrogen Energy 2017;42:25686e94. [6] He WS, Yi JS, Nguyen TV. Two-phase flow model of the cathode of PEM fuel cells using interdigitated flow fields. AIChE J 2000;46:2053e64. [7] Um S, Wang CY. Three dimensional analysis of transport and electrochemical reactions in polymer electrolyte fuel cells. J Power Sources 2004;125:40e51. [8] Wang L, Liu HT. Performance studies of PEM fuel cells with interdigitated flow fields. J Power Sources 2004;134:185e96. [9] Yan WM, Chen CY, Mei SC, Soong CY, Chen FL. Effects of operating conditions on cell performance of PEM fuel cells with conventional or interdigitated flow field. J Power Sources 2006;162:1157e64. [10] Huang YX, Cheng CH, Wang XD, Jang JY. Effects of porosity gradient in gas diffusion layers on performance of proton exchange membrane fuel cells. Energy 2010;35:4786e94. [11] Cooper NJ, Santamaria AD, Becton MK, Park JW. Investigation of the performance improvement in decreasing aspect ratio interdigitated flow field PEMFCs. Energy Convers Manag 2017;136:307e17. [12] Nitta I, Karvonen S, Himanen O, Mikkola M. Modelling the effect of inhomogeneous compression of GDL on local transport phenomena in a PEM fuel cell. Fuel Cell 2008;8:410e21.
n B, Effects of gas diffusion layer deformation on the transport phenomena and perforPlease cite this article in press as: Li S, Sunde mance of PEM fuel cells with interdigitated flow fields, International Journal of Hydrogen Energy (2018), https://doi.org/10.1016/ j.ijhydene.2018.07.064
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 x x x ( 2 0 1 8 ) 1 e1 4
[13] Su ZY, Liu CT, Chang HP, Li CH, Huang KJ, Sui PC. A numerical investigation of the effects of compression force on PEM fuel cell performance. J Power Sources 2008;183:182e92. [14] Taymaz I, Benli M. Numerical study of assembly pressure effect on the performance of proton exchange membrane fuel cell. Energy 2010;35:2134e40. [15] Zhou YB, Jiao K, Du Q, Yin Y, Li XG. Gas diffusion layer deformation and its effect on the transport characteristics and performance of proton exchange membrane fuel cell. Int J Hydrogen Energy 2013;39:12891e903. [16] Cha D, Ahn JH, Kim HS, Kim Y. Effects of clamping force on the water transport and performance of a PEM fuel cell with relative humidity and current density. Energy 2015;93:1338e44. [17] Wang JT, Yuan JL, Yu JS, Sunden B. Investigation of effects of non-homogeneous deformation of gas diffusion layer in a PEM fuel cell. Int J Energy Res 2017;41:2121e37. [18] Toghyani S, Moradi Nafchi F, Afshari E, Hasanpour K, Baniasadi E, Atyabi SA. Thermal and electrochemical performance analysis of a proton exchange membrane fuel cell under assembly pressure on gas diffusion layer. Int J Hydrogen Energy 2018;43:4534e45. [19] Shi ZY, Wang X, Guessous L. Effect of compression on the water management of a proton exchange membrane fuel cell with different gas diffusion layers. J Fuel Cell Sci Technol 2010;7:021012. [20] Mahmoudi AH, Ramiar A, Esmaili Q. Effect of inhomogeneous compression of gas diffusion layer on the performance of PEMFC with interdigitated flow field. Energy Convers Manag 2016;110:78e89. [21] Li SA, Yuan JL, Andersson M, Xie GN, Sunden B. Wavy surface cathode gas flow channel effects on transport processes in a proton exchange membrane fuel cell. ASME J Electrochem En Conv Stor 2017;14:031007. [22] Li SA, Yuan JL, Andersson M, Xie GN, Sunden B. Influence of anisotropic gas diffusion layers on transport phenomena in a proton exchange membrane fuel cell. Int J Energy Res 2017;41:2034e50. [23] Jiao K, Li XG. Water transport in polymer electrolyte membrane fuel cells. Prog Energy Combust Sci 2011;37:221e91. [24] Wang L, Husar A, Zhou TH, Liu HT. A parametric study of PEM fuel cell performance. Int J Hydrogen Energy 2003;28:1263e72. [25] Parthasarathy A, Srinivasan S, Appleby AJ, Martin CR. Temperature dependence of the electrode kinetics of oxygen reduction at the Platinum/Nafion interface-a microelectrode investigation. J Electrochem Soc 1992;139:2530e7. [26] Dalasm NKH, Kermani MJ, Moghaddam DG, Stockie JM. A parametric study of cathode catalyst layer structural parameters on the performance of a PEM fuel cell. Int J Hydrogen Energy 2010;35:2417e27. [27] Bernardi DM, Verbrugge MW. A mathematical model of the solid polymer electrolyte fuel cell. J Electrochem Soc 1992;139:2477e91. [28] Lampinen MJ, Fomino M. Analysis of free energy and entropy changes for half cell reactions. J Electrochem Soc 1993;140(12):3537e46. [29] Yang XG, Ye Q, Cheng P. Matching of water and temperature fields in proton exchange membrane fuel cells with nonuniform distributions. Int J Hydrogen Energy 2011;36:12524e37. [30] Berning T, Lu DM, Djilali N. Three-dimensional computational analysis of transport phenomena in a PEM fuel cell. J Power Sources 2002;106:284e94.
13
[31] Meng H. Numerical studies of liquid water behaviors in PEM fuel cell cathode considering transport across different porous media. Int J Hydrogen Energy 2010;35:5569e79. [32] Springer TE, Zawodzinski TA, Gottesfeld S. Polymer electrolyte fuel cell model. J Electrochem Soc 1991;138(8):2334e42. [33] Zawodzinski TA, Springer TE, Uribe F, Gottesfeld S. Characterization of polymer electrolytes for fuel cell applications. Solid State Ionics 1993;60:199e211. [34] Wilke AJ, Chang P. Correlation of diffusion coefficients in dilute solutions. AIChE J 1995;1:264e70. [35] Gostick JT, Fowler MW, Pritzker MD, Ioannidis MA, Behra LM. In-plane and through-plane gas permeability of carbon fiber electrode backing layers. J Power Sources 2006;162:228e38. [36] Secanell M, Karan K, Suleman A, Djilali N. Multi-variable optimization of PEMFC cathodes using an agglomerate model. Electrochim Acta 2007;52:6318e37. [37] Shah AA, Kim GS, Sui PC, Harvey D. Transient nonisothermal model of a polymer electrolyte fuel cell. J Power Sources 2007;163:793e806. [38] Sun W, Peppley B, Karan K. An improved two-dimensional agglomerate cathode model to study the influence of catalyst layer structural parameters. Electrochim Acta 2005;50:3359e74. [39] Yan QG, Toghiani H, Causey H. Steady state and dynamic performance of proton exchange membrane fuel cells under various operating conditions and load changes. J Power Sources 2006;161:492e502. [40] Li SA, Yuan JL, Xie GN, Sunden B. Effects of agglomerate model parameters on transport characterization and performance of PEM fuel cells. Int J Hydrogen Energy 2018;43:8451e63.
Nomenclature a: Effective surface area, m1/Water activity A: Area, m2 c: Mole concentration, mol m3 CP: Specific heat, J kg1 K1 D: Diffusivity, m2 s1 E: Effectiveness factor f: Platinum mass ratio to Pt/C F: Faraday constant, 96485 C mol1 h: Enthalpy change, J kg1 H: Henry's constant, Pa m3 mol1 i: Exchange current density, A m2 j: Volumetric current density, A m3 k: Thermal conductivity, W m1 K1/Reaction rate constant, s1 K: Permeability, m2 L: Volume fraction m: Mass loading, kg m2 M: Molecular weight, kg mol1 n: Osmotic-drag coefficient P: Pressure, Pa Q: Mass flow rate, kg s1 r: Agglomerate radius, m R: Universal gas constant, 8.314 J mol1 K1 s: Liquid water saturation S: Source term/Entropy, J mol1 K1 t: Thickness, m T: Temperature, K ! u : Velocity vector, m/s V: Voltage, V/Volume, m3 X: Mole fraction Y: Mass fraction
n B, Effects of gas diffusion layer deformation on the transport phenomena and perforPlease cite this article in press as: Li S, Sunde mance of PEM fuel cells with interdigitated flow fields, International Journal of Hydrogen Energy (2018), https://doi.org/10.1016/ j.ijhydene.2018.07.064
14
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 x x x ( 2 0 1 8 ) 1 e1 4
Greek Symbols a: Transfer coefficient b: inertial coefficient, m1 g: Water phase change rate, s1 d: Thickness, m ε: Porosity h: Over-potential, V q: Contact angle l: Water content m: Dynamic viscosity, Pa s x: Stoichiometric ratio r: Density, kg m3 s: Electron/Proton conductivity, S m1/Surface tension, N m2 f: Potential, V F: Theile's modulus Subscripts and Superscripts a: Anode agg: Agglomerate c: Cathode/Carbon con: Condensation d: Osmotic-drag/Dissolved water eff: Effective eq: Equilibrium
evap: Evaporation g: Gas i: i th species/Ionomer j: j th species l: Liquid lg: Liquid and gas water m: Membrane/Mixture mom: Momentum equation oc: Open circuit ref: Reference rl: Liquid relative rg: Gas relative phase: Phase change Pt: Platinum sat: Saturation s: Solid/Liquid water saturation vd: Vapor and dissolved water w: Water wv: Water vapor Abbreviations CC: Current collector CL: Catalyst layer GDL: Gas diffusion layer PEM: Proton exchange membrane
n B, Effects of gas diffusion layer deformation on the transport phenomena and perforPlease cite this article in press as: Li S, Sunde mance of PEM fuel cells with interdigitated flow fields, International Journal of Hydrogen Energy (2018), https://doi.org/10.1016/ j.ijhydene.2018.07.064
Available online at www.sciencedirect.com
ScienceDirect journal homepage: www.elsevier.com/locate/he
Effects of gas diffusion layer deformation on the transport phenomena and performance of PEM fuel cells with interdigitated flow fields Shian Li a,b, Bengt Sunden b,* a b
Marine Engineering College, Dalian Maritime University, Dalian, China Department of Energy Sciences, Lund University, Lund, Sweden
article info
abstract
Article history:
In this study, a three-dimensional, non-isothermal, two-phase flow mathematical model is
Received 7 December 2017
developed and applied to investigate the effect of the GDL deformation on transport phe-
Received in revised form
nomena and performance of proton exchange membrane (PEM) fuel cells with interdigi-
19 June 2018
tated flow fields. The thickness and porosity of the GDL is decreased after compression,
Accepted 8 July 2018
and the corresponding transport parameters (permeability, mass diffusivity, thermal
Available online xxx
conductivity and electrical conductivity) are affected significantly. The alterations in geometry and transport parameters of the GDL are considered in the mathematical model.
Keywords:
The oxygen concentration, temperature, liquid water saturation and volumetric current
PEM fuel cells
density distributions of PEM fuel cells without compression are investigated and then
Interdigitated flow fields
compared to the PEM fuel cells with various assembly forces. The numerical results show
GDL deformation
that the cell performance is considerably improved with increasing assembly forces.
Numerical modeling
However, the pressure drops in the gas flow channels are also substantially increased. It is
Cell performance
concluded that the assembly force should be as small as possible to decrease the parasitic losses with consideration of gas sealing concern. © 2018 Hydrogen Energy Publications LLC. Published by Elsevier Ltd. All rights reserved.
Introduction Proton exchange membrane (PEM) fuel cells are considered as very promising alternative energy conversion devices with a wide range of applications due to their attractive advantages [1]. Although the commercialization of PEM fuel cells is still limited by the cost, durability and stability issues, great efforts have been taken by researchers to improve the cell performance, accelerate and promote the process of commercialization. The gas flow field design of PEM fuel cells has significant impact on the cell performance. The gas flow channels
fabricated on the current collectors provide and guide the reactant gases to the reaction sites in the catalyst layers (CLs). Many configurations of the flow fields in PEM fuel cells have been extensively investigated by researchers, e.g., parallel flow fields, serpentine flow fields, interdigitated flow fields, and pin-type flow fields [2,3]. Recently, the PEM fuel cells with novel flow fields were proposed and investigated, such as the flow field with in-line and staggered blockages [4], and convergent-divergent serpentine flow fields [5]. For the interdigitated flow fields, the interdigitated channel forces the reactants to flow from one channel to adjacent channels through the gas diffusion layers (GDLs) under the ribs. The GDL is one of the fuel cell components, which
* Corresponding author. n). E-mail address: [email protected] (B. Sunde https://doi.org/10.1016/j.ijhydene.2018.07.064 0360-3199/© 2018 Hydrogen Energy Publications LLC. Published by Elsevier Ltd. All rights reserved. n B, Effects of gas diffusion layer deformation on the transport phenomena and perforPlease cite this article in press as: Li S, Sunde mance of PEM fuel cells with interdigitated flow fields, International Journal of Hydrogen Energy (2018), https://doi.org/10.1016/ j.ijhydene.2018.07.064
2
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 x x x ( 2 0 1 8 ) 1 e1 4
provides pathways for reactants/products, conducts electrical current, and provides mechanical support. Numerical and experimental studies related to PEM fuel cells with interdigitated flow fields can be found in the open literature [6e11]. He et al. [6] numerically investigated the two-phase transport characteristics in the cathode GDL of a PEM fuel cell by using a two-dimensional, isothermal and multicomponent transport model. A three-dimensional, isothermal and single-phase model was developed to analyze the transport and electrochemical reactions in PEM fuel cells with conventional and interdigitated flow fields [7]. It has been concluded that the limiting current density of PEM fuel cells with interdigitated flow fields is substantially improved at high current densities as compared to that of the conventional flow fields. Wang et al. [8] experimentally studied the effects of operating conditions including the operating pressure, cell temperature, gas flow rate and humidification on the performance of PEM fuel cells with interdigitated flow fields. A similar work was carried out by Yan et al. [9] to investigate the influences of operating conditions on performance of PEM fuel cells with conventional and interdigitated flow fields. Huang et al. [10] reported influence of a porosity gradient in GDLs on the performance of PEM fuel cells with interdigitated flow fields by using a threedimensional, non-isothermal and two-phase model. It was concluded that the cell performance improvement was not significant, but the pressure drop was remarkably decreased. Recently, Cooper et al. [11] performed an experimental study on the performance improvement of PEM fuel cells with interdigitated flow fields as the ratio of channel length to width was decreased. In order to avoid gas leakage and minimize the contact resistance between different layers, an assembly force is applied on fuel cell stacks. With compression, the deformation of GDL is observed due to its porous feature, and the geometry and transport parameters are changed. Accordingly, the transport processes and cell performance are affected. Hottinen et al. [12] compared effect of homogeneous and inhomogeneous compressions of GDLs on the local transport characteristics and overall performance of PEM fuel cells. In that study, experimentally determined transport properties of GDLs were adopted in a two-dimensional, nonisothermal and single-phase model. Su et al. [13] performed numerical simulations on PEM fuel cells with compressed GDLs by using a three-dimensional, isothermal and singlephase model. The variations of porosity and permeability of GDLs caused by the assembly force were taken into account. A three-dimensional, isothermal and single-phase model was applied to investigate the performance of PEM fuel cells as deformation of GDLs was presented [14]. The changes of contact resistance between the GDLs and current collectors for different assembly forces were considered in the mathematical model, as well as the porosity of GDLs. Zhou et al. [15] reported the effect of assembly force on liquid water saturation in deformed GDLs and the average current density of PEM fuel cells for different types of GDLs and membranes. In that study, the finite element method was used to obtain the geometry of GDLs for different compressions, and a threedimensional, non-isothermal and two-phase model was applied to describe the complicated transport processes within the PEM fuel cells. The water transport and
performance of PEM fuel cells were experimentally investigated for different clamping forces by Cha et al. [16]. It was concluded that the ohmic resistance was decreased and the mass transport resistance was increased as the clamping forces increased. Wang et al. [17] numerically investigated the performance of PEM fuel cells with non-homogeneous deformation of the GDLs by using a three-dimensional, nonisothermal and single-phase model. The cell performance was significantly affected by the decrease in porosity caused by the GDL deformation. Toghyani et al. [18] studied the impact of clamping pressure on the thermal and electrochemical performance of PEM fuel cells with three different widths of channel. Shi et al. [19] studied the influence of compression on the water management of PEM fuel cells with different GDLs by using a two-dimensional, isothermal and two-phase model. Recently, the performance of PEM fuel cells with interdigitated flow fields and inhomogeneous compression of GDL has been numerically analyzed by using a twodimensional, isothermal and two-phase model [20]. The effect of GDL deformation on transport characteristics and cell performance in PEM fuel cells with parallel channels has already been widely and extensively investigated during the past decade [12e18]. However, studies on the effect of GDL deformation on the performance of PEM fuel cells with interdigitated flow fields are still very few [19,20]. A comprehensive study is needed to improve the understanding of transport characteristics of PEM fuel cells with interdigitated flow fields by using a more complete mathematical model. In the present study, a three-dimensional, non-isothermal, twophase flow model based on the finite volume method has been developed and applied to investigate the effect of GDL deformation on transport phenomena of PEM fuel cells with interdigitated flow fields. In addition, the influences of deformation on geometry and transport parameters are taken into account. The oxygen concentration, temperature, liquid water saturation, volumetric current density, and cell performance of PEM fuels with/without compression are presented and compared.
Model description Schematic illustrations of PEM fuel cells using interdigitated flow fields with/without compression are provided in Fig. 1. It is clear that the deformed GDLs intrude into the flow channel, and the cross-sectional area of the flow channel is decreased. Fig. 2 illustrates the three-dimensional computational domain considered in this study. Due to the symmetry of configurations in the y-coordinate direction only, half inlet and outlet channels are included in the computational domain. The detailed information of the fuel cell geometry and operating condition parameters are summarized in Table 1 and Table 2. The fuel cell operating temperature and pressure are 353 K and 1 atm, respectively. In addition, 100% relative humidity is applied for both the anode and cathode reactant gases. In the PEM fuel cell mathematical model, the fluid flow is laminar; ideal gas law is applied for the reactant gases; the CLs are homogeneous and isotropic; the reactant gases cannot diffuse across the membrane; the generated water in the cathode CL is in dissolved phase [23].
n B, Effects of gas diffusion layer deformation on the transport phenomena and perforPlease cite this article in press as: Li S, Sunde mance of PEM fuel cells with interdigitated flow fields, International Journal of Hydrogen Energy (2018), https://doi.org/10.1016/ j.ijhydene.2018.07.064
3
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 x x x ( 2 0 1 8 ) 1 e1 4
Fig. 1 e Schematic illustrations of PEM fuel cells with interdigitated flow fields: (a) without compression; (b) with compression.
Fig. 2 e Schematic illustrations of the computational domain of a PEM fuel cell with interdigitated flow fields: (a) threedimensional computational domain; (b) representative plane and line.
Table 1 e Fuel cell geometric dimensions [15,21,22]. Parameter Fuel cell length Fuel cell width Flow channel width Flow channel height Rib width Current collector height GDL thickness CL thickness Membrane thickness
Table 2 e Fuel cell operating conditions.
Value
Units
50 2 1 1 1 1.5 0.2 0.01 0.05
mm mm mm mm mm mm mm mm mm
Parameter Reactant gas Operating pressure, Pa/Pc Stoichiometric ratio, xa/xc Relative humidity, RHa/RHc Operating temperature, Ta/Tc
Anode
Cathode
Hydrogen/Water 1 atm 1.5 100% 353 K
Air/Water 1 atm 1.5 100% 353 K
n B, Effects of gas diffusion layer deformation on the transport phenomena and perforPlease cite this article in press as: Li S, Sunde mance of PEM fuel cells with interdigitated flow fields, International Journal of Hydrogen Energy (2018), https://doi.org/10.1016/ j.ijhydene.2018.07.064
4
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 x x x ( 2 0 1 8 ) 1 e1 4
With the above mentioned assumptions, the governing equations of mass, momentum, species, energy, charge, liquid water and dissolved water are presented as follows [15,23]: Mass conservation equation:
where Yi and Deff,i are the mass fraction and effective diffusivity for the i th species, respectively. Si represents the amount of consumption or generation of the i th species. The effective diffusivity (Deff,i) is obtained by the Bruggeman correlation (see Table 3) which accounts for the tortuous path in the GDLs and CLs. Energy conservation equation:
V$ðr! u Þ ¼ Smass
u T ¼ V$ keff VT þ ST V$ rCp !
Governing equations
(1)
where r and ! u are the mixture fluid density and superficial velocity, respectively. Smass represents the source term of the mass conservation equation. The mass changes of hydrogen, oxygen and water vapor due to the processes of electrochemical reactions and phase change are completely considered in the corresponding layers. Momentum conservation equation: V$ðr! u! u Þ ¼ V$ðmV! u Þ VP þ Smom
(4)
where Cp and keff are the specific heat and effective thermal conductivity, respectively. The irreversible, reversible, ohmic heat generation and phase change terms are all included in the energy equation source term, ST. Charge conservation equation: V$ seff ;s Vfs þ Ss ¼ 0
(5)
V$ seff ;m Vfm þ Sm ¼ 0
(6)
(2)
where P and m are the mixture pressure and dynamic viscosity, respectively. Smom denotes the source term of the momentum equation. Species conservation equation: V$ðr! u Yi Þ ¼ V$ rDeff ;i VYi þ Si
(3)
where seff,s is the effective electrical conductivity, seff,m the effective protonic conductivity, fs the electrical potential, fm the protonic potential. The Butler-Volmer equation and spherical agglomerate model are adopted to describe the
Table 3 e Parameters used in the mathematical model. Parameter
Value
mg cm2 kg m3 mg cm2 kg m3 kg m3 kg mol1 mm e e A m2
10(0.03741*T16.96)
A m2
ref
0.5 1 56.4
e e mol m3
[27]
3.39
mol m3
4.56 103 0.101325e(666/Tþ14.1) 100/1.7/0.3 0.25 20000/5000/2000 0.104 326.36 2.36 106 3.517 104 0.0625 110 /95 100 100 1 1.0e-13 9.15 105 2.82 105 2.2 105 2.56 105
Pa m3 mol1 Pa m3 mol1 W m1 K1 W m1 K1 S m1 J mol1 K1 J mol1 K1 J kg1 Pa s N m1 e s1 s1 s1 m2 m2 s1 m2 s1 m2 s1 m2 s1
ref
Anode reference exchange current density, ia [22] ref
Cathode reference exchange current density, ic [25] Anode transfer coefficient, aa [26] Cathode transfer coefficient, ac [26] Reference hydrogen concentration, cH2 [27] Reference oxygen concentration,
ref c O2
Units
0.4 2.145 104 0.6 1.8 103 1.98 103 1.1 1 0.4 0.5 100
Platinum loading, mpt [24] Platinum density, rpt [22] Carbon loading, mc [22] Carbon density, rc [22] Dry membrane density, rm [22] Membrane equivalent weight, Mm [22] Radius of agglomerate, ragg [22] Volume fraction of ionomer, Li [22] Volume fraction of ionomer in agglomerate, Li,agg [22]
Hydrogen Henry's constant, HH2 [23] Oxygen Henry's constant, HO2 [23] Thermal conductivity of CC/GDL/CL, kCC/GDL/CL [22] Thermal conductivity of membrane, km [22] Electrical conductivity of CC/GDL/CL, ss,CC/GDL/CL [22] Entropy of hydrogen oxidation, DSa [28] Entropy of oxygen reduction, DSc [28] Latent heat of condensation/evaporation, Dhlg [29] Liquid water viscosity, ml [29] Surface tension, s [29] Contact angle of GDL/CL, qGDL/CL [29] Condensation rate, gcon [23] Evaporation rate, gevap [23] Dissolved water phase change rate, g [23] Permeability of CL, KCL [22] Binary diffusivity, DH2 H2 O [30] Binary diffusivity, DO2 H2 O [30] Binary diffusivity, DO2 N2 [30] Binary diffusivity, DH2 ON2 [30]
n B, Effects of gas diffusion layer deformation on the transport phenomena and perforPlease cite this article in press as: Li S, Sunde mance of PEM fuel cells with interdigitated flow fields, International Journal of Hydrogen Energy (2018), https://doi.org/10.1016/ j.ijhydene.2018.07.064
5
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 x x x ( 2 0 1 8 ) 1 e1 4
hydrogen oxidation reaction and oxygen reduction reaction in the anode and cathode CLs, respectively. Liquid water transport equation: Krl mg ! u ¼ V$ðrl Ds VsÞ þ Sl V$ rl Krg ml
(7)
where Krl and Krg are the relative permeability of liquid and gas phase, respectively. Ds denotes the capillary diffusion coefficient. The source term, Sl, represents the phase change process between water vapor and liquid water, which is determined by the water vapor partial pressure and saturation pressure. Dissolved water transport equation:
V$
n r d sm Vfm ¼ V$ m Dl Vl þ Sd F Mm
(8)
where nd and Dl are the electro-osmotic drag coefficient and water diffusivity in the membrane, respectively. The generated water in the cathode CL is in dissolved phase, and the absorption/desorption processes of ionomer in the anode and cathode CLs take place between water vapor and dissolved water. The parameters and complementary equations used in the mathematical model are summarized in Tables 3 and 4. In addition, the source terms of the governing equations are given in Table 5.
The catalyst layer composition and volume fraction The catalyst layer consists of void space, Pt/C particles, and ionomer phase. The volume fraction of each component is determined by the following expressions [26,36,37].
Table 4 e Complementary equations and definitions. Description
Units 1:5 1:5
m s1 2
Effective mass diffusivity [31]
Deff ;i ¼ ð1 sÞ
Mass diffusivity [17]
1 Xi j¼1;jsi Xj =Di;j 1:5 P0 T Di;j ¼ Di;j ðT0 ; P0 Þ T0 P keff ¼ ð1 εÞks þ εkf !0:5 PH ref ja ¼ ð1 sÞia aeff ref 2 ½eaa Fha =RT eac Fha =RT c H 2 HH 2
Binary mass diffusivity [30] Effective thermal conductivity [21] Electrochemical kinetics [31]
Over-potential [31] Open circuit voltage [22] Active surface area [26] Proton conductivity [32] Effective conductivity [22]
ε
Di;m
m2 s1
Di;m ¼ Pn
m2 s1 W m1 K1 A m3
ha ¼ fs fm ; hc ¼ fs fm Voc
V V
Voc ¼ 1:229 8:456 104 ðT 298:15Þ þ 4:31 105 TLnðPH2 P0:5 O2 Þ mpt 3 2 3 ð227:79f 158:57f 201:53f þ 159:5Þ 10 aeff ¼ tCL 1 1 Þ 1268ð 303 T sm ¼ ð0:514l 0:326Þe
m1 S m1 S m1
seff ;s ¼ ð1 εGDL Þ1:5 ss ; seff ;s ¼ ð1 εCL Li Þ1:5 ss seff ;m ¼
Li1:5 sm
Relative permeability [23]
Krl ¼ Ks3 ; Krg ¼ Kð1 sÞ3
m2
Capillary diffusivity [23]
Ks3 dPc Ds ¼ ml ds ε 0:5 ð1:417s 2:12s2 þ 1:263s3 Þ Pc ¼ scosðqÞ K log10 Psat ¼ 2:1794 þ 0:02953ðT 273:15Þ 9:1837 105 ðT 273:15Þ2 þ
m2 s1
Capillary pressure [23] Saturation pressure [23] Electro-osmotic drag coefficient [32] Dissolved water diffusivity [32]
Equilibrium water content [33] Water activity [23] Oxygen diffusivity in liquid water [34]
1:4454 107 ðT 273:15Þ3 l nd ¼ 2:5 22 8 2:05l 3:25 1 1 < 10 ½2416ð303 TÞ 6:65 1:25l Dl ¼ 10 e : 2:563 0:33l þ 0:0264l2 0:000671l3 1:41 þ 11:3a 18:8a2 þ 16:2a3 ða < 1Þ leq ¼ 10:1 þ 2:94ða 1Þ ða 1Þ Xwv P þ 2s a ¼ Psat DO2 ;w ¼ 7:4 1012
TðjMH2 O Þ0:5 mH2 O VO0:62
Oxygen diffusivity in ionomer [23] DO2 ;i ¼ 2:88 10 Permeability of the GDL [35] Inertial coefficient [35]
K ¼
10
1 T1 ½2933ð313 e
Pa Pa e
ð2≪l < 3Þ ð3≪l < 4Þ ð4 < lÞ
m2 s1
e e m2 s1 m2 s1 m2
d2f ε3 2
16kCK ð1 εÞ
b ¼ 2:88 106
1 ε1:5 K
m1
n B, Effects of gas diffusion layer deformation on the transport phenomena and perforPlease cite this article in press as: Li S, Sunde mance of PEM fuel cells with interdigitated flow fields, International Journal of Hydrogen Energy (2018), https://doi.org/10.1016/ j.ijhydene.2018.07.064
6
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 x x x ( 2 0 1 8 ) 1 e1 4
Table 5 e Source terms in the governing equations [23]. Description
Units kg m3 s1
Smass ¼ SH2 þ Swv Anode CL Smass ¼ SO2 þ Swv Cathode CL Smass ¼ Swv Anode and cathode GDLs m ! Anode and cathode GDLs and CLs u þ brj! u j! u Smom ¼ K ja SH2 ¼ MH2 Anode CL 2F jc SO2 ¼ MO2 Cathode CL 4F Swv ¼ Sl Svd MH2 O Anode and cathode CLs Swv ¼ Sl Anode and cathode GDLs TDSa ST ¼ j a ha ja þ seff ;m kVfm k2 þ seff ;s kVfs k2 þ Sl Dhlg Anode CL 2F TDSc jc þ seff ;m kVfm k2 þ seff ;s kVfs k2 þ Sl Dhlg Cathode CL ST ¼ j c hc 4F ST ¼ seff ;m kVfm k2 membrane ST ¼ seff ;s kVfs k2 þ Sl Dhlg ST ¼ seff ;s kVfs k2
kg m2 s2 kg m3 s1 kg m3 s1 kg m3 s1 W m3
Anode and cathode GDLs
Anode and cathode CCs A m3
Ss ¼ ja Anode CL Ss ¼ þ jc Cathode CL Sm ¼ þ ja Anode CL Sm ¼ jc Cathode CL Sl ¼ Sphase Anode and cathode GDLs and CLs 8 εð1 sÞ > > MH2 O ðPwv Psat Þ Pwv Psat < gcond RT Sphase ¼ > εs > : gevap MH2 O ðPwv Psat Þ Pwv < Psat RT Sd ¼ Svd Anode CL Sd ¼ Svd þ Sl Cathode CL r Svd ¼ g m ðleq lÞ Anode and cathode CLs Mm jc Cathode CL Sl ¼ 2F
A m3 kg m3 s1 kg m3 s1
mol m3 s1
The volume fraction of the void space (εCL) is given by: εCL ¼ 1 LPt=C Li
(9)
where LPt/C and Li are the volume fractions of Pt/C particles and ionomer phase, respectively. The volume fraction of the Pt/C particles (LPt/C) is determined by: LPt=C ¼
mpt tCL
1 1f 1 þ rpt f rc
!
mpt mpt þ mc
(11)
where mc is the carbon loading. The volume fraction of the ionomer phase (Li) is expressed by: Li ¼
r3agg
h
LPt=C r3agg Li;agg þ 1 Li;agg
di ¼ ragg
ffi # "sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 3 Li 1 Li;agg Li;agg þ 1 1 LPt=C
(13)
The total volume of the liquid water is: (10)
where mpt is the platinum loading, tcl the thickness of the CL, rpt the density of platinum, rc the density of carbon. The parameter f is defined as the platinum loading divided by the sum of platinum loading and carbon loading. f¼
thickness di can be determined. It is assumed that the agglomerate is only occupied by Pt/C particles and ionomer phase. The thickness of ionomer film (di) is obtained by:
raggþdi
3
r3agg
i
Vw ¼ sεCL VCL
(14)
where s is the liquid water saturation, VCL is the total volume of the CL. The volume of the liquid water covering the individual agglomerate can be obtained: Vw;i ¼
sεCL N
(15)
where N is the number of agglomerate particles per catalyst layer volume. N¼
3LPt=C 4pr3agg 1 Li;agg
(16)
(12)
The ionomer phase consists of two parts: one part is inside the agglomerate, the other part is the ionomer film. When LPt/C, Li, ragg, and Li,agg are given, a unique value of the ionomer
The thickness of the liquid water film (dw) is obtained by: rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 3 3sεCL 3 ragg þ di þ dw ¼ ragg þ di 4pN
(17)
n B, Effects of gas diffusion layer deformation on the transport phenomena and perforPlease cite this article in press as: Li S, Sunde mance of PEM fuel cells with interdigitated flow fields, International Journal of Hydrogen Energy (2018), https://doi.org/10.1016/ j.ijhydene.2018.07.064
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 x x x ( 2 0 1 8 ) 1 e1 4
The cathode agglomerate model In the present analyses, the spherical agglomerate model [37,38] is used to calculate the volumetric current density. Fig. 3 depicts the structure of the agglomerates within the cathode CL. The identical agglomerates are uniformly distributed in the cathode CL, and the individual agglomerate is evenly covered by ionomer and liquid water films.
jc ¼ 4F
ragg þ di þ dw PO2 1 di dw þ þ HO2 Er kc ð1 εCL Þ ragg aagg;i DO2 ;i aagg;w DO2 ;w
1 (18)
where PO2 is the partial pressure of oxygen, HO2 the Henry's constant, Er the effective factor, kc the reaction rate constant, aagg,i the ionomer effective agglomerate surface area, and aagg,w the liquid water effective agglomerate surface area. The effectiveness factor of the spherical agglomerate is defined as: Er ¼
1 1 1 FL tanhð3FL Þ 3FL
(19)
where the Thiele's modulus FL for chemical reactions is given as: FL ¼
ragg 3
sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi kc L1:5 i;agg DO2
(20)
The reaction rate constant is computed as:
ref
kc ¼
ic aeff 4Fð1
ref εCL ÞcO2
eaa Fhc =RT þ eac Fhc =RT
(21)
The effective agglomerate surface area is obtained from: aagg;i ¼
2 3LPt=C εCL ragg þ di r3agg 1 Li;agg
aagg;w ¼
2 3LPt=C εCL ragg þ di þ dw r3agg 1 Li;agg
(22)
(23)
7
and non-uniformly distributed due to the deformation resulting from the assembly force. It is assumed that the thickness and porosity are decreased only caused by the change in pore size. Therefore, the GDL porosity is determined by the following expression [12,15]: ε ¼ 1 ð1 ε0 Þ
d0 d
(24)
where ε0 is the initial porosity, d0 the initial thickness, d the thickness after compression, ε the porosity after compression. The local porosity is varying with the thickness of the GDL, and the transport parameters associated with porosity are also changed.
Numerical implementation The developed mathematical model for PEM fuel cells is implemented in the commercial software ANSYS FLUENT. The governing equations of charge, liquid water transport, and dissolved water transport are implemented by using user defined scalar (UDS) equations with an under relaxation technique developed by the authors.
Boundary conditions At the inlet of the flow channels, the mass flow rate, temperature, and species mass fractions are prescribed. The mass flow rates of reactants at both anode and cathode sides are calculated at a reference current density of 1.0 A/cm2 and are given by Eqs. (25) and (26) The anode and cathode stoichiometric ratios (xa and xc) are provided in Table 2. In addition, the inlet liquid water saturation is assigned as zero. At the outlet of the flow channels, a pressure-outlet boundary condition is applied. The operating temperature and a constant electric potential, fs ¼ 0, are specified at the anode terminal. At the cathode terminal, the operating temperature and a constant electric potential, fs ¼ Vcell, are applied. A symmetry boundary condition is adopted at the outer surface of the x-z plane.
Determination of porosity of the GDLs after compression
Qa ¼
xa MH2 Iref Am 2FYH2
(25)
For simplicity, a constant and uniform porosity is adopted in the numerical investigation. However, the porosity is varied
Qc ¼
xc MO2 Iref Am 4FYO2
(26)
Fig. 3 e Schematic illustration of the ionomer and liquid water films covering the agglomerate.
n B, Effects of gas diffusion layer deformation on the transport phenomena and perforPlease cite this article in press as: Li S, Sunde mance of PEM fuel cells with interdigitated flow fields, International Journal of Hydrogen Energy (2018), https://doi.org/10.1016/ j.ijhydene.2018.07.064
8
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 x x x ( 2 0 1 8 ) 1 e1 4
The algorithm Semi-Implicit Method for Pressure Linked Equations (SIMPLE) is used to handle the pressure-velocity coupling. The second order upwind discretization scheme is adopted for momentum, species, energy and UDS equations. A grid independence study was carefully performed to balance the accuracy and computational resources. A mesh system with the total number of control volumes of 208,000 was selected in the following simulations. To assess the accuracy of the numerical model, the polarization curve obtained from the model is widely adopted to compare with the experimental data in the open literature. The numerical performance of the developed model are evaluated by comparison with the experimental data reported by Yan et al. [39], as shown in Fig. 4. It is found that the numerical method can reproduce the experimental data well. The information of geometry and operating parameters given in Ref. [39] is relatively complete, and therefore it is selected for model validation. For the experimental study, a 25 cm2 fuel cell with a platinum loading of 0.4 mg/cm2 for both the anode and cathode CLs was tested under the operating temperature of 353 K, operating pressure of 1 atm, and a relatively humidity of 100%. In addition, the present model has been applied to investigate the effects of agglomerate model parameters on local transport characteristics and overall cell performance [40].
because the Young's modulus of the GDLs is much smaller than that of the remaining components and the CLs are very thin. The thickness and porosity of the GDLs under different assembly forces are shown in Fig. 5. As shown in Fig. 5a, the thickness of the GDL under the rib is substantially decreased with increasing assembly forces, while the thickness of the GDL under the channel middle region is not affected by the assembly forces. In addition, the local porosity of the deformed GDLs which is determined by Eq. (24), is given in Fig. 5b. When the assembly force is increased from 0 MPa to 2.0 MPa, the thickness is decreased from 0.28 mm to 0.13 mm and the porosity is decreased from 0.78 to 0.53 at the rib region. Note that the variations of thickness and porosity of GDLs are considered in the y-coordinate direction and the porosity is a function of the y-coordinate. The cell performance is commonly presented and compared in terms of the polarization and power density curves. Fig. 6 shows the effect of assembly forces on the performance of PEM fuel cells with interdigitated flow fields.
Results and discussions In the preceding section, the PEM fuel cell mathematical model and GDL deformation were described in detail. In this section, effect of compression on the thickness and porosity of the GDLs is presented. Then the current densities and power densities of PEM fuel cells with/without compression are investigated and the local transport characteristics at an operating voltage 0.7 V are also analyzed. For simplicity, the configurations of the deformed GDLs under different assembly forces reported by Zhou et al. [15] are employed in the present study. It is assumed that the deformation after compression only takes place in the GDLs,
Fig. 4 e Comparison of the numerical results with the experimental data [39].
Fig. 5 e Distribution of the GDL characteristics under different assembly forces: (a) thickness; (b) porosity.
Fig. 6 e Polarization and power density curves of PEM fuel cells under different assembly forces.
n B, Effects of gas diffusion layer deformation on the transport phenomena and perforPlease cite this article in press as: Li S, Sunde mance of PEM fuel cells with interdigitated flow fields, International Journal of Hydrogen Energy (2018), https://doi.org/10.1016/ j.ijhydene.2018.07.064
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 x x x ( 2 0 1 8 ) 1 e1 4
It can be seen that the current density is increased with increasing assembly force. For an operating voltage 0.7 V, the current densities for four different cases are 0.2383 A/cm2, 0.2422 A/cm2, 0.2565 A/cm2 and 0.2575 A/cm2, respectively. In addition, the maximum power densities are 0.4227 W/ cm2, 0.4263 W/cm2, 0.4464 W/cm2 and 0.4503 W/cm2, respectively. In a similar work carried out by Mahmoudi et al. [20], the cell performance was decreased as the compression ratio was increased. Note that a constant pressure drop between the inlet and outlet channels were adopted in that study. He et al. [6] numerically investigated the influence of pressure drop on the cell performance. It was concluded that the cell performance was substantially improved with increasing pressure drop. The pressure drops between the cathode inlet and outlet flow channels for the four cases at the operating voltage 0.7 V are 146.2 Pa, 428.4 Pa, 1111.3Pa and 3104.1 Pa, respectively. The pressure drop of the PEM fuel cell with 2.0 MPa assembly force is approximately 21.2 times higher than that of the PEM fuel
9
cell without an assembly force. This is because the mass flow rate is constant for all cases in this study, and the thickness, porosity and permeability of the GDLs under the rib regions and the cross sectional area of the flow channels are decreased due to the deformation of the GDLs caused by the assembly force. The distributions of oxygen concentration, temperature and liquid water saturation at the cathode GDL-CL interface plane and the volumetric current density at the middle plane of cathode CL of the PEM fuel cell without compression are given in Fig. 7. It is clear that the oxygen concentration is gradually decreased from the inlet region to the outlet region in Fig. 7a. As shown in Fig. 7b, the temperature is slightly increased along the flow direction, and the minimum temperature is observed under the rib region. In Fig. 7c, the liquid water saturation is gradually increased along the flow direction due to the consumption of oxygen and generation of water. In addition, it is observed that the volumetric current density is gradually decreased along the flow direction and the
Fig. 7 e Transport characteristics of PEM fuel cells without compression: (a) oxygen concentration distribution at the cathode GDL-CL interface; (b) temperature distribution at the cathode GDL-CL interface; (c) liquid water saturation distribution at the cathode GDL-CL interface (c) volumetric current density distribution at the middle plane of cathode CL. n B, Effects of gas diffusion layer deformation on the transport phenomena and perforPlease cite this article in press as: Li S, Sunde mance of PEM fuel cells with interdigitated flow fields, International Journal of Hydrogen Energy (2018), https://doi.org/10.1016/ j.ijhydene.2018.07.064
10
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 x x x ( 2 0 1 8 ) 1 e1 4
Fig. 8 e Oxygen concentration distribution at the middle plane of cathode GDL of PEM fuel cells under different assembly forces: (a) without compression; (b) 1.0 MPa; (c) 1.5 MPa; (d) 2.0 MPa.
Fig. 9 e Temperature distribution at the middle plane of cathode GDL of PEM fuel cells under different assembly forces: (a) without compression; (b) 1.0 MPa; (c) 1.5 MPa; (d) 2.0 MPa. n B, Effects of gas diffusion layer deformation on the transport phenomena and perforPlease cite this article in press as: Li S, Sunde mance of PEM fuel cells with interdigitated flow fields, International Journal of Hydrogen Energy (2018), https://doi.org/10.1016/ j.ijhydene.2018.07.064
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 x x x ( 2 0 1 8 ) 1 e1 4
maximum volumetric current density appears under the rib region and close to the inlet channel side. The oxygen concentration distributions in the middle plane of the GDLs of PEM fuel cells under different assembly
11
forces are shown in Fig. 8. It is observed that the oxygen concentration gradually decreases from the inlet region to the outlet region for all cases due to the consumption of oxygen by the oxygen reduction reaction. The maximum oxygen
Fig. 10 e Liquid water saturation distribution at the middle plane of cathode GDL of PEM fuel cells under different assembly forces: (a) without compression; (b) 1.0 MPa; (c) 1.5 MPa; (d) 2.0 MPa.
Fig. 11 e Volumetric current density of PEM fuel cells under different assembly forces. n B, Effects of gas diffusion layer deformation on the transport phenomena and perforPlease cite this article in press as: Li S, Sunde mance of PEM fuel cells with interdigitated flow fields, International Journal of Hydrogen Energy (2018), https://doi.org/10.1016/ j.ijhydene.2018.07.064
12
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 x x x ( 2 0 1 8 ) 1 e1 4
concentration appearing at the region close to the inlet flow channel is increased when the assembly force is gradually increased. This can be explained by the increase of velocity in the flow channel, which is caused by the decrease of cross sectional area of the flow channel with increasing assembly forces. It is also observed that the region of minimum oxygen concentration in the outlet flow channel is enlarged. The temperature distributions in the middle plane of the GDLs of PEM fuel cells under different assembly forces are presented in Fig. 9. It is clear that the minimum temperature and maximum temperature appear under the rib region and flow channel region, respectively. Due to the presence of rib, the generated heat in CLs and GDLs can be easily and effectively transferred through it. When the assembly force is applied, the thickness and porosity of the GDLs are both decreased and the effective thermal conductivity is consequently increased. Accordingly, the temperature under the rib region is gradually decreased with increasing assembly forces. The liquid water saturation distributions in the middle plane of the GDLs of PEM fuel cells under different assembly forces are shown in Fig. 10. It can be seen that the maximum liquid water saturation appears under the rib regions for all cases, and it is gradually increased with increasing assembly forces. This is because the porosity and permeability of the GDLs are both decreased under a higher assembly force. In addition, a lower local temperature is also observed at this region as shown in Fig. 9. The volumetric current density distributions at the middle plane of CLs (Line-1) are shown in Fig. 11. The minimum and maximum volumetric current densities appear at the outlet channel region and under the rib region for all cases, respectively. The volumetric current density is gradually increased when the assembly force is increased from 0 Mpa to 1.5 Mpa. When the assembly force is increased from 1.5 MPa to 2.0 MPa, the maximum volumetric current density is increased and the minimum volumetric current density is decreased. This indicates that a higher variation of the volumetric current density is obtained under the assembly force of 2.0 MPa.
Conclusions In this study, the effect of GDL deformation caused by the assembly force on transport characteristics and performance of PEM fuel cells with interdigitated flow fields were numerically investigated by using a three-dimensional, nonisothermal, two-phase flow mathematical model based on the finite volume method. The following conclusions were drawn: The thickness and porosity of GDLs are decreased due to the deformation. The corresponding transport parameters are altered and consequently the cell performance is also affected. The oxygen concentration is increased and decreased at the inlet channel and outlet channel regions, respectively. The temperature under the rib region is decreased due to the increased effective thermal conductivity. The liquid water saturation under the rib region is increased due to the decreased local porosity and temperature. In addition, the cell performance improvement is obtained and accompanied with a relatively high pressure drop penalty. It is
concluded that the assembly force should be as small as possible to decrease the parasitic losses with consideration of the gas sealing problem. The predicted results provided detailed transport characteristics of PEM fuel cells with interdigitated flow fields under compression, and may serve as guide lines to improve the cell performance and enable proper design. The present numerical model can be further developed with consideration of the anisotropic properties in GDLs and the contact resistance between different layers. In addition, different types of GDLs can be further studied and compared, as well as the fuel cells with/without micro-porous layer.
Acknowledgments The work was carried out at the Department of Energy Sciences, Lund University. The first author gratefully acknowledges the financial support from China Scholarship Council (CSC).
references
[1] Wang Y, Chen KS, Mishler J, Cho SC, Adroher XC. A review of polymer electrolyte membrane fuel cells: technology, applications, and needs on fundamental research. Appl Energy 2011;88:981e1007. [2] Li XG, Sabir I. Review of bipolar plates in PEM fuel cells: flow field design. Int J Hydrogen Energy 2005;30:359e71. [3] Manso AP, Marzo FF, Barranco J, Garikano X, Mujika MG. Influence of geometric parameters of the flow fields on the performance of a PEM fuel cell. Int J Hydrogen Energy 2012;37:15256e87. [4] Heidary H, Kermani MJ, Advani SG, Prasad AK. Experimental investigation of in-line and staggered blockages in parallel flow field channels of PEM fuel cells. Int J Hydrogen Energy 2016;41:6885e93. [5] Chowdhury MZ, Akansu YE. Novel convergent-divergent serpentine flow fields effect on PEM fuel cell performance. Int J Hydrogen Energy 2017;42:25686e94. [6] He WS, Yi JS, Nguyen TV. Two-phase flow model of the cathode of PEM fuel cells using interdigitated flow fields. AIChE J 2000;46:2053e64. [7] Um S, Wang CY. Three dimensional analysis of transport and electrochemical reactions in polymer electrolyte fuel cells. J Power Sources 2004;125:40e51. [8] Wang L, Liu HT. Performance studies of PEM fuel cells with interdigitated flow fields. J Power Sources 2004;134:185e96. [9] Yan WM, Chen CY, Mei SC, Soong CY, Chen FL. Effects of operating conditions on cell performance of PEM fuel cells with conventional or interdigitated flow field. J Power Sources 2006;162:1157e64. [10] Huang YX, Cheng CH, Wang XD, Jang JY. Effects of porosity gradient in gas diffusion layers on performance of proton exchange membrane fuel cells. Energy 2010;35:4786e94. [11] Cooper NJ, Santamaria AD, Becton MK, Park JW. Investigation of the performance improvement in decreasing aspect ratio interdigitated flow field PEMFCs. Energy Convers Manag 2017;136:307e17. [12] Nitta I, Karvonen S, Himanen O, Mikkola M. Modelling the effect of inhomogeneous compression of GDL on local transport phenomena in a PEM fuel cell. Fuel Cell 2008;8:410e21.
n B, Effects of gas diffusion layer deformation on the transport phenomena and perforPlease cite this article in press as: Li S, Sunde mance of PEM fuel cells with interdigitated flow fields, International Journal of Hydrogen Energy (2018), https://doi.org/10.1016/ j.ijhydene.2018.07.064
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 x x x ( 2 0 1 8 ) 1 e1 4
[13] Su ZY, Liu CT, Chang HP, Li CH, Huang KJ, Sui PC. A numerical investigation of the effects of compression force on PEM fuel cell performance. J Power Sources 2008;183:182e92. [14] Taymaz I, Benli M. Numerical study of assembly pressure effect on the performance of proton exchange membrane fuel cell. Energy 2010;35:2134e40. [15] Zhou YB, Jiao K, Du Q, Yin Y, Li XG. Gas diffusion layer deformation and its effect on the transport characteristics and performance of proton exchange membrane fuel cell. Int J Hydrogen Energy 2013;39:12891e903. [16] Cha D, Ahn JH, Kim HS, Kim Y. Effects of clamping force on the water transport and performance of a PEM fuel cell with relative humidity and current density. Energy 2015;93:1338e44. [17] Wang JT, Yuan JL, Yu JS, Sunden B. Investigation of effects of non-homogeneous deformation of gas diffusion layer in a PEM fuel cell. Int J Energy Res 2017;41:2121e37. [18] Toghyani S, Moradi Nafchi F, Afshari E, Hasanpour K, Baniasadi E, Atyabi SA. Thermal and electrochemical performance analysis of a proton exchange membrane fuel cell under assembly pressure on gas diffusion layer. Int J Hydrogen Energy 2018;43:4534e45. [19] Shi ZY, Wang X, Guessous L. Effect of compression on the water management of a proton exchange membrane fuel cell with different gas diffusion layers. J Fuel Cell Sci Technol 2010;7:021012. [20] Mahmoudi AH, Ramiar A, Esmaili Q. Effect of inhomogeneous compression of gas diffusion layer on the performance of PEMFC with interdigitated flow field. Energy Convers Manag 2016;110:78e89. [21] Li SA, Yuan JL, Andersson M, Xie GN, Sunden B. Wavy surface cathode gas flow channel effects on transport processes in a proton exchange membrane fuel cell. ASME J Electrochem En Conv Stor 2017;14:031007. [22] Li SA, Yuan JL, Andersson M, Xie GN, Sunden B. Influence of anisotropic gas diffusion layers on transport phenomena in a proton exchange membrane fuel cell. Int J Energy Res 2017;41:2034e50. [23] Jiao K, Li XG. Water transport in polymer electrolyte membrane fuel cells. Prog Energy Combust Sci 2011;37:221e91. [24] Wang L, Husar A, Zhou TH, Liu HT. A parametric study of PEM fuel cell performance. Int J Hydrogen Energy 2003;28:1263e72. [25] Parthasarathy A, Srinivasan S, Appleby AJ, Martin CR. Temperature dependence of the electrode kinetics of oxygen reduction at the Platinum/Nafion interface-a microelectrode investigation. J Electrochem Soc 1992;139:2530e7. [26] Dalasm NKH, Kermani MJ, Moghaddam DG, Stockie JM. A parametric study of cathode catalyst layer structural parameters on the performance of a PEM fuel cell. Int J Hydrogen Energy 2010;35:2417e27. [27] Bernardi DM, Verbrugge MW. A mathematical model of the solid polymer electrolyte fuel cell. J Electrochem Soc 1992;139:2477e91. [28] Lampinen MJ, Fomino M. Analysis of free energy and entropy changes for half cell reactions. J Electrochem Soc 1993;140(12):3537e46. [29] Yang XG, Ye Q, Cheng P. Matching of water and temperature fields in proton exchange membrane fuel cells with nonuniform distributions. Int J Hydrogen Energy 2011;36:12524e37. [30] Berning T, Lu DM, Djilali N. Three-dimensional computational analysis of transport phenomena in a PEM fuel cell. J Power Sources 2002;106:284e94.
13
[31] Meng H. Numerical studies of liquid water behaviors in PEM fuel cell cathode considering transport across different porous media. Int J Hydrogen Energy 2010;35:5569e79. [32] Springer TE, Zawodzinski TA, Gottesfeld S. Polymer electrolyte fuel cell model. J Electrochem Soc 1991;138(8):2334e42. [33] Zawodzinski TA, Springer TE, Uribe F, Gottesfeld S. Characterization of polymer electrolytes for fuel cell applications. Solid State Ionics 1993;60:199e211. [34] Wilke AJ, Chang P. Correlation of diffusion coefficients in dilute solutions. AIChE J 1995;1:264e70. [35] Gostick JT, Fowler MW, Pritzker MD, Ioannidis MA, Behra LM. In-plane and through-plane gas permeability of carbon fiber electrode backing layers. J Power Sources 2006;162:228e38. [36] Secanell M, Karan K, Suleman A, Djilali N. Multi-variable optimization of PEMFC cathodes using an agglomerate model. Electrochim Acta 2007;52:6318e37. [37] Shah AA, Kim GS, Sui PC, Harvey D. Transient nonisothermal model of a polymer electrolyte fuel cell. J Power Sources 2007;163:793e806. [38] Sun W, Peppley B, Karan K. An improved two-dimensional agglomerate cathode model to study the influence of catalyst layer structural parameters. Electrochim Acta 2005;50:3359e74. [39] Yan QG, Toghiani H, Causey H. Steady state and dynamic performance of proton exchange membrane fuel cells under various operating conditions and load changes. J Power Sources 2006;161:492e502. [40] Li SA, Yuan JL, Xie GN, Sunden B. Effects of agglomerate model parameters on transport characterization and performance of PEM fuel cells. Int J Hydrogen Energy 2018;43:8451e63.
Nomenclature a: Effective surface area, m1/Water activity A: Area, m2 c: Mole concentration, mol m3 CP: Specific heat, J kg1 K1 D: Diffusivity, m2 s1 E: Effectiveness factor f: Platinum mass ratio to Pt/C F: Faraday constant, 96485 C mol1 h: Enthalpy change, J kg1 H: Henry's constant, Pa m3 mol1 i: Exchange current density, A m2 j: Volumetric current density, A m3 k: Thermal conductivity, W m1 K1/Reaction rate constant, s1 K: Permeability, m2 L: Volume fraction m: Mass loading, kg m2 M: Molecular weight, kg mol1 n: Osmotic-drag coefficient P: Pressure, Pa Q: Mass flow rate, kg s1 r: Agglomerate radius, m R: Universal gas constant, 8.314 J mol1 K1 s: Liquid water saturation S: Source term/Entropy, J mol1 K1 t: Thickness, m T: Temperature, K ! u : Velocity vector, m/s V: Voltage, V/Volume, m3 X: Mole fraction Y: Mass fraction
n B, Effects of gas diffusion layer deformation on the transport phenomena and perforPlease cite this article in press as: Li S, Sunde mance of PEM fuel cells with interdigitated flow fields, International Journal of Hydrogen Energy (2018), https://doi.org/10.1016/ j.ijhydene.2018.07.064
14
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 x x x ( 2 0 1 8 ) 1 e1 4
Greek Symbols a: Transfer coefficient b: inertial coefficient, m1 g: Water phase change rate, s1 d: Thickness, m ε: Porosity h: Over-potential, V q: Contact angle l: Water content m: Dynamic viscosity, Pa s x: Stoichiometric ratio r: Density, kg m3 s: Electron/Proton conductivity, S m1/Surface tension, N m2 f: Potential, V F: Theile's modulus Subscripts and Superscripts a: Anode agg: Agglomerate c: Cathode/Carbon con: Condensation d: Osmotic-drag/Dissolved water eff: Effective eq: Equilibrium
evap: Evaporation g: Gas i: i th species/Ionomer j: j th species l: Liquid lg: Liquid and gas water m: Membrane/Mixture mom: Momentum equation oc: Open circuit ref: Reference rl: Liquid relative rg: Gas relative phase: Phase change Pt: Platinum sat: Saturation s: Solid/Liquid water saturation vd: Vapor and dissolved water w: Water wv: Water vapor Abbreviations CC: Current collector CL: Catalyst layer GDL: Gas diffusion layer PEM: Proton exchange membrane
n B, Effects of gas diffusion layer deformation on the transport phenomena and perforPlease cite this article in press as: Li S, Sunde mance of PEM fuel cells with interdigitated flow fields, International Journal of Hydrogen Energy (2018), https://doi.org/10.1016/ j.ijhydene.2018.07.064