Effects of Cathode Gas Diffusion Layer Configuration on the Performance of Open Cathode Air-Cooled Polymer Electrolyte Membrane Fuel Cell

Abstract

The design of a gas diffusion layer (GDL) is an effective way to manage water transport, thus improving the performance of air-cooled fuel cells. In the present study, three group designs of GDL with polytetrafluoroethylene (PTFE)—uniformly doped, in-planed sandwich doped and through-plane gradient doped—are proposed, and their effects on the performance of air-cooled fuel cells are explored by numerical simulation. The distribution of key physical quantities in the cathode catalyst layer (CCL), current density and the uniformity of current density distribution in the CCL were analyzed in detail. The results show that properly reducing the amount of PTFE in GDL is beneficial to promoting the water retaining capacity of air-cooled fuel cells, and then improving the performance of fuel cells. The performance of the in-plane sandwich GDL design cannot exceed the design with 10% PTFE uniformly doped, and this design will aggravate the uneven distribution of current density in CCL. Compared with the design of GDL with 40% PTFE uniformly doped, the current density can be improved by 22% when operating at 0.6 V by gradient increasing the PTFE content in GDL from the GDL/MPL interface to the gas channel. Furthermore, this design can maintain as good a current density uniformity as uniformly doping schemes.

1. Introduction

In order to overcome the problems of low energy density and long charge time of lithium-ion battery, the polymer electrolyte membrane fuel cell (PEMFC) has been considered an alternative power source for portable and mobile applications such as Unmanned Aerial Vehicles (UAVs). PEMFC is an efficient energy conversion device that converts chemical energy into electrical energy. The reactant is oxygen and hydrogen, while the by-product is only water and waste heat. Due to the merit of high energy density, the endurance of UAVs with PEMFC as the power source can be more than two hours, far beyond that of lithium battery. When applied to UAVs, the miniaturization and compact design of PEMFCs are necessary [1]. Air-cooled fuel cell provides the possibility of miniaturization by removing the bulky balance of plant (BoP) such as air compressors, hydrogen circulation systems, water cooling systems and humidifiers. Therefore, the weight, system complexity and the parasitic power of the fuel cell system are effectively reduced.

During the operation of PEMFCs, more than 50% of the chemical energy is converted into waste heat [2]. Excessive heat will lead to the membrane dehydration thus decreasing the proton conductivity. To tackle the overheating of fuel cell stack, it is necessary to provide excess air to cool the fuel cell and maintain the proper operating temperature [3]. However, high speed dry air will accelerate water evaporation, resulting in membrane dehydration as well. Without the external humidifying, air-cooled fuel cells usually confront the severe situation of membrane dehydration and the performance is vulnerable to the operating state and the change of ambient conditions such as low relative humidity and high temperature [4]. From the above analysis, it can be found that water management significantly affects the performance of air-cooled fuel cells. Maintaining a good hydration state of the membrane is required to increase the proton conductivity [5,6]. To date, many studies have been conducted to improve the performance of air-cooled fuel cells by balancing the water and thermal management.

As one of the key components of PEMFC, gas diffusion layer (GDL) plays an important role in mass and water transport [7,8]. GDL is formed by carbon fibers with diameter of 6–8 μm. The pores between carbon fibers are the bridge between the internal reaction area and the external environment with the pore character diameter of 10–100 μm [9]. Oxygen, liquid water and water vapor transport through pores, while the solid skeleton of carbon fiber constitutes the path for electron conduction and heat transfer. The structural and physical characteristics of GDL, such as porosity distribution, surface contact angle, permeability and diffusivity, have a significant impact on the performance and lifetime of fuel cells [10]. The change of wettability is mainly achieved by adjusting the content of polytetrafluoroethylene (PTFE) in GDL [8], while the effective diffusion coefficient and effective thermal conductivity are mainly controlled by adjusting the porosity of GDL and fiber arrangement. When doped with PTFE, the wettability and porosity of GDL will both be changed. With the increase of PTFE content in the GDL, the material will become more hydrophobic, and the porosity of the GDL will decrease at the same time. In order to keep the water balance and improve the mass transport capacity, the physical characteristics and structural configuration of GDL have been widely studied and optimized. Generally, water-cooled fuel cells operate at high power. It produces more water, and the cathode stoichiometry is less than that of air-cooled fuel cells. If water cannot be expelled quickly, it will accumulate in the gas channel (GC) and the MEA and block the oxygen path to reaction site thus causing water flooding issues [11]. The design of GDL of water-cooled fuel cell usually focuses on alleviating water flooding problems and improving mass transfer. Using the two-phase numerical method, hydrophobic treating GDL surface was found to be beneficial for the water removal performance [12,13]. Chen et al. [14,15] found that the change of GDL structure will significantly affect the transport behavior of liquid water, thus affecting the water management characteristics of fuel cells. Mortazavi et al. [16] experimentally studied the effects of PTFE content on the liquid water dynamic behavior under different gas velocity. Increasing the PTFE content is beneficial to reducing the detachment diameter of droplets. Yu et al. [17] compared the water expelling characteristics of GDL with various through-plane distributing strategies of PTFE. The results indicated that more PTFE distributes near catalyst layer (CL) side in the GDL will reduce the breakthrough time of water and accelerate the process of water removal. Zhou et al. [18] comprehensively investigated the effect of wetting patterns of GDL on the two-phase flow behavior. The compressing deformation of GDL caused by clamping pressure was also considered in the study. As for the through-plane PTFE distribution, the results showed the same trend as Yu et al. [17], namely more PTFE treatment GDL near the CL side is desirable for the water removal. Wang et al. [19] numerically studied the effect of PTFE doping strategies on the water distribution and transport behavior. GDL with higher PTFE content showed a fast water removal performance with the cost of increase of mass transfer resistance. A hybrid wettability GDL along GC was proposed by considering the tradeoff between species transport and water removal and the superiority of the design was proved by numerical simulation. To deepen the understanding of the influence of GDL on flow reaction process, non-uniform pore and wettability distribution of GDL were studied in detail [20]. Then in order to regulate the local saturation of liquid water at the interface between GDL and microporous layer (MPL) to improve the performance of fuel cells, a gradient GDL pore structure was proposed.

However, for air-cooled fuel cells, the main contradiction has shifted from avoiding water flooding to preventing membrane dehydration as discussed above. Through the X-ray tomography and neutron imaging technology, Meyer et al. [21] discovered the water transport mechanism in MEA of air-cooled fuel cell. It was found that the PTFE content and porosity of GDL affect the water transport through the membrane significantly. Pei et al. [22] found that proper PTFE content and thickness of GDL can improve the air-cooled fuel cell stack operation performance. There is a tradeoff between the gas diffusion, contact resistance and water removal ability. Wu et al. [23] experimentally studied effects of GDL wettability on the performance of open cathode fuel cell. It was found that despite the introduction of high contact resistance and gas mass transfer resistance, the fuel cell with hydrophilic diffusion layer can still ensure high performance and stable output at high temperature even at 60 °C. Park [24] numerically studied effects of GDL wettability on the fuel cell performance by modifying the polytetrafluoroethylene (PTFE) content. The results showed that without external humidifying there is an optimum water saturation range between 0.1–0.3 to ensure good hydration of membrane and maintain proper diffusivity at the same time. Zhao et al. [25] studied the effect of PTFE content and GDL thickness on air-cooled fuel cell performance The results suggested that an optimum values of PTFE content for MPL and GDL are 40% and 10% respectively, and the optimum GDL thickness is 200 μm. To enhance the self-humidification effect without introducing additional humidifiers, Kong et al. [6,26] modified the structural design of GDL. Results showed that by placing a thin GDL layer with small porosity near the GC side can effectively improve the water retaining capacity thus promoting fuel cell performance under low humidity conditions.

