T Sousa

M Mamlouk

K Scott
email: k.scott@ncl.ac.uk

Intermediate Temperature Fuel Cell using PBI doped Phosphoric Acid Membranes

Keywords: Numerical Simulation, PEMFC Model, Intermediate Temperature Fuel Cell, Non-isothermal, Polybenzimidazole Polybenzimidazole, Numerical Simulation, Intermediate Temperature Fuel Cell, PEMFC Model, Non-isothermal

A two-dimensional non-isothermal model developed for a single Intermediate Temperature Fuel Cell with a phosphoric acid doped PBI membrane is developed. The model of the experimental cell incorporates the external heaters, and the all body of the fuel cell.

The catalyst layers were treated as spherical catalyst particles agglomerates with porous interagglomerate space. The inter-agglomerate space is filled with a mixture of electrolyte (hot phosphoric acid) and PTFE. All the major transport phenomena are taken into account except the crossover of species through the membrane. This model was used to study the influence of two different geometries (along the channel direction and cross the channel direction) on performance. It became clear, through the performance analyses, that the predictions obtained by the along the channel geometry did not represent the general performance trend, and therefore this geometry is not appropriate for fuel cell simulations. Results also indicate that the catalyst layer was not efficiently used, which leads to large temperatures differences through the MEA.

Introduction

Conventional PEMFCs have typically used Nafion ® or similar perfluorinated polymer as membrane. The proton conductivity of these membranes strongly depends on its water content. This feature leads to temperature restrictions during fuel cell operation, i.e. it can only operate efficiently below the water boiling point, otherwise the membrane dehydrates and its conductivity falls drastically. In order to maintain the correct water balance the inlet streams must be humidified. As a consequence of the low operating temperature (typically 80 °C), fuel cells with Nafion ® membranes exhibit three main problems: inadequate water and heat management, sluggish electrochemical kinetics, and intolerance to impurities, mainly carbon monoxide (CO). To overcome these problems much effort has been made to develop polymer membranes capable of retaining high proton conductivity at higher temperatures without humidification. Poly[2,2-(m-phenylene)-5,5 bibenzimidazole doped with phosphoric acid, generically referred to as polybenzimidazole (PBI), is one membrane that has been proposed in recent years. This kind of membrane has been considered as a promising candidate for use in intermediate temperature PEMFCs [1][2][3][4][5], as although other kinds of polymer and hybrid membranes are under intensive research [6][7][8][9][10][11][12][13][14][15][16][17].

Modelling plays an important role in fuel cell technology because it increases understanding of its behaviour, enables prediction of its performance, and assists with its operational control. Bernardi and Verbrugge [18] were the first to develop a model for a PEMFC running with a Nafion ® membrane. This model was one-dimensional (1D), isothermal and for steady state conditions. It was applied only at the oxygen electrode in order to study the polarization characteristics, water transport, and catalyst utilization. A significant amount of work has been produced since this first model and reviews of the subject have been made [19][20][21].

In comparison with Nafion ® based PEMFCs, only a few publications have considered models for PBI based PEMFCs. Korsgaard et al. [22] developed a simple semi-empirical model to correlate their experimental data. The model showed good agreement with the experimental data, although only predicted the polarisation curves based on the cathode activation overpotential, ohmic and concentration losses. Cheddie and Munroe. [23][24][25][26] developed several papers on this particular issue. In all of them, heat generation due to the electrochemical reaction, and ohmic heating were taken into account. The most recent work [26] was an improvement of the previous models, which accounted for two phase effects due to the transport of dissolved species in phosphoric acid/PBI catalyst layer, and considered aqueous phase electrochemical reactions. This model showed a good fit with the (1.0 mA cm -2 ) was applied by means of 9.0 V battery connected in series with an appropriate resistance.

The temperature of the cell was controlled by thermostatically controlled tubular heaters inserted into the cell bodies. The relative humidity was obtained from an intrinsically safe humidity sensor (Vaisala HUMICAP ® , Finland).

The membranes were prepared from PBI powder dissolved in N,N-dimethylacetamide (DMAc) by casting onto optical glass and kept in an oven at a temperature of 90 to 110 ºC overnight to produce a 40 µm thick membrane. The catalyst ink was prepared either by sonicating the catalyst and PBI in acetone/DMAc, for PBI based MEAs, or the catalyst and PTFE dispersion (60%wt, Aldrich) in water-ethanol mixture, for PTFE based MEAs. Gas diffusion electrodes (carbon cloth) incorporated with wet proofed micro porous layer (H2315 T10AC1) obtained from Freudenberg (FFCCT, Germany) were used as substrates to deposit the catalyst layer for both anode and cathode. Anodes were prepared using 0.2 mg Pt cm -2 20% Pt/C (ETEK) with PBI loading of 0.28 and 0.7 mg cm -2 and fixed doping level of 8 PRU.

MEAs were prepared by hot pressing the GDL supported catalysts onto the PBI membrane at 150 °C and 40 kg cm -2 for 10 min.

For the fuel cell simulations we were particular interested in analysing the effect of different computational domains on model predictions. Therefore, two different computational domains were employed for the simulations (Geometry A and Geometry B).

As can be seen in Figure 1,. Geometry A consisted of a yz-plane which took into account the complete fuel cell body (including: titanium body, PTEF isolation layer, heaters, gas channels and MEA). With this computational domain it was possible to analyse variations along-the-channel (y) and through the membrane (z) due to changes in concentration and pressure. Figure 2 shows details of this computational domain, where it is possible to see with this geometry only two gas channels were modelled, one in the anode side and another in cathode side. Geometry B, also took into account the entire fuel cell body but in a different direction (in this case the xz orientation, as Figure 3 shows). Therefore, this computational domain models an array of channels and considered variations across-the-channel (x) and through the membrane (z) which made possible the analysis of the effect of ribs on performance. With this computational domain 24 channels were modelled (12 in the cathode side and another 12 in the anode side) and the effect of 26 ribs could be analysed. For both computational domains the MEA was divided into 7 different sub-domains: anode and cathode gas diffusion layers (GDL), anode and cathode micro-porous layers (MPL), anode and cathode catalyst layers and the membrane. All the relevant geometric parameters and operational properties are shown in Table 1.

In this model a number of assumptions were adopted:

(i) Steady state operation.

(ii) Single phase flow.

(iii) Ideal gas mixtures.

(iv) Porous media was assumed to be isotropic and macro-homogeneous.

(v) The PBI membrane was assumed to be impermeable to gas flow.

(vi) Negligible contact resistance between components.

(vii) Heat transfer through radiation was neglected.

Governing Equations

This section provides the governing equations which describe the physical, chemical, electrical and electrochemical processes occurring in the PEMFC described above.

Conservation of Momentum

At the gas channels, the gas flow field and pressure was obtained by solving the equation of continuity (Eq. ( 1)) and the weakly compressible form of the Navier Stokes equation (Eq.

(2)).

( ) 0 = ⋅ ∇ v ρ (1)
( )

Γ ⋅ ∇ + -∇ = ∇ ⋅ µ ρ p v v (2)
Where, ρ is the density, v is the velocity vector, p is the fluid pressure, µ is the fluid viscosity and Γ is the deviatoric strain-rate, which is defined by: ( ) ( )

v v v ⋅ ∇ - ∇ + ∇ = Γ 3 2 T ( 3
)
In the porous media Eq. ( 2) had to be corrected for porosity. Therefore, Brinkman's equation (which is a generalization of Darcy's law that facilitates the matching of boundary conditions between porous media and free media [30] ) was used [31]:

