Effect of through-plane and in-plane permeability of gas diffusion layers on the performance of PEM fuel cell

A three-dimensional model with a complete set of governing equations for all components of the PEM fuel cell is developed to investigate the permeability effects on the transport of reactive species and the flooding phenomenon in PEM (Proton Exchange Membrane) fuel cell. This model couple fluid flow distribution, species, electric potential, electric current density, and saturation level in the flow channel, gas diffusion porous layer, catalyst layer, and fuel cell polymer membrane. The electrodes and channels are modeled as two-phase. The effects of isotropic and non-isotropic permeability have been studied and the simulation results show that in order to obtain higher efficiency, the permeability of the gas diffusion layer should be high in the flow direction(in-plane) and perpendicular to the main flow (through plane). High permeability in the direction of the main flow and Low permeability in the direction perpendicular to the main flow will also have a high yield. The current study shows that by optimizing the permeability of the gas diffusion layer, a higher efficiency can be achieved in the fuel cell.


Introduction
Fuel cells are a type of electrochemical energy conversion device that directly converts the chemical energy of the fuel into DC electrical energy through electrochemical reactions.Unlike batteries, which are a source of energy storage, fuel cells are able to operate continuously as long as fuel is available.The combination of high efficiency and environmental benefits make this device a good replacement for polluting energy conversion devices such as internal combustion engines.Fuel cells have a lot of ability and potential to influence future energy production, but before this happens, some technical limitations for fuel cells must be overcome [1].
Observing the phenomenon of transport of chemical species in the fuel cell and measuring it by experiment is a difficult task due to the small dimensions of the cell components.Computational model is a good alternative method in studying the efficiency of fuel cells [2].In the early 1990s, most of the models were one-dimensional and isothermal and focused on the catalyst layer and membrane [3][4][5][6][7][8][9] In the late nineties, the models advanced and included two-dimensional, three-dimensional and two-phase model [10][11][12][13][14][15][16][17][18][19].Many articles about PEM fuel cell and modeling of saturation and flooding phenomenon and the proposed solutions to overcome it have been published in recent years [20][21][22][23][24].
Fuel cell models are dependent on special physical properties that need to be accurately determined in order to properly simulate fuel cell performance.Some of these properties include kinetic parameters such as exchange current density and catalyst surface area per unit volume.Structural information such as porosity and tortuousness of the porous medium related to the gas diffusion layer and polymer properties such as ionic conductivity and diffusion coefficient of water dissolved in it.One of the most important aspects of a fuel cell is water control that describes how water moves inside the cell.In PEM fuel cell, reactant gases have to pass through a porous layer named 'gas diffusion layer' (GDL) for access to the catalyst layer, thereby making the electrochemical reaction possible and the subsequent energy production.Positioned between the bipolar plate and the catalyst layer, the GDL exhibit thermal and electrical conductivity and allows gas and water transport.
Gas transport of reactants across GDLs is both through-plane and in-plane in fuel cells as explained below.Throughplane transport refers to transport perpendicular to the plane of the GDL and the membrane electrode assembly (MEA), from the gas channel to the electrode.In-plane transport corresponds to gas flow parallel to the above planes.The permeability of a porous medium is the ability of the medium to facilitate the flow of fluid through its open pores.Higher permeability of gas across the GDL would lead to a superior reaction at the catalytic site.researchers in studies [25][26][27][28] used computational fluid dynamics to determine the importance of convective transport as a function of GDL permeability.Both isotropic and orthotropic permeability are considered and it is found that the in-plane permeability is the parameter of importance. in study of Zhao [29], the contribution of the catalyst layer to the overall gas permeability of the electrode is experimentally investigated for different catalysts with a range of Pt loadings at various temperatures for air, oxygen and nitrogen gases.Results indicate that the gas permeability of the GDLs can be reduced with the presence of a catalyst layer.Modeling of Belkhiri [30] visualized the effect of permeability on the dynamic behavior of the velocity and pressure profiles in the parallel flow field in the transient model of PEMFC fuel cell.Mokherjee [31], focused on measuring the permeability across gas diffusion layers (GDLs) in PEM fuel cell configuration with varying bi-polar plate designs.The results showed a significant relationship between the equivalent permeability and the fuel cell performance.
Comprehensive study of fuel cell with consideration of all component of fuel cell is necessary.As a novelty of research, this study deals with three-dimensional, two-phase non-isothermal model of PEM fuel cell to determine the effect of isotropic and Non-isotropic gas diffusion medium on saturation, flooding phenomena and overall performance of cell.