To date, considerable studies on the GDL design mainly focus on effects of porosity and wettability caused by modifying the PTFE distribution on the water expelling performance for water-cooled fuel cells. From the above review, it can also be concluded that GDL design of the water-cooled fuel cell has mainly experienced the development of uniform regulation, non-uniform regulation and gradient regulation of porous media properties. However, the water management requirements of air-cooled fuel cells are different from the water-cooled fuel cells. Although there are several experimental studies performed to investigate the effect of GDL characteristics on the air-cooled fuel cell performance, further research about the effects of GDL characteristics on the internal parameter distribution of air-cooled fuel cells is still needed. The purpose of present study is to gain deep understanding of the effects of GDL design on the air-cooled fuel cell performance and to provide a reasonable GDL design strategy. The study is organized as follows: In Section 2, design concept of different PTFE doping strategies is proposed. Various configurations of GDL are designed by modifying the distribution of PTFE content in GDL and can be summarized into three categories, namely, the PTFE uniformly doping design, through-plane gradient doping design and the in-plane sandwich doping design. Numerical method is an effective tool to study the multi-coupling transport process in the air-cooled fuel cell. Therefore, a three-dimensional CFD numerical simulation method, the boundary conditions and grid system adopted in present study are described in Section 3. Then, in Section 4, coupling relationship between various variables and the mechanism of their influence on the performance of air-cooled fuel cells are discussed in detail. Effects of PTFE doping strategies on the distribution of key physical variables, the air-cooled fuel cell performance and the uniformity of current density distribution are comparatively studied. At last, the conclusion of the study is summarized in Section 5.

2. Design of GDL with Different PTFE Doping Strategies

Schematic of the air-cooled fuel cell and its main components are shown in Figure 1, including bipolar plates (BP), gas diffusion layers (GDL) usually with a thinner micro porous layer (MPL), catalyst layers (CL) and proton exchange membrane (PEM) [21,27]. When the metal bipolar plate is stamped, alternating channels are formed on two sides of the metal plate. The air in the cathode gas channel (CGC) can directly contact with the cathode GDL and supply oxygen needed by the electrochemical reaction. On the other hand, the channel isolated from GDL by bipolar plate is cooling channel (CC). The air in the CC is only used as the cooling medium. The configuration of the bipolar plate decouples the heat dissipation and reaction requirements, which can ensure heat dissipation and avoid excessive airflow drying the membrane at the same time [28].

Reasonable designs of GDL can effectively improve the water and heat management and avoid self-weight increase of fuel cell at the same time. In order to address the membrane dehydration of air-cooled fuel cell, several GDL design strategies are proposed. The design concept is shown in Figure 2. These designs of GDL can be summarized into three categories, namely, the uniform PTFE doping scheme as shown in Figure 2a, the through-plane gradient PTFE doping scheme as shown in Figure 2b and the in-plane sandwich PTFE doping scheme as shown in Figure 2c. To realize these designs, PTFE of different content is doped into the different region of GDL. For the through-plane gradient doping scheme as shown in Figure 2b, the PTFE is linearly doped along the thickness direction of GDL. With the increase of PTFE content, the hydrophobicity is increased and the porosity is reduced at corresponding position of GDL. As for the in-plane sandwich scheme as shown in Figure 2c, more PTFE is doped under the GC or under the CC to form two different sandwich designs, respectively. The effects of PTFE doping content on porosity and wettability can be quantitatively measured by mercury intrusion and sessile drop method. The relationship between the characteristics of GDL with different PTFE doping is summarized in

Table 1. Effects of PTFE content on the GDL characteristics.

PTFE Content (wt%)	Contact Angle	Porosity
10%	110°	0.73
20%	130°	0.69
30%	150°	0.65
40%	162°	0.63

[19].

3. Numerical Model

3.1. Fuel Cell Model

To describe the fluid flow, species transfer and heat transport phenomena in PEMFC, a three-dimensional two-phase fuel cell model [29,30] is utilized. The equations are established based on the following assumptions [29]: (1) fuel cell operates in steady state condition, (2) gas mixture obeys laminar flow and ideal gas law at low Reynolds number and pressure, (3) the properties such as permeability and porosity of the porous components of fuel cell are isotropic, (4) water generated by electrochemical reaction is in liquid phase due to the lower operating temperature, (5) negligible contact resistance between different components in PEMFCs. According to the above assumptions, the governing equations of the transport process in the PEMFC are listed in

Table 2. Fuel cell model.

Governing Equations	Expressions	Solution Zones
Mass	$\nabla \left({\rho }_{\mathrm{g}}{u}_{\mathrm{g}}\right)=S_{\mathrm{m}}$	CHs, GDLs, MPLs, CLs
Momentum	$\left(\frac{u_{\mathrm{g}}}{{\epsilon }_{\mathrm{eff}}}\nabla \right)\left({\rho }_{\mathrm{g}}{u}_{\mathrm{g}}\right)=-{\epsilon }_{\mathrm{eff}}\nabla p_{\mathrm{g}}+{\epsilon }_{\mathrm{eff}}{\mu }_{\mathrm{g}}{\nabla }^2\left(\frac{u_{\mathrm{g}}}{{\epsilon }_{\mathrm{eff}}}\right)+S_{\mathrm{u}}$	CHs, GDLs, MPLs, CLs
Species	$\nabla \left({\rho }_{\mathrm{g}}{u}_{\mathrm{g}}{X}_{\mathrm{k}}\right)=\nabla \left({\rho }_{\mathrm{g}}{D}_{\mathrm{k}}^{\mathrm{eff}}\nabla X_{\mathrm{k}}\right)+S_{\mathrm{k}}$	CHs, GDLs, MPLs, CLs
Electron	$0=\nabla \left({\kappa }_{\mathrm{ele}}^{\mathrm{eff}}\nabla {\varphi }_{\mathrm{ele}}\right)+S_{\mathrm{ele}}$	All zones
Proton	$0=\nabla \left({\kappa }_{\mathrm{ion}}^{\mathrm{eff}}\nabla {\varphi }_{\mathrm{ion}}\right)+S_{\mathrm{ion}}$	All zones
Dissolved water	$0=\frac{{\rho }_{\mathrm{mem}}}{EW}\nabla \left(D_{\mathrm{mw}}^{\mathrm{eff}}\nabla \lambda \right)+S_{\mathrm{mw}}$	CLs, PEM
Liquid water	$\nabla \left({\rho }_{\mathrm{l}}\frac{{\mu }_{\mathrm{g}}}{{\mu }_{\mathrm{l}}}\frac{K_{\mathrm{l}}}{K_{\mathrm{g}}}{u}_{\mathrm{g}}\right)=\nabla \left({\rho }_{\mathrm{l}}{D}_{\mathrm{c}}\nabla s\right)+S_{\mathrm{liq}}$	GDLs, MPLs, CLs
Energy	$\nabla \left({\left(\rho c_p\right)}_{\mathrm{f},\mathrm{s}}^{\mathrm{eff}}{u}_{g}T\right)=\nabla \left(k_{\mathrm{f},\mathrm{s}}^{\mathrm{eff}}\nabla T\right)+S_{\mathrm{T}}$	All zones.