Field Code Changed
Field Code Changed Deleted: The intersection of both geometries is a line parallel to the z-axis exactly in the centre of the cell. Geometry A is especially useful to analyse changes along the channel due to changes of concentration and pressure drops, on the other hand, Geometry B will be focus in the effect of the ribs between channels. ¶ Several sub-domains were defined within both geometries in order to successfully implement the numeric solution of the various equations. As can be seen in Figure 2, f

     Γ ⋅ ∇ + -∇ = ε µ µ eff eff p K v (4)
Where, K is the permeability of the porous media, ε is the porosity of the porous media, and eff µ is the effective viscosity

Conservation of Species

The mass conservation for the individual gas phase species was obtained by solving the Stefan-Maxwell equations [27]:

0 ~= + ⋅ ∇ -                 ∇         - + ∇ ⋅ ∇ ∑ s i j j j j j eff ij i St w p p M M x x x D w v ρ ρ (5)
In Eq. ( 5), i w is the mass fraction of the species i , j x is the mole fraction of species j , j M is the molecular weight of species j , M is the average molecular weight, s St is the source term for the conservation of species, and

eff ij D ~
is the effective binary diffusivity of the pair i , j .

The source term accounts for the local rates of oxygen and hydrogen consumption and water vapour production. Therefore, the source term only assumes a value different to zero in the catalyst layers, because the electrochemical reactions only occur in these layers.

F M j St H a H s 2 2 2 , = (6)
F M j St O c O s 4 2 2 , - = (7)
F M j St O H c O H s 2 2 2 , = (8)
Where a j is the anode volumetric current density, c j is the cathode volumetric current density, M is the molecular weight for the different species, and F is the Faraday constant.

Deleted:

Deleted: [25] Deleted:

i n

Kinetics and Thermodynamics

The Butler-Volmer equation was used to describe the kinetics at the anode and cathode.

            -       - = η α η α RT F RT F i a j Ox Rd Pt exp exp 0 (11)
Where, Pt a is the electrocatalytic surface area, α is the transfer coefficient, R is the universal gas constant, T is the temperature, 0 i is the exchange current density which is corrected for the particular conditions of pressure and temperature of the system., and η is the overpotential.

The catalyst layers were treated as spherical agglomerates of catalyst particles with porous inter-agglomerate space, where each agglomerate is covered by a thin film of electrolyte.

Therefore, Eq. ( 11) was modified to account for the geometry of the agglomerate and the transport of the dissolved species through it.

( )

( ) 1 , / 1 1 -         + + - = i liq agg agg agg CAT c r l g i D a r r k E nFC j δ δ ε (12)
Where, c k is the reaction rate constant, l g i C / is the concentration of the reactants at the gas/liquid interface, Er is the effectiveness factor, agg a is the effective agglomerate surface area, agg r is the agglomerate radius, i liq D , is the reactants diffusion in phosphoric acid, δ is the thin film thickness. The reaction rate was obtained by the following equation:

( )             -       - - = η α η α ε RT F RT F C nF i k Ox Rd s l i CAT c exp exp 1 0 (13)
A detailed explanation about the above approach is given in Appendix A The overpotential ,η , is the electrochemical reaction driving force and is defined by the potential difference between the solid and electrolyte phase :