Porosity
The amount of porosity is equal to the ratio of the volume of pores in the porous material to the total volume of the porous material.

Saturation
In two-phase systems, the volume fraction of the total porosity occupied by one phase is called phase saturation.Phase saturation,  j is one of the key parameters in multiphase flow.

Permeability
Permeability is defined as the tendency of a porous medium to pass fluid.Permeability of the porous medium is one of the control factors that determines the rate of fluid penetration inside the porous medium.
Capillary phenomenon in the porous medium of the fuel cell: When two or more fluids are present in the porous medium at the same time, a pressure difference is created at the boundary between the phases due to surface tension, which is caused by the molecular imbalance in the contact line of the phases.The pressure difference at the common boundary of both phases is known as the capillary phenomenon.The capillary pressure between the wetting and non-wetting phases is obtained by the following relationship [15] . c = P l − P g = σ cos(θ c ) (

Water control
The key factor in the efficiency of a polymer membrane fuel cell is the water balance in it.There is a complex relationship between the water content of the cell and its performance.In order for the polymer membrane to have good ionic conductivity, it is necessary for the moisture content inside it to be high, on the other hand, the accumulation of liquid water in the electrodes prevents the reactive gases from passing through the GDL and reaching the catalyst surface.

Net transfer coefficient of water
The net transfer coefficient of water molecules from the membrane is equal to the sum of three water transfer mechanisms.Based on the empirical relationships obtained by Springer [4], the water content in the membrane is obtained by the following relationship.
λ = { 0.043 + 17.81 − 39.85 2 + 36.0 3  0 <  ≤ 1 14 + 1.4( − 1)  1 ≤  ≤ 3 …….(7) In this regard,  is a parameter that shows the activity level of water molecules and actually shows the relative humidity inside the cell and is defined by the following relationship: is the saturated vapor pressure and is expressed by the following equation The vapor pressure of water is defined based on the mole fraction of the vapor and the total pressure:  wv = x H2O P …….. (10) The ionic conductivity of the membrane depends on the water content and the temperature of the fuel cell and is expressed by the following relationship: The electroosmotic coefficient is expressed by the following equation [15]  d = 2.5 ) ……….. ( 14)

Electrochemical reactions
In the catalyst layer, hydrogen is oxidized in the vicinity of platinum particles.The general form of hydrogen oxidation reaction is as follows: In the cathode catalyst layer, oxygen is reduced in the form of water molecules in the vicinity of platinum particles and in the presence of hydrogen ions and electrons.16)

Electrochemical kinetics
Electrochemical kinetics is a complex process and here considers only a simple explanation.The operation of an electrochemical system is a non-equilibrium process that involves the transfer of electrons and reactive species at the electrode surface.The reaction rate is directly dependent on the faradic current through the electrodes.This rate is dependent on three important parameters: exchange current density which is dependent on the catalytic activity on the electrode surface, concentration of reactive species on the electrode surface and activation potential drop on the catalyst surface.Mathematically, these physical quantities are combined in the Butler-Volmer equation to calculate the electric current [15]:

Cathode current
where j cat ref is the base exchange current density of the cathode and α cat is the charge exchange coefficient of the cathode.) …………………….(18) According to the principle of electric current stability, the average current density in the anode is equal to the average current density in the cathode.

Anode current
…………………….(19) which is the total number of cells in the anode and cathode catalyst layer and the cell volume.

Potential drop in fuel cell
When electric current begins to flow in the cell, the cell potential drops below the reversible potential due to several types of potential drops, including activation potential drop, ohmic potential drop, concentration related potential drop, and mass transfer potential.The details of these drops are given in the following section

