Three-Dimensional Transport Modeling for Proton Exchange Membrane(PEM) Fuel Cell with Micro Parallel Flow Field.

Modeling and simulation for heat and mass transport in micro channel are beingused extensively in researches and industrial applications to gain better understanding of thefundamental processes and to optimize fuel cell designs before building a prototype forengineering application. In this study, we used a single-phase, fully three dimensionalsimulation model for PEMFC that can deal with both anode and cathode flow field forexamining the micro flow channel with electrochemical reaction. The results show thathydrogen and oxygen were solely supplied to the membrane by diffusion mechanism ratherthan convection transport, and the higher pressure drop at cathode side is thought to becaused by higher flow rate of oxygen at cathode. And it is found that the amount of water incathode channel was determined by water formation due to electrochemical reaction pluselectro-osmotic mass flux directing toward the cathode side. And it is very important tomodel the back diffusion and electro-osmotic mass flux accurately since the two flux wasclosely correlated each other and greatly influenced for determination of ionic conductivityof the membrane which directly affects the performance of fuel cell.


Introduction
PEMFC (Proton Exchange Membrane Fuel Cell) using very thin polymer membrane has been considered as a promising candidate of future power sources for transportation application and residential power. This type of fuel cell has many advantages such as high efficiency, clean, quiet, low temperature operation, and quick startup capability. However its performance and durability under harsh environments and high cost should be optimized to be competitive to conventional combustion power plants [1][2][3].
Typical schematics of a Proton Exchange Membrane Fuel Cell (PEMFC) is shown in Figure 1-(a) The cell is a sandwich of two graphite bipolar plate with micro flow channels separated by MEA(Membrane Electrode Assembly) which consists of a membrane and two electrodes with dispersed Pt catalyst. The gas diffusion layer (GDL) is porous to supply reactants to the electrodes at unexposed areas of micro flow channel [4]. As shown Figure 1-(b), the typical micro serpentine bipolar flow plates in fuel cells contain many micro flow channels that distribute the reactant gas flow over the catalytic reactive surface. The shape, size and pattern of the micro flow channels were known to greatly affect the fuel cell performance.
Modeling and simulation of the heat and mass transport in micro channel are being used extensively in researches and industrial applications to gain better understanding of the fundamental processes of fuel cell and to optimize fuel cell designs before building a prototype for engineering application.
The first application of CFD (Computational Fluid Dynamic method) to PEMFC focused on one dimensional model by Bernard and Verbrugge [5][6] and Springer [7], who presented several important principles related to the electrochemical reaction process in fuel cell. Recently, a two dimensional model emphasizing the effects of two-phase water and heat transport was developed by Fuller and Newman [8], Nguyen and White [9] and Gurau and Liu [10]. And numerical models incorporating the coupled effects of electrochemical reaction and mass transport in two or three dimensional domains have been developed by Kumar and Reddy [11], Shimpalee [12], Um and Wang [13].
In this study, we attempt to develop a fully three dimensional computational model for PEMFC which can deal with both anode and cathode micro flow fields. A commercial program FLUENT (version 6.3) was modified using User Defined Function to simulate three dimensional micro flow fields with MEA (Membrane Electrode Assembly), which consists of membrane and two electrode with dispersed Pt catalyst. The pressure drops and concentration distributions of reactants and product along the micro channel, water transport though MEA and channel, and current-voltage performance of the fuel cell were calculated through the three dimensional numerical simulation in this study.

Model assumptions
Governing equations for calculating the fully three dimensional flow field are expressed under following assumptions.
1) The gas mixture is incompressible, ideal fluid 2) The flow in the flow channel is laminar 3) Isothermal condition 4) Butler -Volmer kinetics for electrochemical reaction rate where ε is the porosity of the porous media, which is equal to unit for the gas channels, ρ the density, and u v the intrinsic fluid velocity vector, u ε v , reflects the superficial velocity in the porous media. m S denotes source terms corresponding to the consumption of hydrogen and oxygen in the anode and cathode, and the production of water in the cathode.
2) Momentum conservation equation The fluid flow in the fuel cell can be described by the general equation as: where ρ denotes the pressure, μ the effective viscous coefficient. Because the fluid flowing in the channels, gas diffusion layers and catalyst layer, membrane is different μ stand for gas viscous coefficient for gas mixture in the channel and gas diffusion layer, and liquid viscous coefficient for liquid in the catalyst layer and membrane. Furthermore, mass-weighted mixing law gives viscosity of the gaseous mixture. The source terms in the momentum equations are added based on the Darcy's law, representing an extra drag force in the equation as follows: , , 3) Species conservation equation The species conservation equation for the gas mixture is Here, k denotes chemical species that include hydrogen, oxygen, nitrogen and water. eff k D is the effective diffusion coefficient. Source term k S denotes,