ref m s V - - = φ φ η (14)
The reference potential, ref V , on the anode side is zero and is equal to the theoretical cell potential at a given temperature and pressure on the cathode given by Nernst equation.

        +       ∆ - ∆ - = O H O H ref a C C RT nF RT nF S T nF H V 2 2 2 5 . 0 5 . 1 ) ( ln (15)
Where, H ∆ is the enthalpy change, S ∆ is the entropy change,

Conservation of Energy

The conservation of energy for any sub-domain was described as follows:

( )

e p T S T c T k + ∇ ⋅ = ∇ - ⋅ ∇ v ρ (16)
Where, T k is the thermal conductivity, p c is the species heat capacity, and e S is the energy source term.

The energy equation accounted for conduction, convection and heat generation due to Joule heating, as well as the heat of reaction. Joule heating, Joule q , was taken into account in the GDLs, MPL, catalyst layers, and membrane. In the collector plates Joule heating was assumed to be insignificant due to their high conductivity:

k q Joule 2 . φ = (17)
The fuel cell electrochemical reactions are exothermic and heat is generated during the fuel cell operation time. This heat generation is due to entropy changes and irreversibilities An exhaustive study about this issue was presented in [32]. The heat of reaction, . reaction q was calculated from:

c reaction j nF S T q       - ∆ = η . ( 18
)
The heat of reaction is directly proportional to the entropy changes and overpotential.

Because both of these contributions are very small at the anode side it was assumed that the heat of reaction was negligible at the anode side. It is well known that the overpotential at the anode side is very small compared to that at the cathode. From the data provided by Lampinen and Fomino [32] the half cell entropy change at the anode side and at the cathode side at 298.15 K were calculated as: Entropy change at the anode side: -0.104 J mol -1 K -1 Entropy change at the cathode side:-44.3 J mol -1 K -1 Thus the entropy change at the anode side can be ignored when compared to the entropy change at the anode side.

Boundary Conditions

In order to solve the governing equations, boundary conditions were applied for all variables of interest. As mentioned before this model was used in two different geometries; this fact leads to different boundary conditions for each geometry.

Geometry A

Momentum balance

The inlet velocity was set for both the anode and cathode sides. The inlet velocity used in the experimental work carried out in our laboratories was 0.1768 m s -1 and 0.0589 m s -1 for the cathode and anode, respectively. From these values the inlet velocities for each channel could be estimated.

  = = = 0 0 y a x v v v v ( 19
)      = - = = 0 0 y c x v v v v (20)
Where, a v 0 is the anode inlet velocity and c v 0 is the cathode inlet velocity.

In the gas channels, outlets pressure was set for both anode and cathode as follows:

0 p p = (21)
Where 0 p is operational pressure.

Species balance

The inlet mass fraction was set for both the anode and cathode sides. At the anode non humidified hydrogen was used and air was fed to the cathode.

1 2 = H w (22) 011 . 0 23 . 0 2 2 = ∧ = O H O w w ( 23
)
At the gas channels outlets, the convective flux sets the boundary condition:

n ⋅ = v w n H H ρ 2 2 (24) n n ⋅ = ∧ ⋅ = v w n v w n O H O H O O ρ ρ 2 2 2 2 (25)

Energy balance

In the experimental apparatus there was no thermocouple to measure the temperature of the gas feed at the gas channel inlets. As a consequence the inlet temperature was assumed to be equal to the fuel cell operating temperature.

0 T T = (26)
Where, 0 T is the operational temperature.

At the outlets a convective flux as boundary condition was applied. At the exterior boundaries of the domain it was simply required to solve the conservation of energy, therefore, there we only required boundary conditions for this equation. It was assumed that a uniform surface heat flux existed for all exterior walls of the cell. It was also assumed that all exterior walls were subjected to the same exterior air velocity and temperature. To calculate the heat flux, conv q , the following equation was used:

( ) 0 = ∇ - T k T n (27)
( )

air conv T T h q - = (28)
Where, h was the average heat transfer coefficient, and air T was the average temperature of the exterior air.

The average heat transfer coefficient was calculated by the following equation:

L k Nu h Air T = (29)
Where, Nu is the average Nusselt number, , with a uniform surface heat flux, and assuming that the fuel cell can be represented by a rectangular section it can be calculated as follows [33]:

3 1 675 . 0 Pr Re 102 . 0 = Nu ( 30
)
Where, Re and Pr are the Reynolds and Prandtl numbers, respectively.

To simulate the heaters the temperature at the heaters/collector plate boundaries were prescribed to the desired operating temperature.

Charge balance

For the conservation of charge it was necessary to set boundary conditions for both solid and electrolyte phase potentials. The solid phase potential was prescribed at the rib/GDL interface. For the anode side it was set to zero volts, and to the cell potential at the cathode rib/GDL interface. Where, Cell V is the operational potential of the fuel cell.

For the electrolyte phase the potential was only required to prescribe electrical insulation to the exterior boundaries formed by the membrane and catalyst layers.

Geometry B

The boundary conditions for charge and energy balances were the same of Geometry A; however, there were major differences for momentum and species balance. From Fig. 3 it is possible to see that Geometry B did not include the gas channel inlets or outlets, therefore a different approach was applied to solve the momentum and species balance. It was assumed that the species were transported only through diffusion in the z and y directions, so the momentum balance was no longer needed. The species mass fraction was prescribed at the most interior wall of the gas channels. In order to increase the model accuracy these mass fractions should vary with potential. Figure 4 shows the variation of the average species mass fraction with potential at the gas channel\collector plate interface (these data was obtained with a 2D isothermal, along the channel model developed previously and described in [34]).

These data were used to generate a set of equations which were utilized to predict the variation of the species mass fraction with potential. Most of the materials and species properties were treated as variables which depended on temperature, pressure and concentration. The relationship between these properties and the dependent variables are well documented in the literature; therefore this issue will not be subject of attention in this paper. Table 2 show the key parameters that were used for the simulations.

V V V V V V V V w ( 33
)      < + × + + - ≥ + + - = - V 7 . 0 ,
V V V V V V V V w ( 34
)      < × + × - - ≥ + - + - = - - V 7 . 0 ,
V V V V V V V V w ( 35
The species transport through the agglomerates was extremely important in this system.

The diffusion of gases through liquid media had a big impact on performance; in particular at high current densities due to the relative low diffusion coefficient. At high current densities the reaction rate was mainly limited by the sluggish transport of oxygen in the catalyst layer.

Both the Henry's constant and diffusion coefficient for oxygen in phosphoric acid depend on temperature and phosphoric acid concentration. From the experimental data provided by Klinedinst et al. [35], a function to calculate the Henry's constant, , for oxygen in hot phosphoric acid was built: is the phosphoric acid mole fraction.

( ) ( )       - + × + × - = 6
( ) ( )           × + × + × - × + × + × - × + × - = 3 3 2 3 3
To calculate the phosphoric acid mole fraction we assumed that the concentration of phosphoric acid remained constant and the local acid partial pressure was in equilibrium with the acid inside the agglomerate. Therefore, from the experimental data provided by [36] it was possible to generate a set of equations which coupled the concentration of phosphoric acid with the water partial pressure. Therefore, for concentrations higher than 96% the concentration term was kept constant and equal to 0.96.

As result of insufficient studies on hydrogen solubility and diffusion in hot phosphoric acid, the diffusion coefficient and Henry's constant cannot be correlated to experimental data. Thus, it was assumed that these parameters behave in the same way as in water systems [26]:

2 2 4 O H H H = (40) 2 2 2 O liq H liq D D = (41)
Where, is the diffusion coefficient of hydrogen in phosphoric acid.

The proton conductivity has a big impact on fuel cell performance, particularly in intermediate temperatures PEM fuel cells, due to the relative low membrane conductivity.

Bouchet and Siebert [37] and Ma et al. [38] proposed that the conductivity of PBI membranes doped with phosphoric acid obeyed an Arrhenius law:

      - = RT E T B M m exp κ (42)
Where, B is the pre-exponential coefficient and M E is the activation energy.

From the experimental data provided by Scott et al. [27] for a doping level of 5.6 mol of H 3 PO 4 per PBI repeat unit, it was possible to generate two functions that correlated these two parameters with the relative humidity at a constant doping level.

( )R RH E M 47 . 698 86 . 174 + = (43)
( )

6535 . 8 5677 . 0 exp + = RH B ( 44
)
In the catalyst layers the phosphoric acid that comes from the membrane was responsible for its proton conductivity. To calculate the conductivity of phosphoric acid the experimental data was fitted to a polynomial function [36]. The effective conductivity was obtained by the Bruggman model: Isotropic conductivity for all regions where the transport of electrons occurs was assumed.

This means that the model does not distinguish between in-and through plane conductivities.

For the GDL the value reported by the supplier for the through plane conductivity, 1250 S/m was used [START_REF] Brandrup | Polymer handbook[END_REF]. At the MPL and catalyst layer the electron conductivity was calculated by taking into account the volume fraction of the materials responsible for the conductivity. The effective electron conductivity was obtained by the following equation [40]:

Where, s k is the electrical conductivity of the solid phase (carbon black in the MPL and both platinum and carbon black in the catalyst layer), e ε is the volume fraction of the non- conductor materials.

The electrical conductivity of Vulcan XC-72 carbon black was reported by Pantea et al. [41] to be equal to 450 S m -1 . The electrical conductivity of platinum, Pt k , was obtained from the data provided by Powel and Tye [42]. From these data a linear function was generated to correlate the electrical conductivity of platinum ( ) s k ln with temperature ( )

T ln : ( ) ( ) 3 . 21 ln 93 . 0 ln + - = T k Pt ( 47
)
The electric conductivity of the titanium flow fields and carbon paper was assumed constant with temperature and equal to 1.92 x 10 6 S m -1 [43] and 1250 S m -1 [44],

respectively.

The thermal conductivity of the gas mixture fluid T k , was obtained from a semiempiriacal equation developed by Wilke in 1950 [45]:

∑ ∑ = = Φ = n i n j ij j i T i fluid T x k x k 1 1 , , (48)
In which:

1 1 8 1                         +         + = Φ - i j j T i T j i ij M M k k M M (49) ij j i i T j T ji M M k k Φ = Φ , , (50)
i T k , is the thermal conductivity of the species and it was obtained for each species through a linear regression of experimental data [46]

= × + × = - - R T k H T (52) 998 . 0 10 088 . 1 10 453 . 5 2 2 5 , 2 = × + × = - - R T k N T (53) 996 . 0 10 32 . 2 10 189 . 1 2 2 4 , 2 = × + × = - - R T k O H T (54)
For platinum, carbon back, and PTFE the thermal conductivity was also obtained from a polynomial function. In the case of PTFE a linear equation proposed in reference [47]

= + × - × + × - = - - - R T T T k Pt T (56) 998 . 0 979 . 2 10 869 . 2 10 048 . 1 2 3 2 6 , = + × - × = - - R T T k CB T (57)
In the same way as electric conductivity the thermal conductivity of phosphoric acid also depends on temperature and concentration. From the data provided by Turnbull [48] it was possible to generate an equation which couples the concentration of phosphoric acid with the thermal conductivity. At any given concentration a linear expression was obtained between 4 3 , PO H T k and T :

( ) (
)

1
4 3 4 3 - - - - - × × × + × + × - = PO H PO H PO H PO H w w T w λ (58)
The thermal conductivities of the cell body (titanium grade 2), membrane and carbon paper (Toray carbon paper present in the GDL) were assumed to be temperature independent as a result of insufficient data:

Formatted: Line spacing: 1.5 lines

Formatted Table

Formatted: Indent: First line: 0.75 pt, Line spacing: 1.5 lines Formatted: Line spacing: 1.5 lines Formatted: Indent: First line: 0.75 pt Formatted Table Formatted [49,50] -Carbon paper: 1.7 W K -1 m -1 [44] -Titanium grade 2: 16.4 W K -1 m -1 [43] To calculate the effective thermal conductivity of the porous media T k (GDLs, MPLs, and catalyst layers) an equation proposed by Hashin and Shtrikman [51] in 1962 for isotropic medium was used:

( ) ( ) solid T fluid T solid T solid T T k k k k k , , , , 3 1 2 1 2 ε ε - + + + - = (59)
Where, solid T k , is the effective thermal conductivity of the solid phase.

An average of the solid species thermal conductivity was obtained:

MPLs the at , , , CB T CB PTFE T PTFE solid T k k k ε ε + = (60)
layers catalyst the at

4 3 4 3 , , , , Pt Pt PO H T PO H CB T CB PTFE T PTFE solid T k k k k λ ε ε ε ε + + + = (61)
Where, PTFE ε

, CB ε , 4 3 PO H ε
,and Pt ε are the volume fractions of PTFE, carbon black, phosphoric acid and platinum, respectively. In the GDL solid T k , was equal to the conductivity of the carbon paper.

Solution Technique

To solve the governing equations COMSOL Multiphysics version 3.5 was used which is a finite element analyses and solver software. The computations were performed on a 64 bit

Windows platform with 8 GB of RAM, and Intel Core 2 Quad 2.5 GHz processor. Geometry A was divided in 2240 mesh elements, 1536 were triangular and 1004 were quadrilateral.

Geometry B was divided into 4913 mesh elements, 2073 were triangular and 2840 were quadrilateral. In the MEA and gas channels we used a quadrilateral refine mesh, in the other areas of the computational domain a triangular mesh with a growing rate of 1.8. The critical regions of the computational domain were the MEA and gas channels; to these regions was done refinement studies in previous works. With those studies we didn't find any significant numerical error. For example, when was used a much more refined mesh (with more 51% of degrees of freedom) the predicted current density only differed 0.2 mA cm -2 from the base case used in this work.

A direct parametric non-liner solver was applied, and the cell potential, Cell V , was used as the parametric variable. With Geometry B the polarization curve exhibited three different regions, identified in Figure 5. From open circuit to 0.7 V the potential rapidly decreased due to activation losses.

Between 0.7 and 0.2 V the potential fell almost linearly with current density, as a result of Ohmic losses, caused by ionic and electronic resistance. In this particular system the electronic resistance was insignificant when compared with the ionic resistance in the catalyst layers and membrane. Typical values of ionic resistance for PBI doped with H 3 PO 4 can be found in the work developed by Ma et al. [38]. They observed a value of 37.57 Ω cm for a doping level of 630%, at 150 °C, and 5% of relative humidity. Pure H 3 PO 4 (species responsible for ionic conductivity in the catalyst layers) has typically an ionic resistance lower than that of PBI doped with H 3 PO 4 . At 150 °C and 99.72 %wt, MacDonald and Boyack [36] measured a value of 2.053 Ω cm . Both these values are higher than the electric resistance of carbon black, 0.222 Ω cm [41].

At potentials below 0.2 V there was a rapid fall in voltage. This was caused by mass transport limitations, in particular the oxygen diffusion and solubility in hot, concentrated phosphoric acid. The value of oxygen diffusion coefficient, given by Eq. ( 37), was very low

for the entire range of potential, 1.08 x 10 -9 m 2 s -1 . Similarly the solubility (Eq. ( 36)), was The performance predicted with Geometry A dropped shapely and linearly between potentials of 0.7 and 0 V, reaching a maximum of 0.46 A cm -2 . This unrealistic behaviour was caused by an inefficient utilization of the catalyst layers, due to the position of the boundary conditions (Eqs. ( 31)-( 32)) to solve the conservation of charge for the solid phase (Eq. ( 9)). Figure 6 shows the current density profile at the cathode catalyst layer at 0.4 V for Geometry A. As shown the reaction did not occur uniformly through the catalyst layer. The maximum current density was obtained in the adjacent areas of the ribs, 0.725 A cm -2 , whilst, mid way between the ribs the current density was very low, 0.050 A cm -2 . The linear trend of the polarization curve below 0.7 V demonstrates that the electric and protonic resistance had a predominate effect. The electrons could only flow to the external circuit through small regions at the edges of the MEA, therefore, they experienced high electric resistance to flow from the reaction sites to the land thus the rate of the reactions decreased, mainly in the middle of the catalyst layers (farther region from the land).

Figure 7 shows the current density profile at 0.4 V at the cathode catalyst layer for Geometry B. In contrast with Geometry A the reaction rate was more uniform through the catalyst layer, due to the presence of the ribs between channels. The minimum current density was obtained at the edges of the catalyst layer, 0.085 A cm -2 . This low value was not caused by the electric resistance but by mass transport limitations, i.e. oxygen molecules could not reach these areas easily. To analyse to what extent the non-isothermal conditions affected performance, the results predicted by the model were compared with results obtained at isothermal conditions. For both cases the simulations were carried out for a PEM fuel cell utilising 50% Pt/C at the cathode and 20% Pt/C at the anode at 150 °C and at atmospheric pressure. Air was used at the cathode, hydrogen at the anode and Geometry B was used to perform the simulations.

Figure 9 shows the obtained cell voltage vs. current density behaviour. From open circuit potential to 0.4 V the performance was very similar for both isothermal and non-isothermal conditions. However below this point the non-isothermal model predicted a much better performance. The limiting current for the non-isothermal case was 15% higher than the isothermal case. This higher performance, at high currents indicates an improvement in species transport, as can be seen in Table 3. The improvement in species transport can only be explained by a significant temperature increase inside the fuel cell. Figure 10 shows the average temperature of the MEA at different cell potentials. The temperature was practically constant, at 150 0 C until potentials of 0.6 V, however, below this potential the temperature increased until a maximum of 181 °C at 0 V.

Local overpotential distribution

Figure 11 shows the local overpotential distribution predicted for a section of the cathode catalyst layer at three different cell potentials (0.6 V, 0.4 V and 0.2 V). The simulations were carried out for a PEM fuel cell utilising 50% Pt/C at the cathode and 20% Pt/C at the anode at 150 °C and at atmospheric pressure. Air was used at the cathode, hydrogen at the anode and Geometry B was used to perform the simulations. From Figure 11 two general results can be seen. The local overpotential was not uniform within the catalyst layer and the maximum was always reached under the land. At 0.6 V the overpotential variation was nearly constant, the difference between the maximum and minimum value was 2 mV. As the cell potential decreased the overpotential variation increased, at 0.4 V the difference between the maximum and minimum value was 14 mV, and at 0.2 V it reached 47 mV. It can be concluded that the overpotential variation was small for the entire range of cell voltage; however, even small differences are important on performance prediction because it is related exponentially to the current density. This observation is important because it points out one big limitation of models that assume an uniform overpotential distribution. Examples of these models are those that assume the catalyst layer as a reactive boundary. This assumption can lead to an over predicted performance, especially at high current densities where the overpotential variation is higher.

Dead and active reaction zones

The distribution of the reactants reaction rate in the catalyst layers can be an efficient way to measure the effective active area of the catalyst layers. The reaction rate was calculated from the volumetric current density, as is demonstrated by Eq. ( 63).

nF j RR = (62)
Where, RR is the reaction rate.

Figure 12 shows the oxygen reaction rate (ORR) distribution in the catalyst layer for Geometry B at two different positions (catalyst layer/MPL interface and catalyst layer/membrane interface) and at three different cell potentials (0.6 V, 0.4 V and 0.2 V). The simulations were carried out for a PEM fuel cell utilising 50% Pt/C at the cathode and 20%

Pt/C at the anode at 150 °C and at atmospheric pressure. Air was used at the cathode and hydrogen at the anode. From Figure 12 the main feature can be seen was that ORR was nonuniform in both directions, it varied with thickness and length of the catalyst layer.

It is evident that the variation in ORR distribution with the thickness of the catalyst layer, was more significant at high current densities. At 0.6 V the ORRs predicted by the model at the catalyst layer/MPL and catalyst layer/membrane interfaces were very similar, the mean values of both predictions only differed by 5.70 mol m -3 s -1 (2.75 mA cm -2 in terms of current density). On the other hand, at 0.2 V the predictions for these two interfaces were significantly different. The mean value predicted for the catalyst layer/MPL interface (3.49 x 10 3 mol m -3 s -1 ) was 42% higher than the predicted value at the catalyst layer/membrane interface (2.03 x 10 3 mol m -3 s -1 ). These results suggest that the reaction at low potentials occurred mainly in the adjacent region to the MPL; near the membrane the catalyst layer can be considerate a "dead zone". As mentioned before the ORR also varies with the length of the catalyst layer. Two distinct regions can be observed in Figure 12. The first one corresponded to the extremities of catalyst layer. This region was characterized by very low values of ORR rate for all potential ranges. This means that oxygen molecules cannot diffuse easily to these regions, and as consequence the catalyst layer was not efficiently used. The second region is between the extremities of the catalyst layer. In this region the ORR distribution followed two distinct trends. At high potentials the ORR was nearly constant, however at low potentials the trend was a sinusoidal function with a maximum under the gas channel regions and a minimum under the land regions. Sun et al. [52] made a similar study for a Nafion system. At low potentials they found the same sinusoidal trend, with maximum ORR located under the flow channels. At these conditions the electrochemical reaction is very rapid; therefore, oxygen diffusion becomes the limiting factor. At high potentials Sun et el. observed that the maximum ORR moved from regions under the flow channel to regions under the land. This behaviour was caused by the low ORR; the oxygen molecules diffused fast enough to the reaction sites, therefore, the limiting factor was ohmic losses between the regions under the channel and under the land. This phenomenon was not observed in the present study; at high potentials the ORR was uniform. This fact can be explained by the higher solid phase conductivity used in this paper in comparison with the one used by Sun et al. Kulikovsky et al. [53] showed that the most active catalyst area was dependent on the choice of parameters, in particular the solid phase conductivity, s k . They found that an increase of s k leads to a more uniform distribution of the reaction rate. In this paper the average solid phase conductivity in the catalyst layer was 1000 times higher than the value used by Sun et al. due to the introduction of the platinum conductivity in the calculation of s k .

Temperature distribution

The temperature distribution through the MEA is a very important consideration as it can

give important information about the real temperature in the centre of the cell and it can be useful to predict possible problems concerning thermal stresses and hot spots.

In the apparatus used in our laboratory the temperature in the fuel cell was measured by a thermocouple embedded in the collector plates. The data collected by the thermocouple was used by an on/off control system to control the heaters. Figure 13 shows the temperature profile and the heat flux direction for Geometry B at 0.9 and 0 V. The simulations were

Field Code Changed

Field Code Changed Deleted: 10

Deleted: [47] Deleted: [48] Deleted: 11 One main observation can be seen from Figure 13, is that the heat source shifted with potential. At 0.9 V the heat source was located, as expected, at the heaters, however at 0 V the heat source was located at the centre of the domain (MEA region). In addition, the temperature distribution at 0.9 V was more uniform than at 0 V. The difference between the two temperature limits at 0.9 V was 76 K, in contrast, at 0 V it was 116 K. The average temperature of the MEA for the all range of potentials was shown in Figure 10. As can be seen from Figure 8 the average temperature of the MEA at low potentials was much higher than 150 ºC. This meant that the temperature in the middle of the cells was not efficiently controlled with the current location of the thermocouple. This feature can lead to phosphoric acid dehydration and consequently a decrease in proton conductivity.

It is also important to analyse the temperature profile through the membrane in order to detect possible hot spots or abnormal temperature differences. Figure 14 shows the temperature profile predicted as a function of the membrane length at three different cell potentials (0.6 V, 0.4 V and 0.2 V) for Geometry B. The trend in the temperature profile was very similar to the ORR profile shown in Figure 12, which indicates the importance of the reactive heat for the overall heat generation. Once more two distinct regions can be identified. The first one was at the edges of the membrane. This region was characterised by a sharp fall in temperature for the entire potential range. At 0.2 V the temperature fell from 449 K to 432 K, a large difference for such a small distance (8x10 -3 m). This temperature difference could lead to degradation problems due to thermal stress. Between the edges of the membrane, the temperature was constant at high potentials but at low potentials it was nonuniform and followed a sinusoidal trend with hot spots under the channel regions. The temperature difference between the regions under the flow channel and the land was not so marked at that at the membrane edges (only 3 K at 0.2 V), therefore thermal degradation problems in middle regions of the membrane are less likely to occur.

To further illustrate the impact of the reactive heat on the overall heat generation and consequently the temperature profile, the total reactive and ohmic heat were calculated and represented as a function of the cell potential in Figure 15. Both reactive and ohmic heat increased as the potential decreased as a result of a greater current density. However, the Ohmic heat at low potentials was low when compared with reactive heat; for example at 0.2

V it only corresponded to 12.2% of the total heat. This explains why the temperature profile

Conclusions

A two dimensional agglomerate model was developed to simulate the performance of an intermediate temperature PEMFC equipped with a phosphoric acid doped PBI membrane.

The model was capable of predicting the polarization performance, temperature distribution and all transport characteristics.

Two different geometries were analysed and the results compared. It was found that the model predictions are highly dependent on the adopted geometry. It became clear that Geometry A (along the channel direction) was only suitable for fuel cell modelling with the introduction of artificial boundary conditions at the GDL/gas channel interface that simulates the land effect. On the other hand, Geometry B (cross the channel direction) was able to simulate the general trend in fuel cell performance and, in addition, the simulation predictions matched the experimental results. The effect of the non-isothermal conditions on performance was investigated and compared with the results obtained at isothermal conditions. It was concluded that the performance increased when non-isothermal conditions were used in the simulations. This fact was caused by a substantial temperature increase of the MEA at high current densities, leading to an improvement of transport properties.

It was also found that the local overpotential distribution and the oxygen reaction rate were non-uniform through the catalyst layer. The maximum reaction rate was obtained under the channel regions and the lower rate was obtained at the edges of the MEA. This feature, leads to large temperature difference between the central regions and the edges of the MEA.

The temperature difference becomes potentially very problematic at the membrane because it could lead to thermal degradation due to thermal stresses. From the analyses of the oxygen reaction rate it was shown that the catalyst layer was not efficiently used. In particular, the extreme regions of the catalyst layer had a very low reaction rate for the all ranges of potential, this region represents 15% of the total catalyst layer volume. In the electrochemical reactions there are three different participants: electrons, protons and gas reactant/products. The reaction only takes place on the catalyst surface where all the three species have access. The electrons travel through the electric conductive solids, the protons through the electrolyte (phosphoric acid in this particular system) and the gas species should have an easy access to the reaction sites. In order to enhance performance and to reduce platinum loading the catalyst structure should maximise the three phase boundary formed by the electrolyte, solid and void phases. To simulate all the transport and electrochemical processes within the catalyst layer the model treats it as spherical agglomerates of catalyst particles with porous inter-agglomerate space, where each agglomerate is covered by a thin film of electrolyte.

The geometric parameters of the catalyst layer and agglomerate are strongly dependent on the fabrication process. Therefore is not easy to find consensual values in the literature. The thickness of the catalyst layers was measured from the analyses of scanning electron micrograph (SEM) images. For cathodes with 0.4 mg Pt cm -2 , of 50% Pt/C, and 40% wt PTFE the average thickness was 12.5 µm and for anodes with 0.2 mg Pt cm -2 , 20% Pt/C, and 40% wt PTFE the average thickness was 13.8 µm. In this paper an average value of 0.5 µm was adopted for the agglomerate diameter. A wide range of values have been reported in literature, Broka and Ekdunge [28] estimated the agglomerate diameter to be between 1 and 5

µm through the analyses of SEM images. Using the same technique Siegel et al. [54] reported an average value of 6.1 µm. Deleted: [50] Deleted: [47] Deleted: [51] Deleted: 63 is the mass of phosphoric acid per area, A is the catalyst layer area (9 cm2), solids V is the total volume of solids (platinum, carbon black, and PTFE) in the catalyst layer, and agg V is the agglomerate volume. The effective porosity of the catalyst layers will depend on the amount of phosphoric acid inside the catalyst layers.

H r V V V V A L -                               + -         - - = 3 1
Knowing the geometric properties of the catalyst layer and agglomerate the following procedure was adopted to describe the transport and electrochemical reaction processes. [ref, ref] In the catalyst layers the reactants are transported to the catalyst surface as molecular species dissolved in phosphoric acid. The concentration at the gas-liquid interface was obtained by Henry's Law.

i i l g i p H C = / ( 64
)
Where, i H is the Henry constant, and i p is the partial pressure of the reactants.

The reactant mass balance in the catalyst layer was described by the following equivalent equations:

0 ' = + ∇ i agg i N a N (65) 0 = + ∇ i i R N (66)
Where, i N is the overall reactants mole flux, ' i N is the reactants mole flux in the phosphoric acid film, and i R is the reaction rate.

The transport in the liquid phase could be described by Fick's Law:

r C D N i i liq i ∂ ∂ = , ' (67)
Where, i C is the reactants concentration, and r is the radius.

It was assumed that the phosphoric acid film was present in the agglomerate evenly, and the thickness of the film was much smaller than the size of the agglomerate. Therefore the reactant flux can be described by the following equation:

( ) δ s l i l g i i liq i C C D N / / , ' - = (68)
Where, s l i C / is the reactants concentration at the surface of the agglomerate.

The overall oxygen reduction reaction and hydrogen oxidation reaction was assumed to follow first-order kinetics with respect to the reactant concentration [52,57]:

i c i C k R = (

( )

CAT s l i c r i agg C k E R ε - = 1 / , (70)
Where, i agg R , is the reaction rate based on the agglomerate volume, and r E is the effectiveness factor. However to apply this approximation it is necessary to introduce a new factor which takes into account the geometry of the agglomerate and the diffusion resistance of reactant molecules inside the agglomerate. This factor is the effectiveness, which can be described, for spherical agglomerate, by the following equation [START_REF] Thomas | Introduction to the principles of heterogeneous catalysis[END_REF]:

( )         Φ - Φ Φ = 1 tanh 1 3 r E (71)
Where, Φ is an adimensional number called the Thiele's modulus:

eff i liq c agg D k r , = Φ (72)
Where, eff i liq D , is the effective reactants diffusion inside the agglomerate.

The oxygen consumption rate was related to the current density according to:

( ) ( ) s l i CAT c r CAT i agg i C k nFE nFR nFR j / , 1 1 ε ε - = - = = (73)
From Eq. ( 68

( )( ) 1 , / / 1 1 -         + - + = i liq agg agg agg CAT c r l g i s l i D a r r k E C C δ δ ε (76)
The volumetric current density, j , described by Eq. ( 12) was obtained through the combination of Eqs. ( 73)-( 76). The reaction rate (Eq. ( 13)) was obtained from Eqs. ( 11)-(73). Convective heat flux (W m -2 ) q Joulde Joule heating (W m -3 ) Pr Prandtl number q reaction Heat of reaction (W m -3 ) p sat Saturation pressure of water (Pa) R Universal gas constant (J mol ... [2] ... [4] ... [3] ... [12] ... [9] ... [6] ... [10] ... [7] ... [11] ... [8] ... [19] ... [16] ... [13] ... [17] ... [14] ... [18] ... [15] ... [26] ... [23] ... [20] ... [24] ... [21] ... [25] ... [22] ... [33] ... [30] ... [27] ... [31] ... [28] ... [32] ... [29] ... [1] ... [37] ... [34] ... [38] ... [35] ... [START_REF] Brandrup | Polymer handbook[END_REF] ... [36] ...

)+(1-ε)/3k T,solid ) -1 b) [40] W m -1 K -1 a
Calculated from experimental data b ks and k T, solid were treated as dependent variables. Deleted: [38] Deleted: [38] Deleted: [START_REF] Brandrup | Polymer handbook[END_REF] Deleted: (42) Deleted: (42) Deleted: [29] Deleted: [29] Deleted: [40] Deleted: [40] Deleted: [40] Deleted: [41] Deleted: [42] Deleted: [41] Deleted: [43,44] Deleted: [41] Deleted: [41] Deleted: [45] Deleted: [40] Deleted: (45) Deleted: (45) Deleted: [40] 4 Reactants and products mass fraction at gas channel/collector plate interface predicted ad isothermal conditions for a PEMFC utilising 50% Pt/C at the cathode and 20% Pt/C at the anode at 150 °C and at atmospheric pressure. Air was used at the cathode and hydrogen at the anode. (A) ▬▬ H 2 mass fraction, (B) ▬▬ O2 mass fraction, (C) ▬ ▬ H 2 O mass fraction. Fig. 5 Comparison between modelled and experimental data for a PEMFC utilising 50% Pt/C at the cathode and 20% Pt/C at the anode. Air was used at the cathode and hydrogen at the anode at atmospheric pressure. The controller set point was 150 °C. (A) ▬▬ experimental data, (B) ♦ model results using Geometry B, (C) ▲ model results using Geometry A. Fig. 6 Current density profile at the cathode catalyst layer for Geometry A at 0.4 V. The simulation was carried out for PEMFC utilising 50% Pt/C at the cathode and 20% Pt/C at the anode. Air was used at the cathode and hydrogen at the anode at atmospheric pressure. The controller set point was 150 °C. Fig. 7 Current density profile at the cathode catalyst layer for Geometry B at 0.4 V. The simulation was carried out for PEMFC utilising 50% Pt/C at the cathode and 20% Pt/C at the anode. Air was used at the cathode and hydrogen at the anode at atmospheric pressure. The controller set point was 150 °C. Fig. 8 Comparison between modelled and experimental data for a PEMFC utilising 50% Pt/C at the cathode and 20% Pt/C at the anode. Air was used at the cathode and hydrogen at the anode at atmospheric pressure. The controller set point was 150 °C. (A) ▬▬ experimental data, (B) ■ model results using Geometry A with the ribs effect. Fig. 9 Comparison between modelled and experimental data for a PEMFC utilising 50% Pt/C at the cathode and 20% Pt/C at the anode. Air was used at the cathode and hydrogen at the anode at atmospheric pressure. The controller set point was 150 °C. (A) ▬▬ experimental data, (B) ♦ model results using Geometry B at non-isothermal conditions, (C) ▲ model results using Geometry B at isothermal conditions. Fig. 10 Predicted average temperature of the MEA for a PEMFC utilising 50% Pt/C at the cathode and 20% Pt/C at the anode. Air was used at the cathode and hydrogen at the anode at atmospheric pressure. The controller set point was 150 °C. Fig. 11 Local overpotential distribution in a section of the catalyst layer under three different cell potential for a PEMFC utilising 50% Pt/C at the cathode and 20% Pt/C at the anode. Air was used at the cathode and hydrogen at the anode at atmospheric pressure. The controller set point was 150 °C. A) 0.6 V, B) 0.4 V, C) 0.2 V. Fig. 12 O 2 reaction rate distribution of the catalyst layer under three different cell potential and at two different positions of the catalyst layer for a PEMFC utilising 50% Pt/C at the cathode and 20% Pt/C at the anode. Air was used at the cathode and hydrogen at the anode at atmospheric pressure. The controller set point was 150 °C. A) 0.6 V, B) 0.4 V, C) 0.2 V. (A) ▬▬ catalyst layer/MPL interface, ▬▬ catalyst layer/Membrane interface. Fig. 13 Temperature distribution and heat flux direction under two different cell potential for a PEMFC utilising 50% Pt/C at the cathode and 20% Pt/C at the anode. Air was used at the cathode and hydrogen at the anode at atmospheric pressure. The controller set point was 150 °C. A) 0.9 V, B) 0 V. 14 Temperature distribution of the membrane under three different cell potential for a PEMFC utilising 50% Pt/C at the cathode and 20% Pt/C at the anode. Air was used at the cathode and hydrogen at the anode at atmospheric pressure. The controller set point was 150 °C. (A) ▬▬ 0.2 V, (B) ▬▬ 0.4 V, (C) ▬ ▬ 0.6 V. Fig. 15 Total heat released during the operation of a PEMFC utilising 50% Pt/C at the cathode and 20% Pt/C at the anode. Air was used at the cathode and hydrogen at the anode at atmospheric pressure. The controller set point was 150 °C. (A) ▬▬ total heat, (B) ▬▬ reactive heat, (C) ▬ ▬ ohimc heat.

