Mathematical modeling plays an important role in the design of proton exchange membrane fuel cells. In particular, there is a wide variety of continuum macroscopic models that are considered successful as they agree well with the
Constants used in the model are declared at the top.
|
A |
В |
C |
E |
F |
G |
H |
||
|
1 |
Electroclierni |
мі Model of |
=EM Fuel Cel |
3 7 |
||||
|
2 |
||||||||
|
3 |
Eo |
1.229 |
V |
T |
333.15 |
|||
|
4 |
R |
8.314 |
P H2 |
0 26 |
||||
|
5 |
F |
96487 |
P 02 |
1 |
||||
|
6 |
Rs |
0.0035 |
||||||
|
7 |
xl 1 |
-0.948 |
C 02 |
8.7767E-007 |
||||
|
3 |
xi 2 |
3.10E-003 |
n c |
1 |
||||
|
9 |
xl 3 |
7.60E-005 |
||||||
|
10 |
xl 4 |
-1.9ЭЕ-004 |
||||||
|
11 |
в |
0.2 |
||||||
|
12 |
I L |
42 |
||||||
|
13 |
||||||||
|
14 |
||||||||
|
15 |
E N |
V activ |
V cone |
V ohms |
V cell |
|||
|
16 |
1 |
1.20966506 |
0.268339344 |
0.00481951 |
0.0035 |
0.933006226 |
||
|
17 |
2 |
1.20966508 |
0 312907267 |
0.009758033 |
0.007 |
0.87999976 |
||
|
13 |
0.338977862 |
0.014821594 |
0.0105 |
0.845365624 |
||||
|
19 |
20966508 |
0.357475229 |
0.020016692 |
0.014 |
0.818173159 |
|||
|
20 |
1.20966508 |
0.371822902 |
0.025350341 |
0.0175 |
0.794991836 |
The current is in column B, the open circuit voltage (calculated using the Nernst expression) is in column C, the activation voltage drop is in column D, the effects of reduced concentrations at the electrodes are in column E, and the effects of the membrane and electrode resistance are in column F. Combining all of the these effects results in the cell’s voltage in column G.
FIGURE 7.6
Spreadsheet implementation of a simple electrochemical model of a fuel cell.
 |
|
|
|
experimental data. (Although it should be noted that when a model or simulation has a number of free parameters, then agreement with experimental data is often either fortuitous or obtained through fitting the free parameters.) In general, the model or simulation is considered valid when the experimental and modeled polarization curves are in good agreement. However, it is often found that different models and simulations, that are significantly different from one another, are found to agree with the same experimental data. Increasingly complex models are continually being developed, incorporating full dynamic three-dimensional fuel cell models, multiple phase fluid flows, heat transfer and other processes, that capture the important physics of fuel cells. Rather than explore these increasingly complex models and simulations that are emerging within this field, however, we will look at an earlier model which is more suited to spreadsheet modeling.
The model we will consider was originally developed by Springer et al (J. Electrochemical Soc., 1991). The model is a simple one-dimensional model which is capable of capturing the fuel cell current-voltage performance, and was found to agree reasonably well with experimental data. This model is important, not only because it is the cornerstone from which later models were descended, but also because it emphasized the importance of keeping the membrane well hydrated, and the importance of water management. In particular, the model is steady-state, isothermal and one-dimensional, with the fuel cell consisting of five layers (membrane, diffusion media, and catalyst layers), and importantly this model considered the water content in the membrane (allowing for variable properties in the membrane such as conductivity and the water diffusion coefficient).
In the one-dimensional model all quantities are considered to vary only in the direction perpendicular to the anode and cathode interfaces. Recall that at the anode the hydrogen reacts to yield electrons to the platinum and protons to the membrane material
H2 ^ 2H + + 2e – (7.23)
At the cathode, the oxygen reacts with the protons coming through the membrane to form molecules of water
O – + 2H + ^ H2O (7.24)
The crux of this model, therefore, is to capture the water transport throughout the cell. This includes transport through the porous electrodes, and transport through the membrane electrolyte, to obtain the steady-state water distribution throughout the cell.
A schematic of how the model works is shown in Figure 7.8. The interfaces are labeled 1 through 4, as we go from the anode to the cathode, and the water concentrations are obtained at these outer interfaces first. Then the water content is calculated using differential equations (that we will cover shortly) starting from interface 1 and working our way through to interface 3,
 |
 |