Water transport equation
Water management is a critical issue for the performance of a proton electrolyte membrane fuel cell. The transport phenomena of water can be described as follows: the water molecules are transported through the polymer electrolyte membrane by the hydrogen protons and this process is called electroosmotic drag. In addition to the molecular diffusion and electro-osmotic drag, water is generated in the cathode catalyst layer due to electrochemical reaction.

1) Electro-osmotic drag flux
Electro-osmotic water flux through the membrane can be calculated from the proton flux through the membrane, given by the specified current density and Faraday law 2 ( , ) 2 : 2 where d n is Electro-osmotic drag coefficient which depends on water activity as follows.
where λ represents water contend of the membrane described as, 2 3 0.043 17.81 39.85 36.0 , 0 1; where K a , water activity, is expressed as , , ( , ) , where , w K X , sat P are water mole fraction and saturation pressure at each electrode respectively.
2) Back diffusion flux The water formation at the cathode results in a gradient in the water content between the cathode side and anode side of the membrane. For PEMFC, this gradient causes a water flux back to anode side which is superimposed to the electro-osmotic flux. This back diffusion is expressed as following water flux, 2 , , , where , m dry ρ is the dry density of electrolyte, , m dry M is the electrolyte equivalent weight, and z is the direction through the membrane thickness. w D is water diffusion coefficient which is strongly dependent on water content as follows, 3) Current density and membrane ion Conductivity

Numerical algorithm
The numerical simulation model is based on SIMPLE (semi-implicit method for pressure linked equations consistent) algorithm using Segregated solver of Fluent® (version 6.3). Uniform flow velocity and reactant concentration at the inlet are given as the initial condition. Source terms generated by the electrochemical reaction are inserted into mass and species conservation equation using the User Define Function (UDF). The coupled set of equations was solved iteratively until the relative error in each field reached a specific convergent standard (usually 9 10 − ). The flow diagram for this numerical algorithm is represented in Figure 2.

Numerical simulation model
The schematics of the parallel micro flow channel with GDL and catalyst layer is shown in Figure 3. The computational mesh structure which is divided into 40 50 24 × × (total mesh number is 48,000) is presented in Figure 4,. Boundary conditions are set as follows: constant mass flow rate at the inlet and constant pressure condition. Table 1 provides the physical parameters necessary for a numerical calculation.