Table 3 Transport properties at isothermal and non-isothermal conditions

  current effects (the membrane used was ~40µm thick). A small current ~10µA

  electrical conductivity in the porous media, eff elect k is the effective electrolyte conductivity in the porous media, s φ is the solid phase potential, elect φ is the electrolyte phase potential, s St φ and elect St φ are the source terms representing the volumetric transfer current for the solid phase, and the electrolyte phases respectively.

2 HC 2 OC

 22 is the hydrogen concentration (solubility) obtained using Henry's law, is the oxygen concentration (solubility) obtained using Henry's law, and,

k

  is the thermal conductivity of the external air, L is the characteristic length.The Nusselt number is a dimensionless parameter which accounts for the ratio of convective and conductive heat transfer through the boundary.

)

2 HH

 2 is the Henry's constant of hydrogen in phosphoric acid and 2 H liq D

  of phosphoric acid in the catalyst layer.

  75 x 10 -6 mol Pa -1 m -3, which contributed to the low limiting current density of 1.43 A cm -2 .

Figures 5 and 6 5 Deleted: 6 Deleted: 4 Deleted: 5 Deleted

 65645 Figures 5 and 6 demonstrate that Geometry A is not suitable for fuel cell simulations. A realistic performance could only be achieved with the introduction of artificial boundary conditions at the interface GDL/Gas channel, which simulate the ribs effect.Such a simulation was done and the results are shown in Figure8. With these new boundary conditions the area responsible for the electron diffusion from the MEA to the bipolar plates increases. For that reason, the current density was distributed more uniformly in the catalyst layer, leading to more realistic results.

  a PEM fuel cell utilising 50% Pt/C at the cathode and 20% Pt/C at the anode at 150 °C and at atmospheric pressure. Air was used at the cathode and hydrogen at the anode.

  length was similar to the ORR profile through the catalyst layer length.

  can be simplified if the overall reaction rate only depends on the reactant concentration at the outer surface of agglomerate.