Electrochemical modeling
Calculating the rate of hydrogen oxidation and oxygen reduction is important in electrochemical modeling.The electrochemical process occurs on the surface of the catalyst layer on both sides of the membrane.The driving force behind these reactions is the difference between the potential of the solid phase and the electrolyte or membrane phase, which is also called the activation potential drop.can be Therefore, two potential equations are solved in the fuel cell model: one potential equation for electron transport through solid materials and the other equation represents proton transport from the membrane.The equations related to the electric potential in the cathode, anode and membrane are as follows, by solving which the voltage in the battery can be calculated.
For the solid phase potential, there are external boundaries on the anode and cathode sides that are in contact with the external electric circuit, and only through these boundaries does the generated electric current pass.In other external boundaries, zero flux boundary condition is considered for   and in external contact boundaries for   we consider constant values.If the voltage on the anode side is considered zero, the determined value for the cathode side will be the cell voltage.
The exchange current or production term in equation ( 20) and ( 21) is non-zero only inside the catalyst layers.
For the solid phase on the anode side  sol = −R an (< 0) For the solid phase on the cathode side  sol = +R cat (> 0) For the membrane phase on the anode side  mem = +R an (> 0) For the membrane phase on the cathode side  mem = −R cat (< 0)

Potential drop in the anode
The overall potential drop on the anode side is obtained from the following relations:  an = Φ sol − Φ mem …………………….(22) Potential drop in the cathode: The overall potential drop on the cathode side is obtained from the following relations:

Flow in the gas channel
The flow in the channel is modeled as a gas flow in a closed channel.Although laboratory observations show that twophase flow prevails in the gas channel, most models consider single-phase channel flow to simplify.The governing equations are as follows.Here i indicates oxygen or water vapor on the cathode side and hydrogen or water vapor on the anode side.

gas diffusing layer
The transport in the gas diffusing layer is modeled as the transport in the porous medium.) 1.75 …………………….(34) The heat transfer in the GDL layer is modeled by the following equation: k eff which is the effective thermal conductivity coefficient and  h includes the heat source and the sink in the porous medium.The value of  h in the GDL layer is equal to:

Catalyst layer
In most models, the catalyst layer is modeled as a thin dimensionless layer, here the dimensions of the catalyst layer are considered and the flow equations inside this layer are solved.The governing equations in the catalyst layer are similar to the equations of the gas diffusing layer and the differences appear in the properties of the porous layer and the terms of production and consumption.

Mass transfer equation in CL
The chemical reactions that occur inside the catalyst layers are heterogeneous reactions that occur on the surface of the catalyst in the porous material.Therefore, the concentration of hydrogen and oxygen species are surface values.The reactions in the two catalyst layers are surface reactions and it is assumed that the diffusive flux of each reactant is balanced by the production rate.The catalyst layer is the most important part of the computing domain.This area is where the electrochemical reaction takes place and the species in the fuel cell are produced or consumed.Generally, the terms sink and source are used to model the consumption or production of chemical species.In the cathode catalyst layer, oxygen consumption is modeled by the well expression in terms of the local current density in the cathode.
Water production is modeled as a source term and according to the local current density in the cathode, which is closed with the water passing through the membrane due to the effect of electro osmosis.
In the anode catalyst layer, water flows towards the cathode due to electro-osmosis phenomenon and creates a sink expression.
In the anode catalyst layer, hydrogen is consumed to produce electrons and protons.Hydrogen consumption is modeled as a sink term similar to oxygen consumption.
Heat production due to electrochemical reactions on the cathode side is calculated by the following equation.

Membrane
The membrane is assumed to be completely wet.The overall water flux is modeled as an electro-osmotic drag of water from the anode to the cathode and back diffusion from the cathode to the anode, and the role of back diffusion due to the pressure difference can be neglected compared to the electroosmotic phenomenon.
In the membrane, only the energy equation and the electric potential equation are solved.The electric potential equation in the membrane was described in previous section.The ion conductivity coefficient in the electric potential equation in the membrane is dependent on the water content of the membrane.

Two-phase relations in the cathode catalyst layer and the GDL cathode layer:
To obtain the saturation distribution in the cathode catalyst layer and GDL, the saturation equation is solved.In the porous region where there is a high resistance to flow, the capillary diffusion term can be used instead of the convection term in equation .The equation ( 47) acts as a switch and controls the condensation by comparing the concentration of water vapor with the concentration of saturated vapor.

Calculation process and solution algorithm
All the equations are discretized based on the finite volume method and turned into a system of algebraic equations that are solved by the iterative algorithm.The computational domain is divided into a number of control volumes (cells).The SIMPLE algorithm has been used to solve the flow and pressure field.The SIMPLE algorithm was first presented by Patankar [32,33].The equations were solved by FLUENT software.