.

The solution zones of these equations are also presented in Table 2. As for the mass, momentum conservation and species transport equations of gas mixture, the solution zones contain the channel zones (GC and CC) and the porous zones (GDLs, MPLs and CLs). The electron and proton fluxes are controlled by Ohm’s law and the driven forces are the gradient of the electronic and protonic potential respectively. The dissolved water is existed within the PEM and CLs and solved in these limited zones. The energy equation is solved in the entire domain. These equations are coupled with each other. The expressions of source terms [29,31] in governing equations are summarized in

Table 3. Source terms of the governing equation.

Source Terms	Expressions
Mass	$S_{\mathrm{m}}=S_{{\mathrm{H}}_2}+S_{{\mathrm{O}}_2}+S_{\mathrm{vp}}$
Momentum	$\mathrm{For}\mathrm{porous}\mathrm{media},{\mathrm{S}}_{\mathrm{u}}=-{\mu }_{\mathrm{g}}{u}_{\mathrm{g}}/K$
Species	$\mathrm{In}\mathrm{ACL},S_{{\mathrm{H}}_2}=-\frac{j_a}{2F}{\mathrm{M}}_{{\mathrm{H}}_2}$$;\mathrm{In}\mathrm{CCL},{\mathrm{S}}_{{\mathrm{O}}_2}=-\frac{j_c}{4F}{\mathrm{M}}_{\mathrm{O}}{}_{{}_2}$; $\mathrm{In}\mathrm{GDL}/\mathrm{MPL}/\mathrm{CL},{\mathrm{S}}_{\mathrm{vp}}=-{\mathrm{S}}_{\mathrm{vl}}$
Electronic potential	$\mathrm{In}\mathrm{ACL},{\mathrm{S}}_{\mathrm{ele}}=-j_a$$;\mathrm{In}\mathrm{CCL},{\mathrm{S}}_{\mathrm{ele}}=j_c$
Ionic potential	$\mathrm{In}\mathrm{ACL},{\mathrm{S}}_{\mathrm{ion}}=j_{\mathrm{a}}$$;\mathrm{In}\mathrm{CCL},{\mathrm{S}}_{\mathrm{ion}}=-j_c$
Liquid water	$\mathrm{In}\mathrm{ACL},{\mathrm{S}}_{\mathrm{liq}}={\mathrm{S}}_{\mathrm{vl}}-{\mathrm{M}}_{{\mathrm{H}}_2}{}_{\mathrm{O}}{\mathrm{S}}_{\mathrm{wd}}$; $\mathrm{In}\mathrm{CCL},{\mathrm{S}}_{\mathrm{liq}}=\frac{j_c}{4F}{\mathrm{M}}_{{\mathrm{H}}_2}{}_{\mathrm{O}}+{\mathrm{S}}_{\mathrm{vl}}-{\mathrm{M}}_{{\mathrm{H}}_2}{}_{\mathrm{O}}{\mathrm{S}}_{\mathrm{wd}}$ $\mathrm{In}\mathrm{GDL}/\mathrm{MPL},{\mathrm{S}}_{\mathrm{liq}}={\mathrm{S}}_{\mathrm{vl}}$
Dissolved water	$\mathrm{In}\mathrm{CLs},{\mathrm{S}}_{\mathrm{dis}}={\mathrm{S}}_{\mathrm{wd}}+{\mathrm{S}}_{\mathrm{EOD}}$
Energy	$\mathrm{In}\mathrm{BP},{\mathrm{S}}_{\mathrm{T}}={\kappa }_{\mathrm{ele}}^{\mathrm{eff}}{\Vert \nabla {\varphi }_{\mathrm{ele}}\Vert }^2$$;\mathrm{In}\mathrm{GDL},{\mathrm{S}}_{\mathrm{T}}={\kappa }_{\mathrm{ele}}^{\mathrm{eff}}{\Vert \nabla {\varphi }_{\mathrm{ele}}\Vert }^2+{\mathrm{S}}_{\mathrm{pc}}$; $\mathrm{In}\mathrm{ACL},{\mathrm{S}}_{\mathrm{T}}=\frac{j_{a}T\Delta {\mathrm{S}}_{\mathrm{a}}}{2F}+j_a\left|{\eta }_{\mathrm{act}}\right|+{\kappa }_{\mathrm{ele}}^{\mathrm{eff}}{\Vert \nabla {\varphi }_{\mathrm{ele}}\Vert }^2+{\kappa }_{\mathrm{ion}}^{\mathrm{eff}}{\Vert \nabla {\varphi }_{\mathrm{ion}}\Vert }^2+{\mathrm{S}}_{\mathrm{pc}}$; $\mathrm{In}\mathrm{CCL},{\mathrm{S}}_{\mathrm{T}}=\frac{j_{c}T\Delta {\mathrm{S}}_{\mathrm{c}}}{4F}+j_c\left|{\eta }_{\mathrm{act}}\right|+{\kappa }_{\mathrm{ele}}^{\mathrm{eff}}{\Vert \nabla {\varphi }_{\mathrm{ele}}\Vert }^2+{\kappa }_{\mathrm{ion}}^{\mathrm{eff}}{\Vert \nabla {\varphi }_{\mathrm{ion}}\Vert }^2+{\mathrm{S}}_{\mathrm{pc}}$; $\mathrm{In}\mathrm{membrane},{\mathrm{S}}_{\mathrm{T}}={\kappa }_{\mathrm{ion}}^{\mathrm{eff}}{\Vert \nabla {\varphi }_{\mathrm{ion}}\Vert }^2$
	$S_{\mathrm{pc}}=\{\begin{array}{ll}{h}_{\mathrm{con}}{S}_{\mathrm{vl}}& \hspace{1em}\hspace{1em}\mathrm{in}\mathrm{GDLs}/\mathrm{MPLs}\\ h_{\mathrm{con}}{S}_{\mathrm{vl}}-h_{\mathrm{hydr}}{S}_{\mathrm{wd}}{M}_{{\mathrm{H}}_2}{}_{\mathrm{O}}& \hspace{1em}\hspace{1em}\mathrm{in}\mathrm{CLs}\end{array}$

. The parameters related to the porous media, the transport and physical properties are presented in

Table 4. Transport parameters and physical properties.

Description	Value
Porosity	εCLs/MPLs/GDLs: 0.4/0.4/0.78
Permeability	KCLs/MPLs/GDLs: 1 × 10−13/7 × 10−13/2 × 10−12 m2
Contact angle	θCLs/MPLs: 120°/150°
Electronic conductivity,	κele CLs/MPLs/GDLs/BPs: 1000/5000/17500/20000 S m−1
Equivalent weight of membrane	EW: 1.0 kg mol−1
Density of membrane	ρEW: 1980 kg m−3
Heat conductivity	ks MEM/CLs/MPLs/GDLs/BPs: 0.95/1/0.83/1/20 W m−2 K−1
Heat capacity	cp MEM/CLs/MPLs/GDLs/BPs: 833/3300/800/568/1580 J kg−1 K−1
Oxygen diffusivity and hydrogen diffusivity	DO2/DH2: 3.732 × 10−5/5.717 × 10−5 m2 s−1
Water diffusivity	DH2O: 5.717 × 10−5 m2 s−1
Surface tension	σ: 0.0625 N m−1
Vapor condensation latent heat	hcon: 2.308 × 106 J kg−1
Hydronium latent heat	hhy: 3.462 × 106 J kg−1
Entropy change of reaction	ΔS: −149.142 J mol−1 K−1
Dimensionless phase transfer rates	Shcon/evap: 2.04 × 10−3
Specific pore surface area	Apore: 2.0 × 105 m−1
Characteristic length	d: 5.0 × 10−6 m
Phase change rate coefficient	γld: 1.0 s−1
Volume fraction of ionomer in CL	ω: 0.22
Anode exchange current density	$j_{0,\mathrm{a}}^{\mathrm{ref}}$$:1\times {10}^7\exp \left(-1400\left(\frac{1}{T}-\frac{1}{353.15}\right)\right)$ A m−3
Cathode exchange current density	$j_{0,\mathrm{c}}^{\mathrm{ref}}$$:120\exp \left(-\frac{E_c}{R}\left(\frac{1}{T}-\frac{1}{353.15}\right)\right)$ A m−3
Reference H2/O2 concentration	$C_{H_2}^{ref}$$/C_{O_2}^{ref}$: 56.4/3.39 mol m−3
Transfer coefficient	αc/αa: 0.5/0.5
Proton conductivity	$\mathsf{\kappa }{\mathrm{ion}}_{\text{:}}\left(0.5139\lambda -0.326\right)\exp \left[1268\left(\frac{1}{303}-\frac{1}{T}\right)\right]$

.

In the porous zones of GDLs, MPLs and CLs, the flow resistance caused by the porous media is realized by adding viscous resistance S_u described by Darcy’s law to momentum source term. As for the species transport equations, the source terms $S_{{\mathrm{H}}_2}$ and $S_{{\mathrm{O}}_2}$ represent the consumption of H₂ and O₂ in the anode and cathode CLs due to electrochemical reaction respectively and the content is determined by the volumetric current density j_a and j_c. Based on the Butler–Volmer kinetic equation, the electrochemical reaction source terms caused by the hydrogen oxidation reaction and oxygen reduction reaction in the ACL and CCL can be given by

$$j_a={\left(1-s\right)}^5{j}_{0,\mathrm{a}}^{\mathrm{ref}}{\left(\frac{C_{{\mathrm{H}}_2}}{C_{{\mathrm{H}}_2}^{\mathrm{ref}}}\right)}^{0.5}\left[\exp \left(\frac{2{\alpha }_{\mathrm{a}}F}{RT}{\eta }_{\mathrm{act}}\right)-\exp \left(-\frac{2{\alpha }_{\mathrm{c}}F}{RT}{\eta }_{\mathrm{act}}\right)\right]$$

(1)

$$j_c={\left(1-s\right)}^5{j}_{0,\mathrm{c}}^{\mathrm{ref}}\left(\frac{C_{{\mathrm{O}}_2}}{C_{{\mathrm{O}}_2}^{\mathrm{ref}}}\right)\left[-\exp \left(\frac{4{\alpha }_{\mathrm{a}}F}{RT}{\eta }_{\mathrm{act}}\right)+\exp \left(-\frac{4{\alpha }_{\mathrm{c}}F}{RT}{\eta }_{\mathrm{act}}\right)\right]$$

(2)

$${\eta }_{\mathrm{act}}={\varphi }_{\mathrm{ele}}-{\varphi }_{\mathrm{ion}}.$$

(3)

The heat source is mainly composed of entropic heat caused by the electrochemical reaction, Joule heat due to the ion transporting through the PEM and electron conduction in the components of fuel cells and the phase change heat due to the water condensation and evaporation [32,33]. The source term S_liq is composed of the liquid water generation by the electrochemical reaction, phase transition between liquid water and other phases. Phase change is one of the main factors that affects the liquid water profile. The phase change between water vapor and liquid water S_vl can be expressed by [30]:

$${\mathrm{S}}_{\mathrm{vl}}=\{\begin{array}{ll}{A}_{\mathrm{pore}}\frac{{\mathrm{Sh}}_{\mathrm{con}}{D}_{{\mathrm{H}}_2\mathrm{O}}}{d}\epsilon \left(1-s\right)\frac{\left(p_{\mathrm{vp}}-p_{\mathrm{sat}}\right)}{RT}& \hspace{1em}\mathrm{if}{p}_{\mathrm{vp}}>p_{\mathrm{sat}}\\ A_{\mathrm{pore}}\frac{{\mathrm{Sh}}_{\mathrm{evap}}{D}_{{\mathrm{H}}_2\mathrm{O}}}{d}\epsilon s\frac{\left(p_{\mathrm{vp}}-p_{\mathrm{sat}}\right)}{RT}& \hspace{1em}\mathrm{if}{p}_{\mathrm{vp}}<p_{\mathrm{sat}}\end{array}.$$

(4)

Another important factor affecting the distribution of liquid water is the transport characteristics of liquid water. Liquid water transport in the porous electrode including the CL, MPL and GDL can be expressed by the capillary diffusion in the macroscopic models. The diffusivity is given by [30,34]:

$$D_{\mathrm{c}}=-\frac{K_{\mathrm{g}}}{{\mu }_{\mathrm{l}}}\frac{d{p}_{\mathrm{c}}}{ds}$$

(5)

$$p_{\mathrm{c}}=\sigma \cos \theta {\left(\frac{\epsilon }{K_0}\right)}^{0.5}J\left(s\right)$$

(6)

$$J\left(s\right)=\{\begin{array}{ll}1.42\left(1-s\right)-2.12{\left(1-s\right)}^2+1.26{\left(1-s\right)}^3& \hspace{1em}\theta >{90}^{\circ }\\ 1.42s-2.12s^2+1.26s^3& \hspace{1em}\theta <{90}^{\circ }\end{array}.$$

(7)

As for the dissolved water, the source term is caused by the phase change from the water vapor and the transport of dissolved water accompanying with the proton across the PEM due to the electro-osmotic drag effect (EOD).

One factor causing the change of dissolved water content is phase change, and the difference between the equilibrium water content and the local water content is the driven force [35].

$${\mathrm{S}}_{\mathrm{wd}}={\gamma }_{\mathrm{wd}}\frac{{\rho }_{\mathrm{mem}}}{\mathrm{EW}}\left(\lambda -{\lambda }_{\mathrm{equli}}\right)$$

(8)

$${\lambda }_{\mathrm{equli}}=\{\begin{array}{ll}0.043+17.18a-39.85a^2+36a^3& \hspace{1em}\left(0<a<1\right)\\ 14+1.4\left(a-1\right)& \hspace{1em}\left(1\le a\le 3\right)\end{array}$$

(9)

$$a=\frac{p_{\mathrm{vp}}}{p_{\mathrm{sat}}}+2s.$$

(10)

Dissolved water content change due to electroosmotic drag is related to proton conduction and the drag coefficient n_d related to local water content [29]

$${\mathrm{S}}_{\mathrm{EOD}}=\{\begin{array}{l}\frac{{\mathrm{n}}_{\mathrm{d}}{j}_{\mathrm{a}}{M}_{{\mathrm{H}}_2}{}_{\mathrm{O}}}{F}\hspace{1em}\hspace{1em}\mathrm{in}\mathrm{anode}\mathrm{CL}\\ -\frac{n_{\mathrm{d}}{j}_c{M}_{{\mathrm{H}}_2}{}_{\mathrm{O}}}{F}\hspace{1em}\hspace{1em}\mathrm{in}\mathrm{cathode}\mathrm{CL}\end{array}$$

(11)

$${\mathrm{n}}_{\mathrm{d}}=\frac{2.5\lambda }{22}.$$

(12)

As well as the contribution of the above source terms, the diffusion is also an important mechanism that dominates the physical field distribution. To describe the diffusion process in porous media, the Bruggeman correlation is widely used in macroscopic model.

For the diffusion process of gas mixture in porous electrode, the effective diffusivity can be calculated by the following expression [36]:

$$D_{\mathrm{k}}^{\mathrm{eff}}={\epsilon }_{\mathrm{eff}}^{1.5}{D}_{\mathrm{k}}.$$

(13)

Considering the transport path for the electron and proton, the effective conductivity of the electron and proton can be obtained [29].

$$\mathrm{Electron}\mathrm{conductivity}{\kappa }_{\mathrm{ele}}^{\mathrm{eff}}={\left(1-\epsilon -\omega \right)}^{1.5}{\kappa }_{\mathrm{ele}}$$

(14)

$$\mathrm{Proton}\mathrm{conductivity}{\kappa }_{\mathrm{ion}}^{\mathrm{eff}}={\omega }^{1.5}{k}_{\mathrm{ion}}.$$

(15)

The diffusivity of dissolved water can be given by the following correlation [37]:

$$D_{wd}=\{\begin{array}{l}3.1\times {10}^{-7}\lambda \left(\exp \left(0.28\lambda \right)-1\right)\exp \left(-2346/T\right)\hspace{1em}\hspace{1em}0<\lambda <3\\ 4.17\times {10}^{-8}\lambda \left(161\exp \left(-\lambda \right)+1\right)\exp \left(-2346/T\right)\hspace{1em}3\le \lambda <17\\ 4.1\times {10}^{-10}{\left(\lambda /25\right)}^{0.15}\left(1+\tanh \left(\left(\lambda -2.5\right)/1.4\right)\right)\hspace{1em}\lambda \ge 17\end{array}.$$

(16)

3.2. Computational Domain and Boundary Conditions

Due to the periodic geometry configuration as shown in Figure 1, a typical unit is modelled for saving computational resources. Figure 3a shows the computational domain with a cell length of 10 mm. The cathode side contains a whole CC and two half CCs with two symmetrical plane setting at the middle plane of cooling channel. At the cathode inlet side an extended zone with the length of 1 mm is arranged to distribute the air to the GC and CC, respectively. Figure 3b shows the front view of the computational domain marked with the components. The detailed structural and operating parameters for the computational domain are shown in

Table 5. Structural and operating parameters.

Description		Value
GC	Width, WcGC	1 mm
	Depth, dcGC	1 mm
CC	Width, WcCC	1 mm
	Depth, dcCC	1 mm
thickness	BPs, δaBP/δcBP	0.05/0.05 mm
	GDLs, δaGDL/δcGDL	0.5/0.25 mm
	MPLs, δaMPL/δcMPL	0.03/0.03 mm
	CLs, δaCL/δcCL	0.03/0.03 mm
	Membrane, δMEM	0.025 mm
cell length	L	10 mm
Anode side	Inlet temperature, Ta, in	298.15 K
	Inlet relative humidity	0
Cathode side	Inlet temperature, Tc, in	298.15 K
	Inlet relative humidity	0.3
Operating pressure		1 atm

.

The boundary conditions for the anode and cathode flow field are specified with the mass flow rate at the inlet. Such inlet boundary condition can effectively consider the influence of physical property changes caused by environmental parameters on the inlet flow. The outlets of the flow field for both the anode side and cathode side are specified with the absolute pressure of 1 atm to simulate process of discharging exhaust gas to the atmospheric environment. For the anode inlet of air-cooled fuel cell, the mass flow rate based on the reaction demand is provided, and the stoichiometric ratio of 1.4 is given to avoid the shortage of reactants. The concentration C_H2 of the anode inlet can be calculated according to the ideal gas equation.

$${\dot{m}}_a=\frac{{\xi }_a\rho \frac{I_{\mathrm{ref}}}{2F}{A}_{\mathrm{MEM}}}{C_{{\mathrm{H}}_2}{A}_{\mathrm{a},\mathrm{in}}}$$

(17)

$$C_{{\mathrm{H}}_2}=\frac{P_{\mathrm{a},\mathrm{in}}-R{H}_{\mathrm{a},\mathrm{in}}{p}_{\mathrm{sat}}}{R{T}_{\mathrm{a},\mathrm{in}}}.$$

(18)

As for the cathode side, air simultaneously plays the roles of cooling and reaction, namely, the air dual functions [38]. Without a special cooling system, more air is expected to supply the cathode flow field to control the air-cooled fuel cell operating at a proper temperature. The air required for cooling far exceeds the gas required for reaction. Cooling is one of the principal contradictions to be solved by the cathode flow field. Electrochemical reaction in air-cooled fuel cells converts chemical energy into electricity to power loads and the reaction process is inevitably accompanied by heat generation. Assuming that half of the waste heat dissipated by the air forced convection, according to the energy balance, the air flow rate of the cathode side can be given [39]

$${\dot{m}}_{\mathrm{c}}=\frac{0.5\times (U_0-T\frac{d{U}_0}{dT}-V_{\mathrm{cell}})I_{\mathrm{ref}}{A}_{\mathrm{MEM}}}{C_{\mathrm{p}}{A}_{\mathrm{c},\mathrm{in}}(T_{\mathrm{c},\mathrm{out}}-T_{\mathrm{c},\mathrm{in}})}$$

(19)

$$U_0=1.229-0.9\times {10}^{-3}\left(T-298.15\right),$$

(20)

where T_c,in is the actual cathode inlet temperature and T_c,out is a pre-specified cathode outlet set-point temperature. The numerator term in Equation (19) represents half the waste heat produced by the electrochemical reaction. The production of the mass flow rate and the denominator term in Equation (19) is the heat dissipated by the air. By specifying the value of the cathode outlet set-point temperature, the cathode mass flow rate can be obtained. The change of cathode outlet set-point temperature will eventually affect the cathode flow rate. It is worth noting that the actual cathode outlet temperature usually is not equal to the outlet set-point temperature.

Conjugate heat transfer boundary is set for all the fluid and solid interfaces.

$$T_{\mathrm{f}}=T_{\mathrm{s}}.$$

(21)

$${\lambda }_{\mathrm{f}}\frac{\partial T_{\mathrm{f}}}{{\partial }_n}={\lambda }_{\mathrm{s}}\frac{\partial T_{\mathrm{s}}}{{\partial }_n}$$

(22)

The electronic potential for the cathode terminal surface is set to zero and a constant electronic potential V_cell 0.6 V is fixed to the terminal surface of anode terminal. As for the heat transfer the two terminal surfaces of BPs are set as adiabatic boundary condition. A symmetry boundary condition is applied to two end walls in the z-direction of the computational domain as shown in Figure 3b.

3.3. Grid Independent Test

The computational domain is composed of BPs, GCs, CC, GDLs, MPLs, CLs and PEM. Structured grids with hexahedral elements are generated by the software ICEM as depicted in Figure 3a. The governing equations are coupled to each other as discussed in Section 3.1. These governing equations are solved by the commercial software ANSYS FLUENT (ver. 19.2, ANSYS, Inc.: Canonsburg, PA, USA) based on the finite volume method. The SIMPLE algorithm is adopted to tackle with the velocity and pressure coupling. The calculation of conservation equations of electron, proton, liquid water and dissolved water are realized by defining user defined scalars (UDS) in Fluent. The physicochemical parameters, transport coefficients and the boundary conditions are performed by user defined functions (UDF). The second-order differencing scheme is adopted to discretize the convection terms and the diffusion terms. The model was validated against the experimental studies in the previous studies [29,37,40].

In order to ensure the reliability of the results, grid independence test is carried out before the detail research. Figure 3c shows the grid independent test processing. The current density and the cathode outlet temperature under the potential of 0.6 V are selected as the evaluating indicator. With the increase of the gird number, the average current density and cathode outlet temperature increase gradually at first, but when the number of grid system is greater than 270,000, these two parameters hardly change. There is a tiny difference of the current density and cathode outlet temperature can be observed between a grid system with about 270,000 cells and a fine grid system with 400,000 cells. The grid system with 270,000 cell numbers as shown in Figure 3a is employed for the following simulations.

4. Results and Discussion

Although there are many studies focusing on the effect of PTFE content and distribution in GDLs of water-cooled fuel cells, for air-cooled fuel cells the water management requirements are different, and thus GDL design with various PTFE doping strategies for water-cooled fuel cells may be not suitable for air-cooled fuel cells. In this section, effects of the PTFE doping strategies proposed in Section 2 on the distribution of key variables and performance of fuel cells are studied. The cathode outlet set-point temperature is fixed as 315.65 K to maintain the same mass flow rate for all the cases. These designs are compared under the same output voltage of 0.6 V. For the convenience of the description and subsequent comparative studies,

Table 6. Categories of PTFE doping strategies.

Categories	Case No.	PTFE Doping Content (wt%)
Uniform scheme	1	40%
	2	30%
	3	20%
	4	10%
Through-plane Gradient scheme	5	10% to 40% from GC/GDL interface to GDL/MPL interface
	6	40% to 10% from GC/GDL interface to GDL/MPL interface
In-plane Sandwich scheme	7	40% under GC and 10% under CC
	8	10% under GC and 40% under CC

summarized the various designs with different PTFE doping strategies. Cases 1–4 contain four designs with PTFE uniformly doping into the GDL. Cases 5–6 represent two GDL designs with through-plane gradient PTFE doped. These designs are formed by linearly increase or decrease the PTFE content along the GDL thickness direction. Such linear distribution of PTFE makes the porosity and permeability of GDL also linearly distributed along the thickness direction of GDL. Cases 7–8 show two in-plane sandwich doping GDL designs. The difference between the two designs lies in the different PTFE content in GDL under the CC and GC.

4.1. Effects of the PTFE Content

Case 1 is the design with 40% weight percent PTFE doped into GDL among the four uniformly doping designs, and it is used as the baseline design. Figure 4 shows the distribution of several key variables in the central plane along the in-plane direction of cathode flow field for Case 1. As shown in Figure 4a, air is uniformly distributed into the GC and CC, and develops into a parabolic profile. The average velocity in GC and CC in the study is about 3 m/s and the equivalent diameter of the channel is 1 mm, so the Reynolds number is about 200. Therefore, the laminar flow model adopted is reasonable. Along the air flow direction, the pressure drops almost linearly as shown in Figure 4b. The air flowing in the GC and CC absorbs the waste heat generated by the reaction, leading to continuously increasing temperature as shown in Figure 4c. Due to higher thermal conductivity of BP the temperature near the BP is higher than the main flow of CC and GC. Figure 4d shows the distribution of oxygen concentration. The concentration drops along the flow direction. In the CC, air is isolated from the reaction and oxygen is not consumed, and the decrease of oxygen concentration is mainly determined by the temperature increase and the pressure decrease according to the state equation for ideal gas. In the GC, the main reason for concentration reduction is the consumption caused by the electrochemical reaction. Note that in the near CC region of the GC, the decrease of oxygen concentration is more obvious. The reason can be attributed to the higher reaction rate under the CC which will be discussed later. The water vapor concentration distribution in the cathode flow field shows different distribution trend in GC and CC, respectively, as shown in Figure 4e. In the CC, the vapor concentration decreases slightly along the flow direction. Careful observation will find that the water vapor concentration is slightly lower near the rib side. The main reason can be attributed to the increase of the saturated vapor pressure caused by the temperature increase as displayed in Figure 4c and the vapor partial pressure decrease due to the pressure drop of gas mixture as depicted in Figure 4b. In the GC, the water vapor concentration exhibits an increase trend along the flow direction, as water is continuously generated due to the electrochemical reaction. Through evaporation, diffusion and forced convection, water vapor passes through the GDL and continuously accumulates along the GC, resulting in the increase of water vapor content along the flow direction. Similarly, due to the effect of electrochemical reaction, more water is generated near the rib, leading to a higher water vapor concentration near the rib. Although the water vapor accumulates near the rib side in the GC, the relative humidity does show a different trend as shown in Figure 4f. The rise of saturated vapor pressure caused by high temperature near the rib side is the main reason for this phenomenon. It can be found that the above analysis in the flow field shows the complex variation of oxygen concentration, water vapor concentration and the relative humidity affected by electrochemical reaction, temperature and pressure.

The change of the above physical parameters in the cathode flow field is highly related to the complex transport and chemical reaction processes inside the electrodes. Figure 5a–d respectively compare the distribution of oxygen concentration, temperature, dissolved water content and current density in the middle plane of CCL for the four PTFE uniformly doping schemes in Cases 1–4. For all the cases, along the flow direction, the oxygen concentration decreases as it is consumed by the electrochemical reaction. It can also be found that the oxygen concentration is higher under the GC than that under the CC. The higher mass transfer resistance and faster electrochemical reaction rate under the CC are the main reasons for this phenomenon. Additionally, with the decrease of PTFE content in GDL from Case 1 to Case 4, the oxygen concentration in CCL presents a downward trend, which is due to stronger electrochemical reaction as will be explained in Figure 5d. In Figure 5b, the air continuously absorbs waste heat along the flow direction leading to the temperature rise in the MEA from the inlet to the outlet), and Case 4 with the lowest PTFE content presents the highest temperature. Such a high temperature is expected to render more water evaporated resulting in lower dissolved water inside the CCL. However, the dissolved water distribution displayed in Figure 5c shows the opposite trend. The above result results from the fact GDL with less PTFE doping is characterized with lower hydrophobicity and higher porosity, leading to enhanced water retaining ability of the GDL. With more water retained in the electrodes, the proton conductivity increases, and thus the reaction rate is improved due to the reduction of ohmic loss, leading to more water being generated. In turn by the phase change more dissolved water will be retained in the CCL, leading to the positive feedback between higher water content and stronger electrochemical reaction. Due to this positive feedback, the current density distribution in CCL as shown in Figure 5d is positively related with the distribution of dissolved water. The current density is higher under the rib than that under the channel, and it increases with the decrease of PTFE doping content. Careful observation will also find that the current density near the cathode inlet is lower than that at the downstream. The phenomenon is mainly caused by the insufficient hydration of membrane near the inlet due to a large amount of dry air entering the GC. With the water produced by the reaction accumulating at the downstream, the current density is significantly increased.

4.2. Effects of Through-Plane PTFE Gradient Doping

In this section, designs of gradient PTFE doping along the GDL thickness direction are studied. In Case 5, the PTFE amount is linearly increased from GC/GDL interface to the GDL/MPL interface, leading to lower porosity and higher hydrophobic region near the GDL/MPL interface. Note that for water-cooled fuel cell, such PTFE doping strategy is beneficial for the water removal and thus can improve the cell performance [41]. In Case 6, the variation of the PTFE content is the opposite of that in Case 5. To clarify the effects of through-plane PTFE doping strategies on the coupling mechanisms of internal physical quantities and fuel cell performance, distributions of key parameters in the middle plane of CCL as shown in Figure 6 are compared and analyzed.

As depicted in Figure 6a, different designs share the similar features of lower oxygen concentration under the rib and downstream of flow direction. The reasons can be attributed to uneven electrochemical reaction rate distribution and mass transfer resistance under rib and channel as discussed in Section 4.1. Despite the common distribution characteristics, the oxygen concentration in CCL of Case 6 near the outlet of the GC is significantly lower than that in Case 5. Although the inlet flow rate for Case 5 and Case 6 is the same, the temperature distribution is quite different as shown in Figure 6b. The temperature of CCL in Case 6 is significantly higher than that in Case 5 and the hot zone moves towards the cathode inlet. The maximum temperature in Case 6 is nearly 7 K higher than that in Case 5. However, compared with Case 5 there is no obvious dehydration problem in the electrode of Case 6 as shown in Figure 6c. The contradictory is mainly caused by the configuration of GDL. GDL of through-plane PTFE gradient doping with more PTFE doping into the near GC zone adopted in Case 6 shows a good water retaining capacity. Due to the high proton conductivity caused by the high dissolved water content, Case 6 shows a higher current density compared with Case 5. Higher current density in turn will lead to more water production and higher proton conductivity, thus forming positive feedback as well.

To reveal the water retaining mechanism of different designs, the relative humidity and liquid water saturation distribution in the through-plane cross section of the cathode MEA are further discussed as shown in Figure 7. In the in-plane direction, it can be found that the relative humidity and liquid water saturation in the MEA under the CC are higher than that under the GC for both designs. The MEA under the GC is directly exposed to the ambient air with lower relative humidity. Higher relative humidity difference between the GC and MEA intensifies the evaporation of water in the MEA under the GC. In addition, the mass transfer path of water vapor under the CC is longer than that under the GC. Under the combined effects of the above factors, the relative humidity and liquid water saturation under the CC are higher than that under the GC.

Now attention is turned to the through-plane direction. The liquid water saturation shows a decrease trend from the CL to the GDL. For the relative humidity, the distribution is a little complicated. Along the direction from the CCL to GC, the relative humidity presents an increase trend under the CC, however, under the GC the change of relative humidity shows an opposite trend. For the former case, the dominant factor is the local temperature distribution from the CCL to GC. An elevated temperature near the CCL will make the increase of local saturated vapor pressure and reduce the local relative humidity. For the second case, lower humidity in the GC is responsible for the humidity decrease in GDL near the GC.

The above discussion reveals the common distribution features of both Case 5 and Case 6. Comparing the two through-plane PTFE doping schemes, it can be obviously observed that the relative humidity and the liquid saturation of Case 6 are higher than that of Case 5. On the one side, lower hydrophobicity of GDL near the GDL/MPL interface in Case 6 design makes the discharge of liquid water more difficult than Case 5 design. The liquid water content near the CL increases obviously for Case 6 as shown in Figure 7b. Under a higher temperature, the liquid water evaporates into water vapor leading to a higher water vapor content in MEA of Case 6 than that of Case 5. On the other hand, lower porosity of GDL near the GC/GDL interface will weaken the water vapor discharge by the strong forced convection in GC. Under these two favorable factors, more water vapor content finally leads to the increase of relative humidity as shown in Figure 7a. Higher relative humidity and liquid water saturation can improve water activity as shown in Equation (9), thus increasing the equilibrium water content. Large difference between local dissolved water content and equilibrium dissolved water content will drive water vapor transforming to the dissolved water and then improving the dissolved water content in the CL. Based on the above analysis, Case 6 design with a linearly increase PTFE content along the direction from GDL/MPL interface to the GC will strengthen the water retention of air-cooled fuel cells.

4.3. Effects of Sandwich PTFE Doping

In this section, effects of the in-plane sandwich distribution of PTFE in the GDL are studied. Figure 8 compares the local distribution of the key variables for these two different designs in the middle plane of CCL. Along the flow direction, the dissolved water content decreases slightly for both designs as shown in Figure 8a due to the temperature rise. Careful observation will find that the dissolved water content of Case 7 is slightly lower than that of Case 8 near the outlet region. The slightly increased dissolved water content in Case 8 also improves its current density as shown in Figure 8b. In order to reveal the water retaining mechanism of these two PTFE doping strategies, the distribution of relative humidity and liquid water saturation in the middle plane of CCL is analyzed. Case 7 (less PTFE content under the rib) shows a higher liquid water saturation under the rib than Case 8 design as shown in Figure 8c. The lower hydrophobicity under the rib of Case 7 can help to keep more liquid water retained in the MEA. However, more liquid water is expelled under the GC for Case 7 due to more hydrophobicity of GDL thus leading to a more uneven distribution of liquid water saturation than Case 8. Although the maximum liquid water saturation of Case 8 is lower than Case 7, Case 8 shows a more uniform liquid water distribution and a higher liquid water content under the GC. Higher liquid water content under the GC evaporates to more water vapor under high temperature, which slightly improves the relative humidity under the GC to a certain extent as shown in Figure 8d. Due to a more uniform distribution of liquid water and water vapor in the CL of Case 8, the dissolved water content is slightly improved, which improves the current density of the fuel cell as shown in Figure 8b.

4.4. Effects on the Fuel Cell Performance and Current Density Distribution Uniformity

As discussed above, the doping strategies of PTFE have a significant influence on the distribution and coupling mechanisms of internal physical quantities of air-cooled fuel cells. The distribution of these physical quantities will eventually affect the performance of fuel cells. Figure 9 shows the comparison of current density of different doping strategies. In the above analysis, it is found that the dissolved water content in the CL has a significant impact on the performance of air-cooled fuel cell. Therefore, Figure 9 also presents the comparison of the average dissolved water content in the middle plane of the CCL. The air-cooled fuel cell performance shows a downward trend with the increase of the PTFE content in the uniform doping group as shown in Figure 9. The current density of Case 4 with 10% PTFE content in GDL is about 18.2% higher than that of Case 1 at V_cell = 0.6 V. The performance of through-plane PTFE doping designs presents two extreme states. Of various designs in the study, Case 5 has the worst performance, while Case 6 exhibits an outstanding performance. Compared with Case 1, the current density of Case 5 decreases by 6% and that of Case 6 increases by about 22%. The two sandwich doping designs (Case 7–8) show better performance compared with Cases 1–3 designs, but the performance is not as good as that of Case 4. In the study, the air-cooled fuel cell operated at 0.6 V, and the performance of fuel cell is mainly controlled by ohmic polarization. Ohmic polarization is affected by the proton transport process, and the proton conductivity is determined by the dissolved water content in the membrane. It can be found that there is an obvious positive correlation between the current density and the average dissolved water content in the CCL as shown in Figure 9. The design of GDL has an important impact on water retaining capacity, and then affects the performance of air-cooled fuel cells. The gradient PTFE doping strategy of Case 6 shows a higher water retaining capacity compared with the other designs as shown in Figure 9, thus improving the performance of air-cooled fuel cells as discussed above.

As well as the performance of fuel cell, the uniformity of key variables [42,43,44,45,46], such as species concentration, temperature, voltage and current density, is also a concern factor to be studied. To examine effects of PTFE doping strategy on the uniformity of current density distribution in the CL, the following uniformity factor is defined

$$U_i=1-\frac{1}{I_{\mathrm{ave}}}\sqrt{\frac{1}{A_{\mathrm{MEA}}}{\displaystyle \iint {\left(I-I_{\mathrm{ave}}\right)}^{2}dA}},$$

(23)

where I_ave is average current density of the middle plane of CCL and A_MEM is the area of MEA.

Figure 10 shows the comparison of the current density uniformity of various designs. Higher value of uniformity factor means a more uniform distribution of current density. The groups of PTFE uniform PTFE doping (Cases 1–4) and through-plane gradient PTFE doping (Cases 5–6) show relatively high uniformity, while the group of PTFE in-plane sandwich doping (Cases 7–8) shows a slightly lower current density uniformity. In-plane sandwich doping design produces differences of the wettability and porosity under the CC and GC respectively, which intensifies the uneven distribution of dissolved water content. The uneven dissolved water content will aggravate the uneven current density distribution of in-plane sandwich design. Compared with the in-plane sandwich doping designs, the current density distribution uniformity of PTFE uniformly doping design is slightly improved. With the decrease of PTFE content in the group of uniformly doping designs, the uniformity factor shows a decreasing trend. Case 5 shows a higher uniformity among various designs. This is because Case 5 has the worst water retaining capacity, thus leading to a small difference of dissolved water content under the CC and GC. Although the current density uniformity of Case 6 is not as good as that of Case 5, it is equivalent to that of uniform PTFE doping designs. Therefore, Case 6 design not only improves the fuel cell performance, but also maintains a good current density distribution uniformity.

5. Conclusions

In this study, a three-dimensional two-phase fuel cell model was employed for complex transport processes in air-cooled fuel cells. The effects of PTFE doping strategies of GDL on the air-cooled fuel cell performance are comparatively studied. Through comparative study, GDL design with gradient PTFE doping shows a better performance in the water retention of air-cooled fuel cells. The results of the study can be concluded as follows:

(1)

The results show that the GDL design idea of water-cooled fuel cells cannot meet the demand of air-cooled fuel cells. In the PTFE uniformly doping designs, properly reducing the amount of PTFE in GDL is beneficial to promoting the performance of air-cooled fuel cells by improving the water retaining capacity;

(2)

Compared with the design of GDL with 40% PTFE uniformly doped, the in-plane sandwich design of GDL can improve the performance of fuel cells. However, the performance of air-cooled fuel cells with in-plane sandwich GDL design is inferior compared with the uniform PTFE doping design with the least content of PTFE. Between the two sandwich designs, less PTFE doping in GDL under the GC holds a slight advantage for improving the performance of the fuel cell;

(3)

Water retention capacity of the air-cooled fuel cell can be improved by linearly increasing the PTFE content from the GDL/MPL interface to the GC/GDL interface. This design effectively increases the water content of the membrane electrode assembly and, compared with the design of GDL with 40% PTFE uniformly doped, the current density is improved by 22% when operating at 0.6 V;

(4)

Finally, the effects of different GDL designs on the uniformity of current density distribution are compared. Air-cooled fuel cells with in-plane sandwich PTFE doping GDL perform poorly in the uniformity of current density distribution. The through-plane gradient PTFE doping design with the PTFE content linearly decreasing from the GDL/MPL interface to GC/GDL shows the best current density uniformity by sacrificing the performance. In addition to improving performance, a design with linearly increasing PTFE content along the GDL/MPL interface to the GC/GDL interface can also maintain a good uniformity of current density distribution in CCL.

Through numerical simulations, feasible GDL designs to improve the performance of air-cooled fuel cells are proposed. However, these designs inevitably lead to the increase of cost and manufacturing complexity. Future works need to focus on GDL designs that are easier to process and manufacture in order to facilitate the commercialization process.