Fig. 1

 1 Figures Caption

Formatted: 1 Deleted: 1 Deleted

 11 Fig.14 Temperature distribution of the membrane under three different cell potential for a PEMFC utilising 50% Pt/C at the cathode and 20% Pt/C at the anode. Air was used at the cathode and hydrogen at the anode at atmospheric pressure. The controller set point was 150 °C. (A) ▬▬ 0.2 V, (B) ▬▬ 0.4 V, (C) ▬ ▬ 0.6 V. Fig.15Total heat released during the operation of a PEMFC utilising 50% Pt/C at the cathode and 20% Pt/C at the anode. Air was used at the cathode and hydrogen at the anode at atmospheric pressure. The controller set point was 150 °C. (A) ▬▬ total heat, (B) ▬▬ reactive heat, (C) ▬ ▬ ohimc heat.

Table Formatted :

 Formatted Indent: Before: 1.7 pt, First line: 0 pt, After: -1.95 pt

	Formatted Formatted

Table Formatted :

 Formatted Indent: First line: 0 pt, After: 1.7 pt, Tabs: 50.85 pt, Left

	1	2	1	2	1	4	2
		,
		,
	F o r
	P
		e
		e r
					R
						e
							v i e
							w	Formatted Table
							16

  Experimental data for platinum and carbon black were used to generate the polynomial equations computed by the model to calculate the thermal conductivity:

																was
	used:
						, k PTFE T	=	86 . 4	×	10	4 -T	+	253 . 0	(55)
	,	275 . 5	10	4	3	543 . 2	10	5	2	321 . 2	10	2	85 . 77	2	999 . 0