FIGURE 7.8
Schematic of the continuum mathematical model of a fuel cell.
while at the same time the water content is calculated starting from interface 4 and working our way through to interface 3 through the cathode. In this manner we obtain two separate solutions for the water content at interface 3. The ratio of water flux to the molar flux of H2 through the anode is a free parameter, which can be varied until the two solutions agree, giving us a consistent distribution of water content throughout the fuel cell. Note that not only does water enter the fuel cell as humidified fuel and air streams, but it is also generated at the cathode by the oxygen and hydrogen reaction.
The cathode inlet stream is assumed to be saturated with water vapor, and the mole fraction of water is given by
P SAT
xwc _ p (7.25)
where the subscripts W and C represent water and cathode, respectively. PSAT is the saturated pressure at the cathode and Pc is the pressure at the cathode. Essentially this sets up a Dirichlet boundary condition (where we assume that we know the water content at this boundary).
The mole fraction of water and oxygen at interface 4 depend on the mole fraction of water in the cathode inlet stream. In particular, the mole fraction of water at interface 4 is set up as a ratio of the water flow and the total flow.
_ xWcvO + 2(1 + a) (l
vo + (2a + 1) 1 – xWc xon
and, similarly for oxygen
xo4 _ (vo – 1)(1 – xWe) x°N (7.27)
vO + (2a + 1) 1 — xwc xON
where vo is a stochiometric coefficient and the ratio of inlet flux to the total flux across interface 4, and xoN is the inlet dry gas mole fraction of oxygen (and is assumed to be known). The above equations, in a relatively straightforward manner, give us the water and oxygen mole fractions at interface 4.
|
|||
 |
|||
|
|
||
|
|||
|
|||
 |
|||
|
|
||
|
|||
|
|||
|
|
||
 |
|||
 |
|||
 |
|||
|
|||
|
|||
where vH is a stochiometric coefficient and the ratio of inlet hydrogen flux to the total flux across interface 1.
The mole fraction of water at interface 2 can be obtained through integrating the following equation across the anode.
|
|
||
|
|||
where xWa is the water mole fraction in the anode and DWHis the binary diffusion coefficient between water and hydrogen. This differential equation is much easier to solve than Equations 7.28 and 7.29, and a simple analytical equation can be obtained to give us the mole fraction of water at interface 2
|
|||
|
|||
|
|
||
where tA is the anode thickness.
Now we have obtained the water content at interface 2, we must obtain the water content at interface 3, which we will then be able to compare with the solution from earlier which involved the numerical integration of water content from interface 4 to interface 3. When the parameter a is chosen such that these two solutions agree then the simulation is solved.
To get from interface 2 to interface 3 we have to go through the membrane, and for that we must consider A, the membrane water content or the local ratio H2O/SO3 in the membrane. This is calculated directly from the mole fraction of water content. First we calculate the water vapor activity using
xW P
PSAT
where xW is the mole fraction of water, P is pressure and PSaT is the saturation pressure. Next we calculate the water content using a phenomenological equation of the form
( 0.043 + 17.81a – 39.85a2 + 36a3 if a < 1;
14 + 1.4(a – 1) if 1 < a < 3.
It is this water content which we will integrate across the membrane using the differential equation
dA л IM,
dz = [2ndrag – a] PdryD(A) (7.38)
where ndrag is the number of molecules per proton, M is the molecular weight of the membrane, Pdry is the density of dry membrane, and the diffusion coefficient (which depends on water content) is given by
1
(7.39)
From the above system of equations we can obtain the water content at interface 3, A3, both through the integration of water content from the anode and the cathode. We can vary a, the ratio of water flux to the molar flux of H2, until these two solutions converge and we have a complete calculation of the water content throughout the fuel cell. From here we can then calculate the voltage across the fuel cell, and hence obtain the polarization curve.
 |
 |
 |
 |
 |
 |
 |
The conductivity of the membrane as a function of water content, A, has been experimentally determined to be of the form
 |
 |
 |
where T is temperature. From this conductivity, we can calculate the resistance of the membrane using
 |
 |
 |
 |
where the integration is over the entire membrane. To obtain the polarization curve (voltage-current curve) the model accounts for kinetic losses and ohmic losses in the membrane. To calculate kinetic losses, the Tafel equation is used, which relates the current density to the overpotential, q.
where jo is the exchange current density (referenced to pure oxygen at 1 atm at the open cell potential, Voc). Pa is the pressure at the cathode, xliq is the mole fraction of liquid which is generally very small (and often ignored), F is Faraday’s constant, and n is given by
n Voc Vcell JR
where Voc is the open cell voltage, and vcell is the voltage across the cell, which we are trying to obtain in order to generate the polarization curve.
To summarize, the simulation involves calculating the water content at interface 3, A3, from both the anode and cathode sides of the fuel cell such that the two solutions agree. Then we can take the water content throughout the cell and calculate the resistance and voltage across the cell. One of the parameters used throughout these calculations is the current density. Therefore, we can obtain the voltage across the cell as a function of current density, and generate the polarization curve.
Let’s turn our attention to the spreadsheet implementation of the model described above. Figure 7.9 depicts the top of the spreadsheet, where the many constants are contained. Most of these numbers are obtained directly from the literature, but the current density in cell B8 is changed to obtain the polarization (or voltage-current) curve. To the right of Figure 7.9 the
The model uses a large number of constants, most of which are obtained experimentally. The current density in cell В8 is the variable changed when obtaining the voltage-current curves.
from the cathode. When an a is found which results in these two values being close together then the simulation is considered solved.
FIGURE 7.9
Spreadsheet implementation of the continuum mathematical model of a fuel cell.
variable a is determined (a is stored in cell H3). As a function of this variable, the water content at interface 3, A3, is obtained twice (as described in the description of the model above). These two A’s are copied from elsewhere in the spreadsheet model (which we will come back to later) in to the cells H5 and H6. The difference between these two solutions is calculated in H8. The objective, therefore, is to find a value of a which would minimize this difference. While this could be done numerically, using for example the bisection method, in the current implementation this is done by hand (although it doesn’t take too many iterations to converge).
The calculation of the two A3’s, and ultimately the voltage across the cell, is depicted in Figure 7.10. In cells A20 and A23 the mole fractions of water and oxygen, respectively, are calculated at interface 4. In terms of these concentrations, the mole fractions of water and oxygen are calculated by integrating over distance. These integrations (the equations of which are covered above) occur in cells D20:BY20 and D23:BY23, for the mole fraction of water and oxygen, respectively. For example, the code contained in cell E20 is of the form
|
The water and oxygen levels are calculated through the cathode to obtain the membrane water content at the cathode-membrane interface (interface 3)
The water content is calculated through the anode and through the membrane to again obtain the membrane content at the cathode-membrane interface. |
FIGURE 7.10
Spreadsheet implementation of the continuum mathematical model of a fuel cell (continued).
( (1-D20-D23)*(1+$H$3)/$B$12 + (0.5*D20+D23*(1+$H$3))/$B$13)
and corresponds with the calculation of the mole fraction of water in Equation 7.28 above, while for the mole fraction of oxygen we discretize Equation 7.29 and the code in cell E23 is of the form
=D23-$E$15*($B$6*$B$7*$B$10/$B$11) *
((D23*(1-$H$3)+0.5*D20)/$B$14 + (1-D20-D23)/$B$13)
In cell BY20, at the end of the spreadsheet to the right (not shown), the mole fraction of water at interface 3 is obtained. From this it is easy to obtain the water content, and this is then displayed in cell H5 (Figure 7.9).
In cell A28 the mole fraction of water at interface 1 is calculated, and the mole fraction of water at interface 2 is determined in cell C28. The code in cell C28 is
=(A28-(H3)/(1+H3))*EXP(B6*B7*B10*E7/(E5*E6)) + (H3)/(1+H3)
and the integration of mole fraction of water is obtained analytically. To obtain the water content at interface 3, we first calculate the water content at interface 2 and integrate this through the membrane to interface 3. The water vapor activity and water content at interface 2 are calculated in cells E28 and
F28, respectively. The diffusion coefficient of water is calculated from the water content in cells I27:AR27. For example, the diffusion coefficient of water calculated in cell J27 uses the code
(2.563-0.33*J28+0.0264*J28*J28-0.000671*J28*J28*J28)
The water content throughout the membrane is now calculated in cells I28:AR28. For example, the code in cell J28 is
=I28+$E$15*(2*2.5*128/(22) – $H$3)*($B$10*$E$10/($E$11*I27))
which is a simple finite difference approximation of the differential equation for water content (Equation 7.38). This ultimately yields the water content at interface 3, in cell AR28 (not shown) and this is then displayed in cell H6 (Figure 7.9). By comparing these two A’s (in cells H5 and H6) we can find a value of a such that the two solutions reasonably agree, and the water content throughout the membrane.
Once we have the water content throughout the membrane, we can determine the resistance of the membrane and finally the voltage across the cell. Cells I30:AR30 calculate the conductivity of the membrane (as a function of the water content). This is then integrated in cells I31:AR31 to obtain the resistance of the membrane in cell I33. In other words, the numerical integration (essentially a summation) of the resistance of small portions of the membrane are added together to give the total resistance of the membrane. From this the voltage across the cell is obtained in cell I34. By systematically varying the current density in cell B8, and solving the above simulation to obtain the voltage across the cell, the polarization curve can be obtained. The polarization cell, obtained from the spreadsheet model, is depicted in Figure 7.11. For high currents, the concentration losses play a more important role and these are not taken into consideration here. Therefore the model is more accurate for regions of the plot with low to intermediate current. The curve, however, is just as we might expect, and this model has been successfully fitted to experimental data, validating the model to some degree.