Discussions
To validate the numerical simulation model used in this study, the performance curves of voltage and current density compared with the experimental data of the fuel cell with a serpentine flow channel under the same conditions in Figure 5. The computed polarization curve is in favorable agreement with the experimental polarization curve [14]. But at high current density, the discrepancy between the computational results and experimental data is a little large, and the model always over predicted the current density. In the high current density region, the low current density of the experimental results may be caused by the presence of liquid water in the catalyst layers and the gas diffusion layers. Due to the presence of liquid water, the effective porosity of the gas diffusion layers and catalyst layers decreased, and the mass transfer resistance increased. Since the current model neglected the above two-phase effect, the predicted current density at high current density was always higher than the experimental density. But the comparison of relative performance of the fuel cell with different flow configurations using single phase model may give the data necessary to design of the flow channel.
Using this numerical simulation model and conditions listed in Table 1, calculations were carried out to examine parallel micro flow fields.
It is very important to know how reactants are transported into the membrane. Figure 6 shows the velocity vectors of hydrogen and oxygen at the anode and cathode sides of flow channel at the operation condition of 0.5V. These velocity vectors indicated the mass transport at each section of the fuel cell very well. Figure 6 shows that almost all velocity vectors are directed to the membrane side in the Y direction, which meant that hydrogen and oxygen were supplied to the membrane by the diffusion mechanism. And the velocity vector pattern was symmetric pattern to the Y axis center line of flow field.  Parasite loss such as pressure drop in micro flow channel should be reduced in order for a fuel cell to produce optimum performance. Figure 7 presents the pressure distribution along the parallel flow channel at the condition of 0.5V. Pressure drop of 11 Pa at the anode and 22.5 Pa at the cathode flow channel were obtained from calculation. The higher pressure drop at the cathode side than that at the anode side is thought to be caused by the higher flow rate of oxygen at the cathode. Oxygen, water concentration, current density band activation over potential at the catalyst layer surface is presented in Figure 8 for two different cell voltages. The figure shows that oxygen, current density and over potential are high at the entrance of flow channel and then decrease along the channel due to electrochemical reaction as shown in Figure 8-(a), (c), (d). On the other hand, water concentration exhibits the gradual increase along the channel length as presented in Figure 8-(b).
The contour of current density is found to correlate with oxygen and activation over potential contour [15]. The current density is maximum on the land areas, but the oxygen concentration maxima appears at the flow channel area. And the maxima of current density on the land area was attributed to the higher ohmic loss on the flow channel area due to the longer electron path on the flow area than that on the land area.
Average value of current density is higher for low cell voltage that that for high cell voltage as shown in Figure 8-(d). In a lower voltage region, the spatial gradient of the local current density varied much more than that for the high voltage.
Oxygen It is very important to know the water transport along the cathode and anode channel since the water concentration is closely associated with the ionic conductivity of the membrane. Figure 9 shows three dimensional contour of the molar concentration of water at each channel.
In Figure 9-(a) showing the water distribution at the fully humidified anode side, the water concentration decreased along the flow channel. This decrease of water concentration along the anode channel was caused by electro-osmotic drag directing from the anode to the cathode. But water concentration at the cathode channel increased along the flow channel. This increase of water concentration was associated with the fact that the water was formed by electrochemical reaction along the channel and water was transported from anode side by electro-osmotic drag simultaneously (a) (b) Figure 9. Three dimensional distribution of moral water concentration at anode side (a) and parallel micro flow field at cathode side (b).
The oxygen concentration along the anode and cathode channels is presented in Figure 10 with water concentration. The oxygen concentration at the cathode decreased along the channel due to oxygen consumption by water formation.
Therefore, the water transport through the membrane is very crucial for simulating the water concentration and determining of water content in the membrane. Figure 11 presents the electro-osmotic mass flux and back diffusion flux across the membrane along each channel. The net flux across the membrane is expressed by the following equation It is found that the water molecules at the inlet of the anode channel are transported mostly to the cathode by electro-osmotic drag but the electro-osmotic mass flux decreased along the channel.
As shown the figure, the back diffusion took place at about 30mm away from the inlet and its amount was much smaller than the electro-osmotic mass flux. Therefore the net water mass across the membrane is directed from anode to cathode side.
The current density on the catalyst layer is presented in Figure 12. As expected, the current density at 1.0 V in Figure 12-(a) is much smaller than that at 0.3V. The current density at the inlet was the highest and decreased along the channel. The highest value of current density at inlet is probably because of the high concentration of hydrogen and oxygen and high electro-osmotic mass flux at the inlet region. Figure 12-(e) presents the current density (i) -voltage(V) polarization curve for the fuel cell with parallel micro flow channel. The current density increased with decreasing cell voltage. The voltage for the maximum power density as found to be 0.5V was optimal operation condition for fuel cell in this study.

Conclusions
Using a single-phase, fully three dimensional simulation model for PEMFC which can deal with both anode and cathode flow, the following conclusions for PEMFC with parallel micro flow channel were obtained. 1) Hydrogen and oxygen were supplied to the membrane by diffusion mechanism rather than convection transport. Pressure drop was 11 Pa at anode and 22.5 Pa at the cathode flow channel for 0.5V, and a higher pressure drop at the cathode side is thought to be caused by the higher flow rate of oxygen at the cathode.
2) The amount of water in cathode channel was determined by the water formation due to the electrochemical reaction plus the electro-osmotic mass flux directing toward the cathode side.
3) It is very important to model the back diffusion and electro-osmotic mass flux accurately since the two flux were closely correlated with each other and have greatly influenced the determination of ionic conductivity of the membrane, which directly affected the performance of the fuel cell.
4) Since our model was only a single-phase model, it has some limitations in representing water flooding phenomena properly. However, a comparison of the relative performances of fuel cells with different flow configurations may give the data necessary for designing the micro flow channel.