:

  After: 1.7 pt

	Formatted

Table Formatted Table Formatted :

 FormattedFormatted Left

	-Membrane: 0.41 W K -1 m -1
	F o r
	P
	e
	e r
	R
	e
	v i e	Formatted Table
	w
	17

  The membrane PA doping level was 630 %. Air was used at the cathode and hydrogen at the anode. The open circuit potential (OCP) was obtained from Eq.(15), and was the same for both geometries, 1.16 V. Until 0.7 V a good agreement between experimental data and model results was obtained for both geometries. This region is characterized by low current densities and a sharp decrease in voltage due to activation losses, mainly in the cathode side. Below this potential the performance predicted for each geometry was different. The results predicted by the model agreed with the experimental data when geometry B was used. On the other hand, with geometry A the predicted performance was much lower than the experimental performance.

	3 Results and Discussion	Formatted: Portuguese Portugal
	3.1 Polarisation curve
	Figure 5 compares model predictions with the corresponding experimental data for a PEM	Formatted Table
	fuel cell utilising 50% Pt/C at the cathode and 20% Pt/C at the anode at 150 °C and at
	atmospheric pressure.
		Formatted: Left
		Formatted Table
	F o r
	P
	e
	e r
	R
	e
	v i e
	w
	18

Table 2 Physical and electrochemical parameters for base case

 2

	Titanium
	Electric conductivity	7.407x 10 5 [59]	S m -1
	Thermal conductivity	7.5 [59]		W m -1 K -1
	PTFE
	Thermal conductivity	2.5 x10 -1 [47]		W m -1 K -1
	Membrane
	Protonic conductivity	Eq. (42)		S m -1
	Thermal conductivity	40 [26]		W m -1 K -1
	GDL
	Porosity	0.78 [39]
	Permeability	3.01 x 10 -11a)		m 2
	Electric conductivity	1250 [39]		S m -1
	Thermal conductivity	-2k T,solid +(ε/(2k T,solid +k T,fluid )+(1-ε)/3k T,solid ) -1 b) [40]	W m -1 K -1
	MPL
	Porosity	0.6
	Permeability	1.97 x 10 -12a)		m 2
	Electric conductivity	k S (2-2ε)/(ε+2)	[40]	S m -1
	Thermal conductivity	-2k T,solid +(ε/(2k T,solid +k T,fluid )+(1-ε)/3k T,solid ) -1 b) [40]	W m -1 K -1
	Catalyst layers
		Anode	Cathode
	Pt loading	2 x 10 -3	4 x 10 -3	Kg m -2
	Exchange current density	1.44 x10 3 [60]	2.63 x 10 -4 [61]	A m Pt	-2
	Gama	1 [60]	1 [62, 63]
	Transfer coefficient	0.5 [60]	0.74
	Activation energy	16.9 x 10 3 [60]	72.4 x 10 3[64]	J mol -1
	Permeability	2.96 x10 -13a)		m 2
	Electric conductivity	ks(2-ε)/(ε+2) b) [40]	S m -1
	Protonic conductivity	Eq. (45)		S m -1
	Thermal conductivity	-2k

