Outlines

The scope of this dissertation is mainly focused on four parts. Chapter 2

introduces the homogeneous model of the multiphase mixture formulation (M² model)

development as it is particularly suitable for two-phase flow modeling and applying

computational fluid dynamics (CFD) to solve the complete set of transport equations

governing mass, momentum, species, energy, and charge conservation in a PEM fuel

cell. Furthermore, the proposed model incorporates the basic assumptions is also

presented. Chapter 3 investigates the effects of cathode humidification and gas

diffusion layer porosity on the location of the interface when the liquid water began to

condense along the flow channel of a PEM fuel cell; its contribution to the cell

performance is also discussed. Additionally, the resulting oxygen and water fraction

distributions and liquid water saturation fields at fixed cathode humidity are obtained

to validate the simulation results. Chapter 4 reports the formation and influence of the

gas-liquid interface location along the flow channel direction at various cell

temperatures and humidification temperatures with a non-isothermal,

three-dimensional, multi-component, and two-phase model of a PEM fuel cell. The

temperature distribution in the cell domain membrane and the distributions of

temperature and liquid water saturation in the cross-section of the cathode gas

diffusion layer in the inlet region are investigated. In Chapter 5, the conclusions of

this investigation are drawn and suggestions for future study are proposed.

Figure 1.1 Basic components of a single PEMFC End plate

Current collector

Bipolar plate Diffusion media Membrane

Current collector

End plate

Current collector

Bipolar plate Diffusion layer

Catalyst layer

Membrane

Diffusion layer

Bipolar plate

Current collector End plate End plate

Figure 1.2 Operating principle of a single PEMFC Load

e

^-

H

⁺

H

2

e

-e

-Anode Cathode

O

₂

H

₂

Electron Proton

O

2

H

₂

O

Figure 1.3 The dominant mechanisms of cell performance Open-circuit potential

Electrode Kinetic Dominate

Ohmic Polarization Dominate

Mass Transport Dominate

Cell Voltage

Current Density

CHAPTER 2

MATHEMATICAL MODELING

Traditionally, macroscopic problems of two-phase flow and transport in

porous media have been modeled using a two-fluid approach. However, this approach

results in a large number of primary variables for each phase, and highly nonlinear

equations. Therefore, exact solutions of two-phase problems with two-fluid models

are limited to a very limited number of problems with many simplifying assumptions.

Furthermore, the two-fluid models require explicitly tracking the irregular and moving

interface between two phases, increasing the numerical complexity of the problem.

Particularly in PEM fuel cells, the gas-liquid interfaces, i.e., the condensation and

evaporation fronts, are expected as well as the coexistence of single- and two-phase

regions. Therefore, a convenient model capable of describing both single- and

two-phase regions without a need to track the irregular, a priori unknown interface is

required. For these reasons, the multiphase mixture (M²) formulation of Wang and

Cheng [12] is particularly suitable for two-phase PEM fuel cell modeling. The

multiphase mixture model is a mathematical reformulation of the classical two-phase

model that views the multiphase system as a chemical mixture. With this approach,

the multiphase flow is then described in terms of a mass-averaged mixture velocity

and diffusive flux, representing the difference between the mixture velocity and

individual phase velocity. One major advantage of the M² model over the classical

two-fluid models is that it eliminates the need for tracking phase interfaces, thus

simplifying the numerical complexity of two-phase flow and transport modeling.

Another salient feature of the M² model for PEM fuel cell is that all model equations

are valid in all three types of regions possible in a PEM fuel cell: single-phase (gas),

liquid-gas (two-phase), and single-phase (liquid). Finally, the M² model is

mathematically equivalent to two-fluid models without invoking any additional

approximations. These aforementioned advantages render the M² model to be a

suitable and widely adopted two-phase flow and transport modeling framework for

PEM fuel cells [11, 13, 15].

In this chapter, a three-dimensional, two-phase, multi-component model of a

PEM fuel cell is presented and the present analytical study focus on two-phase

transport taking place in the cathode due to the production of water in the cathode

catalyst layer.

2.1 Model Description