Computing domain
The computational domain of the cell includes the cathode current channel, the anode current channel, the gas diffusion layer on the cathode side, the gas diffusion layer on the anode side, the catalyst layer on the cathode side, the catalyst layer on the anode side, and the polymer membrane, which is shown in Figure (1).

Boundary conditions
Boundary conditions are applied on the outer surfaces of the computational domain.
The inlet velocity is calculated as follows.where n e is the number of electrons per mole of reactive species, for oxygen in the cathode ne=4 and for hydrogen in the anode ne=4, ξ is the stoichiometric ratio of current.P i is the total gas pressure and x i is mole fraction of hydrogen or oxygen at the anode or cathode.

Outlet
Newman's boundary condition is used for other variables in the outlet, i.e. zero normal gradient for velocity, temperature, species mass fraction and potential field.Mathematically, this boundary condition is expressed as follows.

External surfaces of the cell
A constant temperature condition has been applied in the upper and lower surface of cell. = 353

Independence from the number of computational grids
In order for the numerical results to be independent of the number of computing grids, the solution was done with several different meshes and it was observed that the electric current (for the voltage of 0.7) reaches convergence in the meshes size of 40x60x70 and the meshes refining after this meshing does not affect the results.Simple square meshes are used.The size of the meshes is the same in the x direction, but it is different in the y direction depending on the size of the different components of cell.Limited laboratory results can be found in the articles that fully describe the geometry of the pile and its parameters and test conditions.The dimensions used in the computing domain in this model are given in table (1).The working conditions of the fuel cell are given in table (2).physical properties of electrodes are given in table (3).

Comparison with laboratory results
To verify the accuracy of computer and simulated models, the polarization curve of the model was compared with the polarization curve obtained from laboratory results prepared by Ticianelli [34] and the accuracy of this simulation was confirmed.Figure (2) shows acceptable agreement between the model polarization curve and the laboratory polarization curve.

Results
In addition to the polarization curve, the model allows obtaining important information about the details of the transfer phenomenon inside the fuel cell.These transport phenomena include local distribution of reactant concentration, temperature field, local current density, and potential field.

Permeability effects Table 4 Permeability compounds for the model
In-plane permeability Through-plane permeability The permeability of the gas diffusion layer is one of the parameters that affects the efficiency of the fuel cell.the importance of permeability in the main flow direction of the channel has been studied in some articles and it has been shown that permeability values are important in the flow direction.The experimental evidence of the permeability of the gas diffusion porous layer in different materials has been investigated in the direction of the channel flow (z axis) [30].In this research, a comparison has been made between isotropic and non-isotropic diffusion layer and its effects on fuel cell efficiency.The effects of isotropic and non-isotropic permeability have been studied using the values given in table (4).
Items P11, P22, P33 indicate isotropic permeability.Therefore, cases P11, P12, P13, P23 are closer to the permeability of GDLs used in real fuel cells.The simulation was done for voltages of 0.5V and 0.75V and the calculated current density is shown in Tables (5) and (6).As shown in table (5) at 0.75V voltage, the permeability effect is almost negligible because at this voltage, the fuel cell works in single-phase mode and saturation and condensation start from 0.65V.At a voltage of 0.5 V, permeability effects appear.For isotropic cases P11, P22, and P33, the decrease in permeability reduces the efficiency of the fuel cell.For the non-isotropic cases P11, P12 and P13, where the permeability in the direction of the channel flow is kept at a constant value and the permeability is variable in the direction perpendicular to the flow, decreasing the permeability in the directions perpendicular to the main flow increases the current density from 1.77 to 1.93 It can be For cases P11, P21, P31, where the permeability perpendicular to the flow is kept at a high value and the permeability changes in the direction of the main flow of the channel, the effects are significant and current density decreases from 1.77 to 1.37.These results show that the permeability of GDL must be isotropic and high value, or must be non-isotropic and its permeability be high in the direction of channel flow.Non-isotropic GDL with high permeability in the direction perpendicular to the flow and lower permeability in the direction of the flow causes The efficiency of the fuel cell decreases.Comparison of cases P13 and P31, which have current densities of 1.93 and 1.37, respectively, also shows the advantage of large permeability in the direction of flow channel over permeability in the direction perpendicular to the channel.
Taking a precise attention at the permeation effects on the fuel cell according to the liquid water profile (saturation) at the GDL/CL interface at the outlet of the cell Figure (3).we can get significant results.For the permeability case p13, we have the lowest saturation and for the P31 case, we have the highest saturation.It can be concluded that different permeability compounds have different effects on the accumulation of liquid water at the end of the fuel cell, so the P13 case causes a better discharge of the accumulated or saturated water and thus facilitates the presence of oxygen in this area to continue the chemical reaction.This can be seen from the amount of electric current density, the current density of P13 case is equal to 1.93 while for P31 case this value is 1.37.Therefore, there are two possible reasons behind these phenomena: First, less oxygen penetrates due to less permeability.Second, the presence of liquid water prevents oxygen from reaching the reaction surface.The water removal rate is better in P11, P12 and P13 cases.As can be seen in figure (3), the saturation level in the GDL and CL part of the cathode for the P13 case is much lower than the P31.cases P11, P12 and P13 are more optimal than others, which is also evident from the current and voltage values.