T,solid +(ε/(2k T,solid +k T,fluid

Wiley-VCH Fuel Cells

Acknowledgements

The authors acknowledge the support of EPSRC through the Supergen; fuel cell consortium award.

Font: 12 pt Page 30: [62] Formatted 4/7/2010 11:07:00 AM Font: 12 pt, Bold Page 30: [62] Formatted 4/7/2010 11:07:00 AM Font: 12 pt Page 30: [62] Formatted 4/7/2010 11:07:00 AM Font: 12 pt, Italic Page 30: [62] Formatted 4/7/2010 11:07:00 AM Font: 12 pt Page 30: [63] Formatted 4/7/2010 11:07:00 AM Font: 12 pt, Italic Page 30: [63] Formatted 4/7/2010 11:07:00 AM Font: 12 pt Page 30: [63] Formatted 4/7/2010 11:07:00 AM Font: 12 pt, Bold Page 30: [63] Formatted 4/7/2010 11:07:00 AM Font: 12 pt Page 30: [63] Formatted 4/7/2010 11:07:00 AM Font: 12 pt, Italic Page 30: [63] Formatted 4/7/2010 11:07:00 AM Font: 12 pt Page 30: [64] Formatted 4/7/2010 11:07:00 AM Font: 12 pt, Italic Page 30: [64] Formatted 4/7/2010 11:07:00 AM Font: 12 pt Page 30: [64] Formatted 4/7/2010 11:07:00 AM Font: 12 pt, Bold Page 30: [64] Formatted 4/7/2010 11:07:00 AM Font: 12 pt Page 30: [64] Formatted 4/7/2010 11:07:00 AM Font: 12 pt, Italic Page 30: [64] Formatted 4/7/2010 11:07:00 AM Font: 12 pt Page 30: [65] Deleted 3/26/2010 2:50:00 PM