The physical model employed in this study consists of nine parts, namely

bipolar plates, gas flow channels, gas diffusion layers, and catalyst layers on both the

anode and cathode sides, and a polymer membrane that is sandwiched between them,

as illustrated in Fig. 2.1. Symmetry is assumed, and a single straight flow channel is

therefore considered herein.

2.2 Basic Assumption

Utilizing the homogeneous model of M² formulation for two-phase transport,

the proposed model incorporates the following assumptions. These basic assumptions

are made to simplify actual cell conditions in the theoretical model and thus facilitate

the modeling approach of transport component influence on transport phenomena and

cell performance.

z The gaseous phase of the working fluid behaves as an ideal gas and the liquid

water is incompressible;

z The Reynolds number of the fluid is below 100 because the velocity of the

mixture is low, and the flow is considered to be laminar;

z The properties of the porous medium are isotropic and homogeneous.

z The system operates in a steady state;

z The thermal conductivities of the cell components are constant.

2.3 Governing Equations

The three-dimensional, multi-component, and two-phase model of a PEM fuel

cell includes five nonlinear coupled conservation equations of mass, momentum,

energy, species, and charge, which are described as follows.

Mass conservation equation

Since the gaseous and liquid water are present simultaneously in the control

volume, the conservation equation of mass for a multiphase mixture is:

( )

⁼⁰

⋅

∇ ε_effρuK (2.1)

where ρ represents the density of the mixture and is defined as the volume-weighted

average of the phase mass concentration in two-phase flow [11]. When the liquid

water is present, the effective porosity is given byε_eff = 1ε

(

−s

)

. Momentum conservation equation

The general form of the Navier–Stokes equation is used with source terms that

describe the drag forces in the porous medium. The equation is:

(

effuu

)

= − eff∇P+∇⋅

(

eff eff∇u

)

+ eff kg+Sm

⋅

∇

ρε

^KK

ε ε μ

^K

ε ρ

(2.2)

where μ_eff is the effective viscosity of the mixture, ρ_k is the kinetic density, and

Sm is the sum of the Darcy and the Forchheimer drag forces,

)