Discussion
3D computational fluid dynamics model for PEM fuel cell was used in this study.This model provides useful information about transport phenomena in the fuel cell, such as reactant concentration distribution and saturation distribution.In addition to revealing basic information about the transfer phenomena inside the fuel cell, this simulation is able to predict the effects of working and design parameters on the PEM fuel cell efficiency.The results obtained in this research showed that the permeability of the gas diffusion layer is one of the parameters that affects the efficiency of the fuel cell.
The effects of isotropic and non-isotropic permeability were studied.These results showed that the permeability of GDL should be isotropic and its value should be high or it should be non-isotropic and its permeability should be high in the direction of channel flow.Non-isotropic GDL with high permeability in the direction perpendicular to the flow and lower permeability in the direction of the flow causes a decrease in the efficiency of the fuel cell.The comparison of the cases shows the advantage of the large permeability in the direction of the channel flow over the large permeability in the direction perpendicular to the main flow.
Taking a precise look at the permeation effects on the fuel cell according to the liquid water profile (saturation) at the GDL/CL interface at the end of the cell, we achieve significant results.It can be concluded that different permeability value has different effects on the accumulation of liquid water at the end of the fuel cell, so that the large permeability of the diffusion layer in the flow direction compared to the permeability in the direction perpendicular to the channel causes a better discharge of saturated water and oxygen in this area and facilitates the continuation of the chemical reaction.This can also be seen from the electric current density.As seen, the use of materials with appropriate permeability can reduce the phenomenon of flooding in PEM fuel cells and increase its efficiency.

Conclusion
In this study, a numerical investigation has been carried out to determine the effect of through-plane and in-plane gas permeability of the gas diffusion layers in PEM fuel cells.These results highlight the importance of non-isotropic permeability of gas diffusion layers, hence the large permeability of the diffusion layer in the flow direction compared to the permeability in the direction perpendicular to the channel causes a better discharge of saturated water and oxygen in this area and facilitates the continuation of the chemical reaction.The permeability analysis reported in this study can be used for understanding of the physical processes involved in an operating PEM fuel cell, and will be also useful for fuel cell industry in their design calculations of PEM fuel cells meeting specific performance targets.

Disclosure of conflict of interest
No conflict of interest to be disclosed.
.surf − y i.cent )r = S k …………………….(40)The left side of the above equation represents the infiltration flux at the reaction surface and the right side represents the mass production rate.The above equation is used to obtain surface values of hydrogen and oxygen concentration.

Figure 2
Figure2Comparing the polarization curve obtained from the simulation with the results of[34]

Figure 3 Figure 4 Figure 5 Figure 6
Figure 3 Saturation value at the end of the cell (in the width of the cell) for different permeability

Table 3
Physical properties of electrodes and parameters of electrochemical reactionsCorrect selection of model parameters is an important issue for comparing model results with laboratory results.

Table 5
Electric current density in isotropic and non-isotropic GDL conditions for 0.75 V voltage

Table 6
Electric current density and saturation in isotropic and non-isotropic GDL conditions for 0.5 V voltage