( _{Darcy Forch}

m F F

S G G

+

−

= (2.3)

which are given by:

Darcy drag force =

K F_Darcy ^{eff eff}u

G ε ²μ G

= (2.4)

Forchheimer drag force = uu

K

FG_Forch ε^eff3C^Fρ G G

= _(2.5)

The parameters CF andK represent the quadratic drag factor and the permeability.

Species conservation equation

The species conservation equation for the gas mixture is

(

effuCk

)

=∇⋅

(

Dkeff∇Ck

)

+Sc

⋅

∇

ε

G _, (2.6)

where k represents the chemical species, including hydrogen, oxygen, nitrogen, and

water. Dk,eff =D_k

ε

^τrepresents the effective diffusion coefficient of the k-th component of the fuel reactant [17]. The exponent τon the porosity ε is the tortuosity of the

porous medium. The source term Sc defines the production or consumption of the k-th

species in the gas phase and is given by:

a

( )

The above conservation equations of mass, momentum, and species are derived on

the basis of the M² model. The constitutive relationships of mixture parameter and

variables are all dependent on the liquid saturation, defined as the ratio of the liquid

volume to the pore volume:

Effective viscosity:

( )

(

rl l

) (

rg g

)

Relative permeability:

( )

^phase

phase

Relative mobility:

( )

( ) ( )

^phase

During the operation of a fuel cell, liquid water forms if the partial pressure of

water vapor exceeds the saturation vapor pressure. The liquid water thus formed may

occupy the pores and thereby prevent the diffusion of fuel, causing mass transport

overpotential in the porous medium. Hence, the effect of liquid water is taken into

account. Additionally, capillary forces dominate the transport of liquid water on the

hydrophilic surfaces because the pores in the porous medium are extremely small.

Therefore, the generalized Richards equation, developed by Wang and Beckermann

[12, 40] to elucidate the two-phase flow transport in capillary porous media, is

applied:

whereD^candP^crepresent the capillary diffusion coefficient and the capillary pressure,

respectively:

and where J(s) is the Leverett function, which takes the following form [41, 42]:

( ) ( ) ( ) ( )

In Eq. (2.22), the contact angle, θc, of the gas diffusion layer is dependent upon the

hydrophilic (0° < θc < 90°) or hydrophobic (90° < θc < 180°) nature of this layer, and

varies with the Teflon content. We assume here that the gas diffusion layer is a

hydrophilic medium. Further, the surface tension, ζ, for the liquid water-air system is

taken as 0.0625 N/m [41].

The source term S_l is a simplified switch function between condensation

and/or evaporation of liquid water under these non-equilibrium conditions [15]. When

the partial pressure of water vapor exceeds the saturation pressure of water, liquid

water may form and occupy the pores in the porous medium. Conversely, the liquid

water will evaporate if the partial pressure of water vapor is less than the saturation

pressure of water:

where rcon and reva are the condensation and evaporation rate constants, respectively;

x is the molar fraction of water vapor, and w P_satis the saturation pressure of water, which varies with the temperature [4].

3

The heat generation sources in a PEM fuel cell account for the irreversible

heat and entropic heat that is generated by electrochemical reactions, Joule heating

that arises from proton/electronic resistance, and the latent heat of water condensation

and/or evaporation. A generalized energy conservation equation is:

( )

T

where k represents the effective thermal conductivity. On the right-hand side of _eff

the equation, the first two terms represent the conduction energy and the reactant

enthalpy flux; the third and fourth terms represent electrical-related thermal effects,

and the last term is a source term, associated with the phase change.

Charge conservation equation

In a fuel cell, the potential gradient effect causes electrons and protons to

move along individual paths. Solid phase potential controls the movement of electrons.

Electron transport generally occurs only in the bipolar plates, the diffusion layers, and

the catalyst layers. However, ionomer phase potential controls the motion of protons,

which occurs in the catalyst layer and the membrane. Potential fields in these two

media are described as follows.

For electrons:

(

^∇

)

^{+ =0}

S_φp denote electron conductivity, proton conductivity, electronic phase potential, electrolyte phase potential, and consumption

rates of charge and product in the electrochemical reaction in the catalyst layer,

respectively.

The membrane conductivity is strongly related to the temperature and the

water content λ . It is defined as the ratio of the number of water molecules to the

number of charge sites [4]:

( )

_⎥

The water content of the membrane surface depends on the activity of the water vapor,

which also depends on the partial pressure of the water. Therefore, the empirical

relationship between them can be applied:

At anode At cathode

At anode At cathode

⎩⎨

2.4 Boundary Conditions

Boundary conditions are necessary and crucial for solving the above equations.

They describe the operating conditions as well as the model geometry characteristic of

the PEM fuel cell. Due to the single-domain formulation, boundary conditions are

required only at the external surfaces of the computational domain. The most

important ones are as follows.

Inlet boundaries: The fuel and oxidant flow rates along the flow channel can be

described by a stoichiometric flow ratio, ξ, which is defined as the ratio of the

amount of reactant supplied to the amount of reaction to generate the specified

reference current density Iref. [37]. At the gas channel, the temperature and gas species

concentrations are assumed to be uniform. The inlet velocities are specified by

ch

Where ξ_aandξ_care the reactant stoichiometric flow ratio of anode and cathode,

respectively. Am is the geometrical area of the membrane and Ach is the cross-sectional

area of the gas channel.

Outlet boundaries: Fully developed flow is applied. At the outlets, both anode and

cathode channels assumed sufficiently long so that velocity and species concentration

fields are fully developed.

Walls:Neumann conditions and no-slip conditions are applied.

Symmetric boundaries:Mass flux or momentum flux have zero gradients.

Electronic phase potential boundaries: Fixed total cell overpotential at the outer

boundary of the cathode is specified.

plate

The cell potential can be obtained from the following expression

tot l oc

ce V

V = +η (2-37)

where η_tot andV_ocare total cell overpotential and pen circuit voltage, respectively

[17].

2.5 Numerical Procedures

2.5.1 Numerical Method

The solution to the governing equations is performed by employing a finite

volume scheme with the model domain divided into a number of cells as control

volumes. The governing equations are numerically integrated over each of these

computational cells or control volumes. The method exploits a collocated

cell-centered variable arrangement with the local or cell-averaged values of the

physical quantities evaluated and stored at each cell center.

The governing equations can be expressed in the form of a generalized

convection-diffusion type of transport equation:

(

ρuϕ−Γϕ∇ϕ

)

=Sϕ

⋅

∇ (2-39)

where ϕ denotes the general dependent variable, Γ the exchange coefficient, _ϕ S _ϕ

the source term, uG

velocity vector, and ρ the density. With the discretization of

the governing equations, the coupled finite-difference equations can be expressed in

the form of

ϕ φ

ϕ ϕ

ϕ

ϕ a a a a S

a_{P P} = _{E E} + _{W W} + _{N N} + _{S S} + (2-40)

where ϕ is the value of ϕ at the current point P, _P ϕ …_E ϕ stand for the values of _S

the grid points adjacent to the point P, and a_P…a_S are known as the link

coefficients.

2.5.2 Calculation Procedure

The governing equations with their related boundary conditions are solved

using a commercial code based on the SIMPLE algorithm for convection-diffusion

problems. The numerical flow diagram of present investigation is shown in Fig. 2.2.

As a convergence criterion it is imposed that the normalized residual for each model

variable is smaller than 10⁻⁴ [43].

2.5.3 Model Validation

In the simulation, a uniform grids distribution is used to calculate the complex

electrochemical reaction and physical phenomena in the fuel cell. Three mesh

systems- 41 x 13 x 47, 51 x 16 x 58, and 61 x 21 x 67 are constructed to explore

numerical result dependence on computational cell numbers. Table 2.1 presents

geometrical and operating parameters of the base model in the PEM fuel cell. The

results of the polarization curve by the base model under different grid systems are

shown in Fig. 2.3. Considering both accuracy and economics, the grid system of 51 in

the z-direction, 16 in the x-direction, and 58 in the y-direction was selected for present

research.

To further check the adequacy of the numerical scheme, it is clearly seen from

Fig. 2.4 that the present predictions agree reasonably with the experimental data of

Squadrito et al. [44]. The above preliminary runs confirm that the present model and

the numerical method used are generally appropriate in analysis of the present

problem.

3 [14]

Tortuosity of the diffusion and catalyst layers

0.5 / 1.5 [14] Relative humidity of anode/cathode inlet

1.5 / 3 [44]

Stoichiometry, at / at 1.0

0.28 [45]

Porosity of membrane

0.4 / 0.4 [45]

Porosity of diffusion and catalyst layers

[14]

Permeability of membrane

/ [14]

Permeability of diffusion and catalyst layers

3 [14]

Tortuosity of the diffusion and catalyst layers

0.5 / 1.5 [14] Relative humidity of anode/cathode inlet

1.5 / 3 [44]

Stoichiometry, at / at 1.0

0.28 [45]

Porosity of membrane

0.4 / 0.4 [45]

Porosity of diffusion and catalyst layers

[14]

Permeability of membrane

/ [14]

Permeability of diffusion and catalyst layers

Table 2.1. Geometrical and operating parameters

Quantity Value

Cathode BP Anode BP

Cathode Channel

Anode Channel

Cathode GDL Cathode CL

Membrane Anode CL Anode GDL

X Y

Z

Cathode BP Anode BP

Cathode Channel

Anode Channel

Cathode GDL Cathode CL

Membrane Anode CL Anode GDL

X Y

Z

Figure 2.1 Physical and computational domains considered in this study

Begin

Solve momentum eq. obtain u*, v*, w*

Obtain pressure and velocity correction Guess P*

Set P*=P Stop

Solve for scalar equations of species, enthalpy, potential field Obtain mass imbalances

Correct pressure and velocity fields

Converged ?

Solve momentum eq. obtain u*, v*, w*

Obtain pressure and velocity correction Guess P*

Set P*=P Stop

Stop

Solve for scalar equations of species, enthalpy, potential field Obtain mass imbalances

Correct pressure and velocity fields

Converged ?

Figure 2.2 Numerical flow diagram of the solution procedure.

0 0.2 0.4 0.6 0.8 1 1.2

0 0.2 0.4 0.6 0.8 1 1.2

41 x 13 x 47 51 x 16 x 58 61 x 21 x 67

Vo lt a g e ( V )

Current Density ( A/cm

²

)

Figure 2.3 Comparison of predictions on the three different grid systems.

0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8

Squadrito et al. [44]

Present result

Vo lt a ge ( V )

Current Density ( A/cm

²

)

Figure 2.4 Comparison of the predicted I–V curve and the experimental data of Squadrito et al. [44].

CHAPTER 3

EFFECTS OF CATHODE HUMIDIFICATION AND POROSITY OF THE GDL ON THE GAS–LIQUID INTERFACE LOCATION IN A PEM FUEL CELL

3.1 Introduction

A PEM fuel cell is prone to gas-liquid two-phase formation due to its low

operating temperature, particularly under highly humidified or high current density

conditions. When the gas diffusion layer and the catalyst layer become saturated with

water vapor, the product water starts to condense and block open pores, reducing the

available paths for oxygen transport. This phenomenon is termed “flooding” and

becomes a major limiting factor of PEM fuel cell performance. Hence, it is critical to

understand the two-phase flow and transport in a PEM fuel cell, and a mathematical

model is useful to improve this understanding. In practice, humidification of anode

fuels and/or cathode oxidants is often used to provide sufficient membrane hydration.

Since water is generated in the cathode catalyst layer from electrochemical reaction

and it also tends to migrate from the anode side to the cathode side by the

electro-osmotic drag, it becomes a key issue in the design and operation of PEM fuel

cells to avoid the flooding phenomenon in the cathode. In this chapter, effects of

cathode humidification and gas diffusion layer porosity on the gas-liquid interface

location in a PEM fuel cell are investigated and described in detail. For all of the

calculations carried out in this study, the fuel flows were co-flows and inlet

stoichiometric ratios of 1.5 and 3 are used for the anode and cathode sides based on a

reference current density of 1 A/cm². The inlet fuel, hydrogen, is assumed to be fully

humidified on the anode-side; the cell temperature and humidification temperature are

343 K and the operating pressure is 1 atm.

3.2 Effects of cathode humidification scheme

In this study, the location of the gas-liquid interface along the flow channel

direction at various cathode humidity conditions and its effect on cell performance are

elucidated by modeling a three-dimensional PEM fuel cell system using CFD. Figure

3.1 shows the effect of the relative humidity of the cathode on the location of the

interface where the liquid water begins to condense along flow channel at a cell

operating voltage of 0.7 V. The interface is defined as the location where liquid water

begins to condense. The horizontal dotted line indicates the interface between the flow

channel and the gas diffusion layer. A higher cathode relative humidity corresponds to

a smaller distance between the gas-liquid interface and the gas inlet, as clearly

displayed in Fig. 3.1 (a). Furthermore, the gas-liquid interface location moves to the

flow channel inlet region as the relative humidity of the cathode increases. This

phenomenon is caused by: (i) the decrease in the amount of evaporated water through

the flowing gas stream; (ii) the increase in the partial pressure of water and the ability

to reach the saturation pressure of the vapor water relatively quickly to form liquid

water earlier as the relative humidity of the cathode increases. On the contrary, the

gas-liquid interface location moves closer to the cathode catalyst layer when the

relative humidity of the cathode is less than 60%, as shown in Fig. 3.1 (b). In

conclusion, increasing the relative humidity of the cathode can reposition the

gas-liquid interface and cause liquid water to appear, affecting the performance of the

cell, because liquid water may occupy the pores in the porous media, reducing the

cell, because liquid water may occupy the pores in the porous media, reducing the

在文檔中 質子交換膜燃料電池中陰極層汽液界面位置對性能之影響 (頁 35-0)