Multiple Transport Processes in Solid Oxide Fuel Cells Pei-Wen Li, Laura Schaefer and Minking K. Chyu Department of Mechanical Engineering University of Pittsburgh Pittsburgh, PA 15261, USA Abstract There are three important issues discussed in this topic concerning about the theoretical fundamentals and practical operation of a solid oxide fuel cell. The thermodynamic and electrochemistry fundamentals involved in a fuel cell are reviewed in the first section, which discusses the ideal efficiency and energy distribution, when a fuel cell converts the chemical energy in the oxidization of a fuel directly into electrical energy. Issues of chemical equilibrium for the solid oxide fuel cell with internal reforming and shift reactions, in case of methane or natural gas being used as the fuel, are also discussed in details in this section. The losses of electrical potential in the practical operation of a fuel cell is elucidated in the second section, which includes the discussions about the overpotentials from activation polarization, ohmic loss, and the loss due to mass transport resistance. In the third section, the coupling processes of flow, heat/mass transfer, chemical reaction, and electrochemical reaction, which influence the output performance of a fuel cell, are analyzed. The issues in a numerical computation for the fields of flow, temperature, and species concentration, which collectively determine the local and overall electromotive force in a solid oxide fuel cell, are discussed in details finally.
1 Introduction A fuel cell is a device that converts the chemical energy in a fuel oxidation reaction directly into electricity; it is substantially different from a conventional thermal power plant, where the fuel is oxidized in a combustion process and thermal-mechanical-electrical energy conversion process is employed. Therefore, unlike heat engines, fuel cells are not necessarily subject to the Carnot cycle efficiency limitation; its energy conversion efficiency is generally higher than
that of heat engines [1]. Ideally, the Gibbs free energy change in fuel oxidation is directly converted into electricity [1, 2] in a fuel cell. Common as any kinds of fuel cells, the core component of a solid oxide fuel cell (SOFC) is a thin gas-tight ion conducting electrolyte layer sandwiched by porous anode and cathode as shown in Figure.1; the special point is that this electrolyte is solid oxide material that only allows the passage of charge-carrying oxide ions. To produce useful electrical work, free electrons released in the oxidation of a fuel at the anode must travel to the cathode through an external load/circuit. Therefore, the electrolyte must has function to conduct ions while preventing electrons released at the anode from returning back to the cathode by the same rout. The oxide ions are driven across the electrolyte by the chemical potential difference on the two sides of the electrolyte due to the oxidation of fuel at the anode; such a difference of chemical potential is proportional to the electromotive force across the electrolyte, which, therefore, sets up a terminal voltage across the external load/circuit.
Figure 1: Principle of operation of an SOFC Only at high temperatures from 800-1000 oC, does the solid oxide electrolyte have sufficient ion conductivity. The high operating temperature of SOFC also ensures rapid fuel-side reaction kinetics without requiring expensive catalyst. In addition, the high temperature exhaust from SOFC can be directed to a gas turbine (GT); thus, using SOFC-GT hybrid system, one can achieve an efficiency of at least 66.3% based on low heat value (LHV, means the electrochemical product, water, is in the gas state) of the SOFC [3-6]. Since it operates via transport of oxide ions rather than that of fuel-derived ions, in principle, an SOFC can be used to oxidize a number of gaseous fuels. In particular, an SOFC can consume CO as well as hydrogen as its fuel, therefore it can be fueled with reformer gas containing a mix of CO and H2 [7, 8]. Since an SOFC operates under high temperature, its energy conversion efficiency and component safety are both highly concerned in industry. In the following sections, the issues to be discussed will include: (1) the thermodynamic and electrochemistry fundamentals concerning to the energy conversion and species variation, (2) the potential losses in practical operation, (3) the influence of flow, heat and mass transfer to the operational efficiency and safety, (4) numerical model simulating the output performance, the fields of flow, temperature, species concentration, etc.
2
Thermodynamics and electrochemistry fundamentals for solid oxide fuel cells
To study the energy conversion efficiency and energy distribution in the energy conversion process in a fuel cell, one has to understand the basic principles in the set-up of the chemical potential and thereof electromotive force across the electrolyte, which has multiple disciplinary of thermodynamic, electrochemistry, ion/electron conduction, and heat/mass transfer being involved in. In this section, the fundamentals of thermodynamics and electrochemistry for solid oxide fuel cell system are reviewed. Considering the isothermal oxidation of a fuel A with oxidant B, we express it by the following equation aA + bB + ⋅ ⋅ ⋅ ⋅ ⋅⋅ → xX + yY + ⋅ ⋅ ⋅ ⋅ ⋅ ⋅
(1)
The systematic change of enthalpy, Gibbs free energy, and entropy production in the reaction has the relationship of
ΔH = TΔS + ΔG
(2)
In a solid oxide fuel cell, the operating temperature is from 800 oC to 1000 oC, thus the gas species of reactants and products can be treated as ideal gases, which allows the chemical entropy change being expressed in the form of
ΔH = ( xh X + yhY + ⋅ ⋅ ⋅ ⋅ ⋅⋅ ) − ( ah A + bh B + ⋅ ⋅ ⋅ ⋅ ⋅⋅ )
(3)
where the h is the enthalpy. When a gas is pure, ideal, and at 1 atm, it is said to be in its standard state. The standard state is usually designated by writing a superscript 0 after the symbol of interest [9]. The Gibbs free energy pertains to one mole of a chemical species is called the chemical potential. For an ideal gas at temperature, T, and pressure of p, the chemical potential is expressed as g& = g& 0 + RT ln
p p0
(4)
where R is the gas constant; the p0 is the standard pressure of 1atm. One may omit the p0 in the denominator of the logarithm in Eq.(4), as it is 1 atm, however, in such a case, the pressure p must be measured in atm. The systematic change of Gibbs free energy in Eq.(1) can be expressed in terms of the standard state Gibbs free energy and the partial pressures of the reactants and products.
ΔG = ( xg& X + yg& Y + ⋅ ⋅ ⋅ ⋅ ⋅⋅ ) − ( ag& A + bg& B + ⋅ ⋅ ⋅ ⋅ ⋅⋅ )
= [ xg& X0 + yg&Y0 + ⋅ ⋅ ⋅ ⋅ ⋅⋅] − [ag& A0 + bg& B0 + ⋅ ⋅ ⋅ ⋅ ⋅⋅] + RT ln
= ΔG 0 + RT ln
( p X / p 0 ) x ( pY / p 0 ) y ⋅ ⋅ ⋅ ⋅ ⋅ ⋅ ( p A / p 0 ) a ( pB / p 0 ) b ⋅ ⋅ ⋅ ⋅ ⋅ ⋅
( p X / p 0 ) x ( pY / p 0 ) y ⋅ ⋅ ⋅ ⋅ ⋅ ⋅
(5)
( p A / p 0 ) a ( pB / p 0 )b ⋅⋅⋅⋅⋅⋅
where ΔG0 = ( xg& 0X + yg&Y0 + ⋅ ⋅ ⋅ ⋅ ⋅⋅ ) − ( ag& 0A + bg& B0 + ⋅ ⋅ ⋅ ⋅ ⋅⋅ )
(6)
which is the Gibbs free energy change of the standard reaction at temperature T (i.e., with each reactant supplied and each product removed at p0=1 standard atmosphere pressure). The theoretical electromotive force (EMF) induced from the chemical potential ( ΔG ) is the Nernst potential, which is in the form of E=
( p A / p 0 )a ( pB / p 0 )b ⋅⋅⋅⋅⋅⋅ − ΔG − ΔG 0 RT = + ln ne F ne F n e F ( p X / p 0 ) x ( pY / p 0 ) y ⋅ ⋅ ⋅ ⋅ ⋅ ⋅
(7)
where F (=96486.7 C/mol) is the Farady’s constant; the first part on the righthand side by the Gibbs free energy change of the standard reaction is also called as ideal potential, which is E0 =
− ΔG 0 ne F
(8)
where the ne is the number of electrons derived from per a molecules of the fuel, when the fuel is oxidized in the reaction of Eq.(1). While the Gibbs free energy change, − ΔG , converts to electrical power; the entropy production, −TΔS , is the thermal energy to be released in the process of electrochemical oxidation of the fuel. Both the h and g& 0 are solely functions of temperature for ideal gases, which are given in Table 1 for the gas species involved in the reactions of an SOFC. While the electromotive force in a fuel cell is determinable from the chemical potentials as discussed above; the current that a fuel cell takes out, denoted by I , is directly related to the molar consuming rate of fuel and oxidizer through the following expressions: m fuel =
where
I nefuel F
;
mO2 =
I
(9)
O2 e
n F
neO2 is for oxygen, which is the number of electrons per b molecules of
oxygen obtained in the electrochemical reaction; the n efuel is for fuel, which is the number of electrons derived from per a molecules of the fuel.
2.1 Operation on hydrogen fuel If an SOFC operates based on the fuel of hydrogen gas, the oxidation of hydrogen is the only electrochemical reaction in the fuel cell, which may be expressed as the following chemical equation H2 +
1 O2 = H 2 O( gas ) 2
(10)
The Nernst potential from this electrochemical reaction will be E( H 2 +1 / 2O2 = H 2O ) =
− ΔG(0H 2 +1 / 2O2 = H 2O ) 2F
+
RT [ln( pH 2 / p 0 ) anode + ln( pO2 / p 0 )0.5cathode 2F
− ln( p H 2O / p 0 ) anode ]
(11)
The ideal chemical potentials at temperature T(K) can be calculated from the data given by handbooks [10]. As a convenient reference, Table 1 gives the data of ideal chemical potentials and enthalpy for the gas species, which an SOFC will involve. Recognizing the fact of electrochemical equilibrium in the anode gas mixture, there is
ΔG = −ΔG(0H 2 +1 / 2O2 = H 2O ) + RT [ln( p H 2 / p 0 ) anode + ln( pO2 / p 0 ) 0.5 anode
− ln( p H 2O / p 0 ) anode ] = 0
(12)
Substitute Eq.(12) into Eq.(11), we obtain the Nernst potential in another form for the electrochemical reaction of Eq.(10)
E( H 2 +1 / 2O2 = H 2O ) =
RT [ln( pO2 / p 0 ) 0.5 cathode − ln( pO2 / p 0 ) 0.5 anode ] 2F
(13)
Since the oxygen partial pressure at the anode is very low (of the order of 10-22 bar) due to the anode reaction [1], it does not cause appreciable effect to the partial pressures of the major species in anode flow, which therefore gives sufficient accuracy when using Eq.(11) to determine the Nernst potential for the electrochemical reaction of Eq.(10). Representatively, for an electrochemical reaction by Eq.(1) hereafter, the common practice in determining the Nernst potential is to use Eq.(7), in which the partial pressure of oxygen on the cathode side and that of the fuel and product species on the anode side are used. The molar consumption rate of hydrogen and oxygen in the electrochemical reaction of Eq.(10) can be easily derived from Eq.(9) as: mH2 =
I ; 2F
mO2 =
I 4F
(14)
Table 1a: Enthalpy and standard state Gibbs free energy of species
Table 1b: Enthalpy and standard state Gibbs free energy of species
Table 1c: Change of enthalpy and standard state Gibbs free energy of reactions
2.2 Operation with methane reforming and shift involved internally It is necessary to have high operating temperature in a solid oxide fuel cell to maintain sufficient ionic conductivity for the solid oxide electrolyte [2, 4]. This gives a quite favorable environment for the reforming of hydrocarbon fuels, like methane. In fact, a solid oxide fuel cell operates based on the transport of oxide ion through the electrolyte layer from cathode side to anode side, the reforming products of hydrogen and carbon monoxide in fuel channel can both serve as fuel. Bearing this advantage, solid oxide fuel cells can directly utilize hydrocarbon fuels or, at least, the methane pre-reformed or partly reformed gas with components of CH4, CO, CO2, H2 and H2O is applicable. If this is the situation, fuel reforming and shift will happen in the fuel channel in a solid oxide fuel cell. The anode is, in fact, a good material serving as the catalyst for such chemical reactions at the high temperature in an SOFC, which therefore has no need of noble metals for catalyst [11]. If there are five gas species, CH4, CO, CO2, H2, and H2O, in the fuel channel, the solid oxide fuel cell will operate with internal reforming and shift reactions; therefore, we need to consider the electrochemical reaction and the coexisting chemical reactions of reforming and shift for determining the species’ molar fractions that are crucial to the electromotive forces in the fuel cell. Reforming:
CH 4 + H 2O ↔ CO + 3 H 2
(15)
Shift:
CO + H 2O ↔ CO2 + H 2
(16)
Since the high operating temperature of SOFC ensures rapid fuel reaction kinetics, it is a common practice to assume that the reforming and shift reactions are in chemical equilibrium [4] when determining the molar fractions of the species, which makes the computation significantly convenient. From the concept of chemical equilibrium, reactants and products species must satisfy the condition of ΔG = 0 . Therefore, the mole fractions or partial pressures of the five gas species in fuel stream will relate through the following two equations at the same time [12]:
pH 2
pCO ) 0 ΔGreforming p p0 K PR = = exp(− ) pCH 4 pH 2O RT ( 0 )( 0 ) p p pCO pH ( 0 2 )( 02 ) 0 ΔGshift p p K PS = ) = exp(− pH 2 O p RT ( CO )( ) p0 p0 (
0
)3 (
(17)
(18)
The dominant electrochemical reaction was reported to be the oxidation of H2 [11], which is primarily responsible to the setting-up of the electromotive force.
However, at the meantime, the electrochemical oxidation of the CO is also possible, and likely occurs to some extent in the solid oxide fuel cell. It has been reported that fuel cells operated by using mixtures of CO and CO2 have shown the electrochemical oxidation of CO to be an order of magnitude slower than that of hydrogen [13]. Nevertheless, there is no necessity to distinguish whether the electrochemical oxidation process involves H2 or CO in order to formulate the electromotive force. This is going to be discussed in the following paragraph. When the shift reaction of Eq.(16) in anodic gas is at the chemical equilibrium, there is 0 ΔG = [ g& CO + RT ln( 2
pCO2 p0
) + g& H0 2 + RT ln(
pH 2 p0
0 )] − [ g& CO + RT ln(
pH O pCO ) + g& H0 2 O + RT ln( 02 )] = 0 p0 p
(19) Rearrange this equation, there is 0 g& CO + RT ln( 2
p CO2 p0
0 ) − [ g& CO + RT ln(
p CO p0
)] = g& H0 2O + RT ln(
p H 2O p0
) − [ g& H0 2 + RT ln(
pH2 p0
)]
(20) /2 To substrate a term of [(1 / 2) g& O0 + RT ln( p O / p 0 )1cathode ] from both sides of Eq.(20), we have 2
2
pO / 2 p CO p 1 0 0 ) − RT ln( 02 )1cathode g& CO − g& CO − g& O0 2 + RT ln( 0 2 ) − RT ln( CO 0 2 2 p p p pH O pH pO / 2 1 = g& H0 2O − g& H0 2 − g& O0 2 + RT ln( 02 ) − RT ln( 02 ) − RT ln( 02 ) 1cathode 2 p p p
(21)
It is easy to find that the left-hand side of Eq.(21) is the Gibbs free energy change of the electrochemical oxidation of CO and the right-hand side is that for H2. Divided by (2F) on both sides, Eq.(21) is further reduced to E( H 2 +1 / 2O2 = H 2 O ) = E( CO +1 / 2O2 =CO2 )
(22)
where E( H 2 +1 / 2O2 = H 2 O ) is given in Eq.(11), while the Nernst potential from the electrochemical oxidation of CO is E (CO +1 / 2 O2 =CO2 ) =
− ΔG (0CO +1 / 2O2 =CO2 )
− ln( p CO2 / p 0 ) anode ]
2F
+
RT [ln( p CO / p 0 ) anode + ln( p O2 / p 0 ) 0.5 cathode 2F (23)
It is preferable that the EMF of an internal reforming SOFC is calculated from the electrochemical oxidation of H2; however, the species’ consumption and production are the results collectively determined from the reactions of
Eqs.(10), (15) and (16). For convenience, the mole flow rates of CH4, CO and H2 are denoted by their formula. Assuming that, x , y , and z are the mole flow rates respectively for CH4, CO, and H2 that are consumed in the three reactions by Eq.(15), Eq.(16) and Eq.(10) in the fuel channel, the coupled variations of the five species are in the following forms; when an interested section of fuel channel is considered from its inlet to outlet [8,14]:
CH4
out
= CH4 − x in
(24)
CO out = CO in + x − y
CO2 H2
out
out
(25)
= CO2 + y in
(26)
= H 2 + 3x + y − z in
(27)
H 2 O out = H 2 O in − x − y + z
(28)
The overall mole flow rate of fuel, denoted by Mf , will vary from the inlet to the outlet of the interested section of fuel channel in the form of
Mf
out
=Mf
in
+ 2x
(29)
Meanwhile, the partial pressures of species, proportional to its molar fraction, must satisfy the Eqs.(17) and (18) at the outlet of the interested section, which thus gives: 3
K PR
⎛ CO in + x − y ⎞⎛ H 2 in + 3 x + y − z ⎞ p 2 ⎜ ⎟⎜ ⎟ ( ) in ⎜ M in + 2 x ⎟⎜ ⎟ p0 M f + 2x f ⎝ ⎠ ⎝ ⎠ = in in ⎛ CH 4 − x ⎞⎛ H 2O − x − y + z ⎞ ⎜ ⎟⎜ ⎟ in ⎜ M in + 2 x ⎟⎜ ⎟ M f + 2x ⎝ f ⎠⎝ ⎠
K PS
⎛ H 2in + 3 x + y − z ⎞ ⎛ CO2in + y ⎞ ⎟ ⎜ ⎟⎜ in ⎜ ⎟ ⎜ M in + 2 x ⎟ M f + 2x f ⎝ ⎠ ⎝ ⎠ = ⎛ CO in + x − y ⎞⎛ H 2O in − x − y + z ⎞ ⎟ ⎜ ⎟⎜ in ⎟ ⎜ M in + 2 x ⎟⎜ M f + 2x f ⎠ ⎝ ⎠⎝
(30)
(31)
where the p is the overall pressure of the fuel flow in the interested section. Since, the electrochemical reaction is responsible by the oxidation of H2 as discussed in the preceding section, the consumption of hydrogen is directly related to the charge transfer rate, or ‘current’, I, across the electrolyte layer, which gives
z = I /( 2 F )
(32)
From the electrochemical reaction, the molar consumption of oxygen on cathode side will be calculated by using Eq.(14). By conducting a solution for the Eqs. (32), (31) and (30), we can determine the species variations, x , y and z . Finally, with the reacted mole numbers of CH4 and CO determined, the heat, absorbed from reforming reaction and released from shift reaction can be obtained, respectively as: Q Reforming = ΔH Reforming ⋅ x
(33)
QShift = ΔH Shift ⋅ y
(34)
Nevertheless, prior to the solution for the Eqs.(30) and (31), the electric current in Eq.(32) from the fuel cell must be known. This features the processes in SOFC a character of strong coupling that the species molar variation and the electromotive force, ion conduction as well as current flow are of interdependency. The ion transfer rate or current conduction in a SOFC will be discussed in the following section.
3 Electrical potential losses The ideal efficiency is never attained in practical operation for any fuel cells. In fact, there are three potential drops in a fuel cell that makes the actual output potential to be lower than the ideal electromotive forces of the electrochemical reaction. The nature of the fuel cell performance in response to loading condition can be realized by its polarization curve, typically as shown in Figure. 2. 1.2 Voltage Loss Due to Ohmic Resistance
Cell Potential (V)
1 0.8 Voltage Loss Due to Activation Resistance
0.6
Voltage Loss Due to Mass Transport Resistance
0.4 0.2 0 0
0.2
0.4
0.6
0.8
1
1.2
Current Density (A/cm2)
Figure 2: Over-potential in the operation of a fuel cell With an increase in current density, the cell potential experiences three kinds of potential losses due to different dominant resistances. The potential drop due
to the activation resistance, which is activation polarization, associates with the electrochemical reactions in the system. Another potential drop comes from the ohmic resistance in the fuel cell components, when ions and electrons conduct in electrolyte and electrodes, respectively. The third drop, which might be very sharp at high current densities, is attributable to the mass transport resistance, or concentration polarization, in the flow of fuel and oxidizer. It is known from observing the Nernst equation that the electromotive force of a fuel cell is a function of the temperature and the involved gas species’ partial pressures at the electrolyte/electrode interfaces, which are actually right proportional to their molar fractions. It is worth to note that, in the fuel stream, fuel must transport or diffuse from the core region of the stream to the anode surface, and also, the product of the electrochemical reaction must conversely transport or diffuse from the reaction site to the core region of the fuel flow. On the cathode side, oxygen must transport and diffuse from the core region of airflow to the cathode surface. Along with the fuel and air streams, consumption of reactants, or coming out of products will make the molar fractions of reactants to decrease and product to increase. Due to these resistances in the mass transport process, the feeding of reactants and removing of products to/from the reaction site only can proceed under large concentration gradient between the bulk flow and the electrode surface when the current density is high, which therefore induces the sharp drop of the fuel cell potential. As a consequence of all the above-mentioned potential drops, extra thermal energy will be generated as a heat release together with the theoretical heat of −TΔS from entropy production. The heat transfer issues in a solid oxide fuel cell will be considered later in section 4. 3.1 Activation polarization The activation polarization is the electronic barrier to be overcome prior to current and ion flow in the fuel cell. Chemical reactions, including electrochemical reactions, also involve energy barriers, which must be overcome by the reacting species. The activation polarization may also be viewed as the extra potential necessary to overcome the energy barrier of the rate-determining step of the reaction to a value such that the electrode reaction might proceeds at a desired rate [15]. The Butler-Volmer equation is the well-used expression for the activation polarization, η Act : ⎧ n Fη n Fη ⎤ ⎫ ⎡ i = i0 ⎨exp( β e Act ) − exp ⎢ − (1 − β ) e Act ⎥ ⎬ RT RT ⎦ ⎭ ⎣ ⎩
(35)
where β is usually 0.5 for the fuel cell application, which is the transfer coefficient; i is the actual current density in fuel cell, and io is the exchange current density. The transfer coefficient is considered to be the fraction of the change in polarization that leads to change in the reaction-rate constant. The
exchange current, io , is the forward and reverse electrode reaction rate at the equilibrium potential. A high exchange current density means that a high electrochemical reaction rate and good fuel cell performance can be expected. The ne in Eq.(35) is the electrons transferred per reaction, which is 2 for the reaction of Eq.(10). Substituting the value of β = 0.5 into the Eq.(35), one can obtain a new expression as follows ⎛ n Fη ⎞ i = 2i0 sinh⎜ e Act ⎟ ⎝ 2 RT ⎠
(36)
from which, the activation polarization is expressed as η Act =
⎛ i ⎞ 2 RT ⎟ sinh −1 ⎜⎜ ⎟ ne F ⎝ 2io ⎠
or
η Act =
⎞ i 2 RT ⎛⎜ i ln ( ) + ( ) 2 + 1 ⎟ ⎟ ne F ⎜⎝ 2i0 2io ⎠
(37)
For high activation polarization, Eq. (37) can be approximated as the following simple and well-known Tafel equation [15] η Act =
2 RT ⎛ i ⎞ ln ⎜ ⎟ ne F ⎜⎝ i0 ⎟⎠
(38)
On the other hand, if the activation polarization is small, Eq. (37) can be approximated as the linear current-potential expression [15] η Act =
2 RT i ne F i 0
(39)
Nevertheless, Eq. (37) is recommended for its integrity and accuracy to calculate the activation polarization. The value of exchange current density io is different for anode and cathode and is also depend on the electrochemical reaction, as well as the electrode materials. In the work by Chan et al. [15], io of 5300 A/m2 for the anode and 2000 A/m2 for the cathode for a SOFC were used; however there was no comment on setting these values in the paper. On the other hand, Keegan et al. [16] adjusted the io so as to obtain a simulation result to satisfy their experimental data; no report, however, is about the adjusted
io in their paper.
2
The present authors used a higher value of 6300 A/m and 3000 A/m2 [17, 18] respectively for the io of anode and cathode, and very good agreement between our numerical simulated cell terminal voltage and the experimental results from different researchers [19-22] has been obtained. Nevertheless, the io might be more or less different as reported in literatures for SOFCs.
3.2 Ohmic loss The ohmic loss comes from the electrical and contact resistance of electrodes and the current collecting components, as well as by the ionic conduction resistance of the electrolyte layer. Therefore, the conductivity of the materials for cell components and the current collecting pathway are the two factors most influential to the overall ohmic loss of an SOFC. In the state-of-the-art SOFC technology, lanthanum manganite suitably doped with alkaline and rare earth elements is used for cathode (air electrode) [21]; yttria stabilized zirconia (YSZ) has been most successfully employed for electrolyte; nickel/YSZ is applied over the electrolyte to form the anode. Temperature could affect the conductivity of SOFC materials more or less significantly; especially for the electrolyte, for example, the resistivity could be two orders smaller if its temperature increases from 700oC to 1000oC. The equations of resistivity for SOFC components suggested in literature are collected in Table 2. Table 2: Data and equations for resistivity of SOFC components
Bessete et al. [23] Ahmed et al. [24] Nagata et al. [25] Ferguson et al. [26]
Cathode Ω ⋅ cm
Electrolyte Ω ⋅ cm
0.008114 e500 / T
0.00294e10350 / T
*0.0014
0.3685 + 0.002838e
*0.1
10.0e[10092 (1 / T −1 / 1273) *0.013
T 4.2 × 105 o
Anode Ω ⋅ cm
0 .00298 e −1392 / T 10300 / T
1200 / T
e
1 e10300 / T 3.34 × 10 2
Interconnect Ω ⋅ cm
*0.0186 T 9.5 × 10 5
*0.5 *0.5
1150 / T
e
T 9.3 × 10 4
e1100 / T
* At temperature of 1000 C. A careful check for the equations in Table 2 was conducted. The expressions by Bessete et al. [23] was found reliable, which gives quite identical prediction as that by Ahmed et al. [24] and Nagata et al.[25]; the predicted data for anode resistivity by the equation of Ferguson et al.[26] has significant discrepancy with the predictions by other equations.
Figure 3: Schematic of a planar type SOFC
It is rational to assume that the passage of the charge-carrying species through the electrolyte, or the ion conduction through the electrolyte, is a charge transfer, like a current flow. In a planar type SOFC as shown in Figure 3, the current collect through the channel walls, or called ribs, after it perpendicularly across the electrolyte layer. The network circuit for current flow modeled by Iwata et al. [27] considers the channel walls as current collection pathway in a planar SOFC. However, the height and the width of the gas channel are both small (less than 3 mm), and the electric resistance through the channel wall might be negligible [24]. This simplification leads to the consideration that the current is almost perpendicularly collected, which means that the current flows normally to the trilayer of cathode, electrolyte and anode. When calculating the local current density, the ohmic loss is thus simply accounted for in the following way [24]: I = ΔA ⋅
a c ( EMF − η Act − η Act ) − V cell
(δ a ρ ea + δ e ρ ee + δ c ρ ec )
(40)
where ΔA is an unit area on the anode/electrolyte/cathode tri-layer, through which the current I passes; δ is the thickness of the individual layers; ρ e is for the resistivity of electrodes and electrolyte; Vcell is the cell terminal voltage; the denominator of the right-hand side term is the summation of the resistance of the tri-layer. The Joule heating due to current flow in the volume of ΔA × δ is expressed for anode, in the form of a QOhmic = I 2 ⋅ (δ a ρ ea / ΔA)
(41)
which is also applicable to electrolyte and cathode, by replacing the thickness and resistivity accordingly.
Figure 4: Schematic of a tubular SOFC In case the current pathway is relatively long in a fuel cell, for example, a tubular type SOFC, as shown in Figure 4, the current collects circumferentially,
which leads to a much longer pathway [28] compared to that of a planar type SOFC. In order to account for the ohmic loss and the Joule heating of the current flow in the circumferential pathway, a network circuit [17, 19, 29, 30] for current flow may be adopted as shown in Figure.5. Because the current collection is symmetric in peripheral direction in the cell components, only half of the tube shell is deployed with mesh in the analysis. The local current routing from the anode to cathode through electrolyte is determinable based on the local electromotive force, EMF , the local potentials in anode and cathode, as well as the ionic resistance of the electrolyte layer, which yields the expression: I =
EMF − η aAct − η cAct − (V c − V a ) Re
(42)
where V a and V c are the potentials in anode and cathode, respectively; the R e is the ionic resistance of the electrolyte layer in the thickness of δ e and a unit area of ΔA , which is in the form of R e = ρ ee ⋅ δ e / ΔA
(43)
where the ρ ee is the ionic resistivity of electrolyte.
Figure 5: Ion/electron conduction network in a tubular SOFC In order to obtain the local current across the electrolyte by using the Eq.(42), supplemental equations for V a and V c are necessary. Applying the Kirchhoff’s law of current to any grid located in anode, we obtain the equation associating the potentials of the interested grid P with the potentials of its neighboring points, east, west, north, south and the corresponding grid P in the cathode. (
VEa − VPa Rea
+
VWa − VPa Rwa
)+(
VNa − VPa Rna
+
VSa − VPa Rsa
)+[
P VPc − VPa − ( E P − η Act )
RPe
]=0
(44)
In the same way for a grid P in the cathode, there is (
V Ec − V Pc R ec
+
VWc − V Pc R wc
)+(
V Nc − V Pc R nc
+
V Sc − V Pc R sc
)+[
P V Pa − V Pc + ( E P − η Act )
R Pe
]=0
(45)
where R a and R c are the discretized resistance in anode and cathode respectively, which are determined according to the resistivity, the length of P current path and the area the current act on; the η Act is the total activation polarization, including both from anode side and cathode side. With all the equations for the discretized grids in both cathode and anode enclosed, a matrix represented by the pair of Eqs.(44) and (45) is formatted. When conducting a solution for such a matrix equation for potentials, the following approximations are useful: (1) at the two ends of the cell tube there is no longitudinal current flow, and therefore, insulation condition is applicable. (2) at the symmetric plane A-A as shown in Figure.4 and Figure.5, there is no peripheral current in the cathode and anode, unless the cathode or anode is in contact with nickel felt, through which the current flows in or out. (3) the potentials on the nickel felts are assumed uniform due to their high electric conductivities. (4) since the potential difference between the two nickel-felts is the cell terminal voltage, the potential at the nickel-felt in contact to the anode layer can be assumed as zero; thus the potential at the nickel felt in contact to the cathode will be the terminal voltage of the fuel cell. Once we obtained all the local electromotive forces, the only unknown condition for the equation matrix is either the total current flowing out from the cell or the potential at the nickel felt in contact to the cathode. This originates two approaches when predicting the performance of an SOFC. In case the total current taken out from the cell is prescribed as the initial condition, the terminal voltage can be predicted. On the other hand, one can prescribe the terminal voltage and predict the total current, i.e., the summation of local current I across the entire electrolyte layer. With the potentials being obtained in an electrode layer, the volumetric Joule heating in the electrode for an interested volume, controlled by P, will be q& Pa =
1 ⎡ ( VEa − VPa )2 ( VWa − VPa )2 ( VNa − VPa )2 ( VSa − VPa )2 ⎤ a a + + + ⎥ /( ΔxP ⋅ r ⋅ Δθ P ⋅ δ ) ⎢ 2 ⎢⎣ Rea Rwa Rna Rsa ⎥⎦
(46)
q& Pc =
1 ⎡ ( VEc − VPc )2 ( VWc − VPc )2 ( VNc − VPc )2 ( VSc − VPc )2 ⎤ c c + + + ⎢ ⎥ /( ΔxP ⋅ r ⋅ Δθ P ⋅ δ ) 2 ⎣⎢ Rec Rwc Rnc Rsc ⎦⎥
(47)
P ⎡ ( E − η Act − V Pc + V Pa ) 2 ⎤ e e q& Pe = ⎢ P ⎥ /( Δx P ⋅ r ⋅ Δθ P ⋅ δ ) R Pe ⎣⎢ ⎦⎥
(48)
where the r and δ with corresponding superscripts of a, c and e are the average radius and thickness respectively for the anode, cathode and electrolyte and, ΔxP and Δθ P are the P-controlled mesh size in axial and peripheral directions as shown in Figure.3. The volumetric heat induced from activation polarization in anode and cathode is in the form of: P ,a P ,a q& Act = I P ⋅η Act /( Δx P ⋅ r a ⋅ Δθ P ⋅ δ a )
(49)
P ,c P ,c q& Act = I P ⋅η Act /( Δx P ⋅ r c ⋅ Δθ P ⋅ δ c )
(50)
The thermodynamic heat generation occurring at the anode/electrolyte interface in the interested area, controlled by P, is: Q PR = ( ΔH − ΔG ) ⋅ I P /( 2 F )
(51)
3.3 Mass transport and concentration polarization Due to the gradual consumption, the reactants and oxidizer will decrease in their fractions respectively in the fuel and air streams, which will make the electromotive force to decrease gradually along the flow stream. On the other hand, due to the mass transport resistance, the concentration of gas species will have polarization in between the core flow region and the electrode surface, which will results in lower partial pressures for reactants but higher partial pressures for products at the electrode surfaces. Therefore, the fuel cell terminal voltage will be lower than the ideal value as is indicated by the Nernst equation. When at high cell current density, the increased requirement for feeding of reactants and removing of products can make the concentration polarization higher and thus the cell output potential will sharply decrease. In order to take the concentration polarization into account when calculating the electromotive force, the local partial pressures of reactants and products at the electrode surface are used. However, this needs the solution of the concentration fields for gas species in fuel and oxidizer channels, which might be simply based on a one-dimensional [31-33], or based on complicated two- or three-dimensional solution for mass conservation governing equations [17, 3436]. In fact, the concentration fields strongly coupled with the gas flow, temperature, and the distribution of the electromotive force in the ways as indicated in Fig. 6. First, the gas species mass fraction determines gas properties in the flow field; while the flow fields affect the gas species concentration distribution and temperature. Second, the gas species concentration field and temperature distribution determines the electromotive forces; while the ion/electron conduction due to the electromotive force determines the mass variation and heat generations in the fuel cell. The inter-dependency of these
parameters will be discussed in details in the following section when modeling an SOFC for predicting a fuel cell performance and the detailed distributions of temperature, gas species concentration, and flow field.
Figure 6: The inter-relation amongst concentration and other parameters
4
Introduction of computer modeling based on a tubular SOFC
An operation curve for SOFC, characterizing the average current density versus the terminal voltage, is very important when designing a hybrid system of SOFC/GT [37-40]. Other information like the temperature and concentration fields in an SOFC is also the high concerns for the safe operation of SOFC and downstream facilities in a hybrid system. Although there have been some experimental data about the operation performance and temperature of SOFC [ 19-22], an experimental test for SOFC is still rather tough because of its high operating temperature. This makes the modeling for SOFC very necessary. The purpose of computer simulation for an SOFC is to predict the operation characteristic in terms of average current density versus the terminal voltage, based on prescribed operation condition. The operation condition of an SOFC is solely fixed by fixing the flow rates and thermodynamic state of fuel and oxidant, as well as one of the following load conditions, namely, terminal voltage, taking-out current, and the external load [41]. The flow rates and thermodynamic conditions of fuel and oxidant may be called as internal conditions; the terminal voltage, taking-out current and external load may be called as external condition. Like any kinds of “battery”, the external load condition of SOFC determines the consumption of fuel/oxidizer and generation of products in the electrochemical reaction [42]; the only difference in a fuel cell
is its continuous feeding of fuel/oxidant and removing of products as well as waste species. According to the different way prescribing the external parameters, the following three schemes might be designed when constructing a numerical model for the SOFC in order to predict the other unknown parameters: (1) Use internal condition and terminal voltage to predict the taking-out current. (2) Use internal condition and taking-out current to predict the terminal voltage. (3) Use internal condition and external load to predict terminal voltage and current density. The cost of iteration computation using the three schemes is quite different. In the first scheme, the cell terminal voltage is known and thus the local current can be obtained, for example, by using Eq.(42) and solving Eqs.(44)-(45) respectively for a planar and tubular type SOFC, once the temperature and partial pressure fields of gas species being available. The integration of the local current will be the total taking-out current from the SOFC. In the second scheme, however, it needs to assume the terminal voltage and check the total current integrated from the local current until that the calculated total current agrees with the prescribed value; in this computation process, it needs a proper method to find the best-fit terminal voltage iteratively. The third scheme is resemble to the second scheme; one needs to assume a terminal voltage and get the total current; only when the voltage-current ratio equals to the prescribed load, can the computation be stopped. With the understanding of the principles of the energy conversion, chemical equilibrium, potential loss and the operation of an SOFC, we are now ready to construct a computer model for an SOFC. Since the manufacturing and stack operation for tubular SOFC is relatively mature at the present stage [21]. The following analysis takes tubular SOFC as an example, and there are three issues to be addressed to construct a numerical model. 4.1 Outline of a computation domain In a practical tubular SOFC stack, many tubular cells mount in a container to be a cell bundle, as shown in Figure 7. For most of the tubular cells in the bundle, everyone of them is surrounded by four others in the same style. Therefore, it is possible and very likely that most of the single tubular SOFCs work under the same environment of temperature and concentrations of gas species. This allows us to define a controllable domain in the cross-section, which pertains to one interested single cell, as outlined by the dashed-line square in Figure7. In the meantime, it is understandable that there must be no flow velocity and fluxes of heat and mass across the outline. This will significantly simplify the analysis for a cell stack. In other words, through analysis of the heat/mass transfer and the chemical/electrochemical performance for the single cell and its controllable area, one can obtain results very useful for evaluating the performance of a cell stack.
Figure 7: Schematic of a tubular SOFC in a cell stack With the longitudinal direction also being considered, the heat and mass transfer in the above outlined square area enclosing the tubular SOFC are in three dimensions, in nature. For a solution of the three-dimensional governing equations of momentum, energy, and species conservation, a large number of discretized mesh is necessary, which will make the computation too heavy to be conductible through a personal computer. In order to reduce computation cost, the square area enclosing the tubular SOFC is approximated to be an equivalent circular area; therefore, the domain enclosing the interested single tubular SOFC is viewed as a 2-dimensional axial symmetric one, as seen in Figure 7. However, the zero-flux, or insulation of heat and mass transfer at the boundary remains unchanged, even though the geometric approximation is applied to. From the preceding discussions, an axial-symmetrical two-dimensional (x-r) computation domain is profiled as shown in Figure 8, which includes two flow streams and solid area that locates the cell tube and air-inducing tube.
Figure 8: Computation domain for a tubular SOFC 4.2 Governing equations and boundary conditions As the mass fractions of species vary in the flow field, all the thermal and
transport properties of fluids are local functions of the species concentration, temperature, and pressure; therefore, the governing equations for momentum, energy, and species conservation (based on mass fraction) have variable thermal and transport properties being considered. ∂ ( ρu ) 1 ∂ (rρv) + =0 ∂x r ∂r
(52)
∂( ρuu ) 1 ∂(rρvu ) ∂p ∂ ∂u 1 ∂ ∂u ∂ ∂u 1 ∂ ∂v + =− + (μ ) + (rμ ) + ( μ ) + (rμ ) ∂x ∂x ∂x r ∂r ∂x r ∂r ∂r ∂x ∂x r ∂r ∂x
(53)
∂v 1 ∂ ∂v ∂ ∂u 1 ∂ ∂v 2μv ∂ ( ρuv) 1 ∂ (rρvv) ∂p ∂ (rμ ) + ( μ ) + (rμ ) − 2 + = − + (μ ) + ∂x r ∂r ∂r ∂x ∂x r ∂r ∂r ∂x ∂r r ∂r ∂r r
(54)
∂ ( ρCpuT ) 1 ∂ (rρCpvT ) ∂ 1 ∂ ∂T ∂T + = (λ ) + (rλ ) + q& ∂x ∂r r ∂x ∂x r ∂r ∂r
(55)
∂( ρuYJ )
(56)
∂x
+
∂Y ∂Y 1 ∂( rρvYJ ) ∂ 1 ∂ = ( ρD J , m J ) + ( rρD J ,m J ) + S m r ∂r ∂x ∂x r ∂r ∂r
These equations are applied universally to the entire computation domain; however, zero velocities will be assigned to solid area in the numerical computation. In the energy conservation equation, thermal energy from electrochemical reaction and Joule heating, represented by q& , are introduced as source terms in the proper locations in the fuel cell; some terms due to energy diffusion driven by the concentration diffusion of the gas species is very small and thus neglected [43, 44]. The boundary conditions for the momentum, heat and mass conservation equations are as follows: On the symmetrical axis, or at r = 0 : there is v = 0 and ∂φ / ∂r = 0 , where φ represents general variables except v . At the outmost boundary of r = r fo : there is thermally adiabatic conditions; impermeability for species and non-chemical reaction are assumed, which gives v = 0 and ∂φ / ∂r = 0 , where φ represents general variables except v . At x = 0 : the fuel inlet has prescribed uniform velocity, temperature, species mass fraction; the solid part has u = 0 , v = 0 , ∂T / ∂x = 0 ∂YJ / ∂x = 0 . At x = L : the air inlet has prescribed uniform velocity, temperature species mass fraction; the gas exit part has v = 0 , ∂u / ∂x = 0 , ∂T / ∂x = 0 ∂YJ / ∂x = 0 ; all the tube-end solid part has u = 0 , v = 0 , ∂T / ∂x = 0 ∂ YJ / ∂ x = 0 .
and and and and and
At the interfaces of air/solid, r=rair , and fuel/anode, r=rf : u = 0 is assumed. In the fuel flow passage, mass flow rate increases along x direction due to the transferring-in of oxide ion. Similarly, reduction of air flow rate occurs in the air flow passage, due to the ionization of oxygen and transferring of the oxide
ion to fuel side. Therefore, radial velocity exists at r=rair and r=rf , which are:
vf =
v air =
∑ m& xfuel ,species ρ xfuel
∑ m& xair ,species ρ xair
r =rf
(57)
r = rair
(58)
where m& [kg/(m2s)] is mass flux of gas species at the interface of electrodes and fluid, which arises from the electrochemical reaction in fuel cell. As important information, the mass fraction of all participating chemical components at the boundaries of r=rair and r=rf is calculated with consideration of both diffusion and convection effect [ 45, 46]:
m& xJ ,air = − D J ,air ρ xair
∂YJ + ρ xair YJ v air ∂r
m& xJ , fuel = − D J , fuel ρ xfuel
∂YJ + ρ xfuel YJ v f ∂r
(59)
(60)
It is worth to point out that the mass fluxes for species in the above equations, Eqs. (57)-(60), strongly relate to the ion/electron conduction; the determination of mass variation and related mass flux that arise from electrochemical reaction have been discussed in Section 2. As a consequence, the mass/mole fraction at the solid/fluid interface, derived from Eqs. (59) and (60), will be used for the determination of the partial pressures and thereof the local electromotive forces by Eq. (11). Table 3: Properties of SOFC materials Cp Density Thermal conductivity W/(m K) J/(kg K) kg/m3 d b a 11; c2.0; b2.0 623 4930 Cathode d c b b a 2.7; 2.7; 2.0 623 5710 Electrolyte c b a Anode 11.0; d6.0; b2.0 623 4460 c Support tube 1.0 Air-inducing tube c1.0 b b a Interconnector 13; c2.0; d6.0 800 6320; b7700 a b c Ahmed et al. [24]; Recknagle et al. [35]; Nagata et al. [25]; dIwata et al.[27]
The properties of solid materials in an SOFC are given in Table 3, which have some discrepancy from different literature. The single gas properties are available from reference [10] and [47]. For gas mixtures, equations from
reference [10] are selected for calculating the properties; they are listed in the following section. The mixing rule for viscosity is X i μi
n
μm = ∑
i =1
n
∑ X jφij
φij =
;
1 81 / 2
j =1
⎡ ⎛ M μ (1 + i ) −1 / 2 ⎢1 + ⎜ i ⎢ ⎜ μj Mj ⎣ ⎝
1/ 2
⎞ ⎟ ⎟ ⎠
⎛Mj ⎜ ⎜M ⎝ i
1/ 4 ⎤
⎞ ⎟ ⎟ ⎠
2
⎥ ⎥ ⎦
(61)
where μ m (Pa sec) is viscosity for mixture, and μ i (Pa sec) or μ j is viscosity of species (Pa sec); M i or M j is molecular weight of species; X i or X j is molar fraction; when i=j, φ ij = 1 . The mixing rule for thermal conductivity of gases at atmosphere pressure or less is: X i ki
n
km = ∑
i =1
n
∑ X j Aij
j =1
;
⎧ 1 ⎪ ⎡μ Aij = ⎨1 + ⎢ i 4 ⎪ ⎢μj ⎩ ⎣
⎛Mj ⎜ ⎜M ⎝ i
⎞ ⎟ ⎟ ⎠
3/ 4
⎛ T + Si ⎜ ⎜T + Sj ⎝
1/ 2 ⎫
⎞⎤ ⎟⎥ ⎟⎥ ⎠⎦
⎪ ⎬ ⎪ ⎭
2
⎛ T + Sij ⎜ ⎜T +S i ⎝
⎞ ⎟ (62) ⎟ ⎠
where k m [W/(m K)] and k i [W/(m K)] are thermal conductivities of mixture and species; Sij = C ( Si S j )1 / 2 , where C=1.0, but when either or both components i and j are very polar, C=0.73; for helium, hydrogen, and neon, Si or S j is 79K; otherwise, Si = 1.5Tbi , S j = 1.5Tbj , where Tb is boiling point temperature; the unit of T is K. When the gas mixture is above atmosphere pressure, the following correction is applied to the km obtained above: A × 10−4 (e Bρ r + C ) ⎛ Tc1 / 6 M 1 / 2 ⎞ 5 ⎜ ⎟Z ⎜ P 2/3 ⎟ c c ⎝ ⎠ ρr < 0.5 A=2.702 B=0.535 C=-1.000 0.5 < ρ r < 2.0 A=2.528 B=0.670 C=-1.069 2.0 < ρ r < 2.8 A=0.574 B=1.155 C=2.016 k = k' +
(63)
where k [W/(m K)] is the gas thermal conductivity at temperature T(K) and pressure P of interest of the mixture; k’ [W/(m K)] is the thermal conductivity at T and atmosphere pressure obtained by Eq.(62); ρ r = Vc / V is reduced density; Vc (m3/kmol) is critical molar volume; V ( m3/kmol) is molar volume at T and P; Tc (K) is critical temperature; M is molecular weight; Pc (M Pa) is critical pressure; Z c = PcVc /( RTc ) is critical compressibility factor, R is gas constant
which is 0.008314 [MPa m3 /(kmol K)] here. The mixture critical properties are obtained via the following equations: ⎡ ⎛n ⎞⎤ ⎡ Tcm − Tpc ⎤ Pcm = Ppc + Ppc ⎢5.808 + 4.93⎜⎜ ∑ X iωi ⎟⎟⎥ ⎢ ⎥ ⎝ i =1 ⎠⎦ ⎢⎣ Tpc ⎥⎦ ⎣
(64)
⎛ ⎞ ⎜ ⎟ n X jVcj Tcj ⎟ Tcm = ∑ ⎜ n ⎜ ⎟ j =1 ⎜ ∑ X iVci ⎟ ⎝ i =1 ⎠
(65)
Vcm = ∑ ∑ φiφ j vij (i ≠ j ) i
(66)
j
where φ j =
X jVcj2 / 3 n
∑ X iVci2 / 3
⎛ V −V ; vij = Vij (Vci + Vcj ) ; Vij = −1.4684⎜ ci cj ⎜ Vci + Vcj 2.0 ⎝
⎞ ⎟ + C ; C is zero ⎟ ⎠
i =1
for hydrocarbon systems and is 0.1559 for systems containing a nonhydrocarbon gas. In all the above equations, Xi or Xj is mole fraction of a species in the mixture; ωi is the acentric factor of a species; Pcm, Tcm and Vcm are mixture critical properties; Ppc and Tpc are pseudocritical properties of mixture, which are expressed as: n
n
i =1
i =1
Tpc = ∑ X iTci ; Ppc = ∑ X i Pci
(67)
Gas diffusivity of one species against the remaining species of a mixture is expressed in the form of
Dim
1− X i = ∑ X j / Dij j j ≠i
0.01013T 1.75 (
; D = ij
1 1 0.5 + ) Mi M j
P[(∑ vi )1 / 3 + (∑ v j )1 / 3 ] 2
(68)
where units of T, P and D are K, Pa, and m2/sec, respectively; Mi or Mj is molecular weight; all vi or vj are group contribution values for the subscript component summed over atoms, groups and structural features, which are listed in Table 4 for gases that an SOFC might involves. Table 4: Atomic diffusion volumes for use in Eq.(68) Atomic and structural diffusion-volume increments v [10] C 16.5 7.07 H2 17.9 H 1.98 N2 16.6 O2 O 5.481 CO 18.9 N 5.69 26.9 CO2 Aromatic ring -20.2 12.7 H2O Heterocyclic ring -20.2
4.3 Numerical computation
In order to conduct a numerical computation for flow, temperature, and concentration fields in an SOFC, a mesh system with sufficient grid number both in r and x direction must be deployed to the computation domain. All the governing equations may be discretized by using the finite volume approach, and the SIMPLE algorithm can be adopted to treat the coupling of the velocity and pressure fields [48, 49]. The temperature difference between the cell tube and the air-inducing tube might be large enough to have radiation heat exchange; therefore, numerical treatment based on the method introduced in literature [50] can be used to consider the radiation heat exchange. As has been discussed at the beginning of Section 4, the computation may base on the internal conditions and taking-out current; and as a consequence of simulation, the terminal voltage will be given out together with other interested details. The convenience of using this procedure in simulation is discussed next. It is quite common in practice that the taking-out current is prescribed in terms of the average current density of the fuel cell. Also, instead of the flow rates of fuel and air, the stoichiometric data are prescribed in terms of the utilization percentage of hydrogen and oxygen. This kind of designation of the operating conditions gives the convenience for comparing the fuel cell performance based on the same level of average current density and the hydrogen and oxygen utilization percentage. The inlet velocities of fuel and air are, then, obtainable in the forms of: u fuel = (
RT f Acell i cell ) 2 FU H 2 X H 2 A fuel P f
(69)
u air = (
A cell i cell RT air ) 4 FU O 2 X O 2 A air Pair
(70)
where Acell is the outside surface area of the fuel cell; A fuel and Aair are the cross-sectional inlet areas of the fuel and air; Pf , Pair and Tf , Tair are the inlet pressure and temperature of the fuel and air flows respectively; X H 2 and X O2 are the molar fractions of hydrogen in the fuel and oxygen in the air, respectively; U H 2 and U O2 are the utilization percentage for hydrogen and oxygen. The computation process is highly iterative and interdigitated, in nature. As the first step, the latest local temperature, pressure, species’ mass fraction are used in the network circuit analysis to obtain the cell terminal voltage, local current across the electrolyte and thus the local species’ transfer fluxes and local heat sources. In the second step, the local temperature, pressure and species’ mass fractions are, in turn, obtained through solution of the governing equations
under new boundary conditions determined by the latest-available species’ flux and heat sources. The two steps iterate interactively until convergence obtained. 4.4 Typical results from numerical computation for tubular SOFCs
The present authors have conducted numerical computation for three tubular SOFCs [17], which have been tested by Hagiwara et al [20], Hirano et al [19], Singhal [21], and Tomlins et al. [22]. The fuel tested by Hirano et al. [19] had components of H2, H2O, CO and CO2; therefore, there is water-shift of the carbon monoxide in the fuel cell to be considered together with the electrochemical reaction. The fuel used by other researchers [20-22] had components of H2 and H2O, and there is no chemical reaction except the electrochemical reaction in the fuel channel. The dimensions of the three different solid oxide fuel cells tested in their studies are summarized in Table 5, in which the mesh size adopted in our numerical computation is also given. Listed in Table 6 are the operating conditions, including species mole fraction and temperature of fuel and air in those test, which are the prescribed conditions for numerical computation. In the experimental work by Singhal [21], test of the pressure effect was also conducted by varying fuel and air pressure from 1 atm to 15 atm. It is expected that the experimental data for these SOFCs in different dimensions and operating conditions will facilitate a wide range benchmark for validation of the numerical modeling work. Table 5: Example SOFCs with test data available Data sequence: Outer diameter (mm ) / Thickness (mm) / Length (mm) Hagiwara et al. [20] Hirano et al. [19] Singhal [21] Tomlins et al. [22] 6.00 /1.00/290 12.00/1.00/1450 Air-inducing tube 7.00 /1.00 /485 13.00/1.50/300 Support tube 15.72/2.20/500 14.40/0.70/300 21.72/2.20/1500 Cathode 15.80/0.04/500 14.48/0.04/300 21.80/0.04/1500 Electrolyte 16.00/0.10/500 14.68/0.10/300 22.00/0.10/1500 Anode 18.10/ - /500 16.61/ - /300 24.87/ - /1500 Fuel boundary 66x602 66x1602 Grid number (rxx) 66x602 Table 6: Species’ mole fraction, utilization percentage, and temperature Air Fuel O2 %– U O 2 / N2 %/T(oC) H2 %– U H 2 /H2O%/CH4 %/CO%/CO2% /T(oC) I II
21.00-17.00/79.00/600.0 *21.00-25.00/79.00/600.0 **21.00-25.00/79.00/400.0 III 21.00-17.00/79.00/600.0
98.64-85.00/1.36 /0/ 0 / 0 / 900.0 55.70-80.00/27.70/0/10.80/5.80/800.0 55.70-80.00/27.70/0/10.80/5.80/800.0 98.64-85.00/1.36 /0/ 0 / 0 / 800.0
*Current density= 185 mA/cm2; ** Current density= 370 mA/cm2.
I- Tested by Hagiwara et al. [20] II- Tested by Hirano et al. [19] III- Tested by Singhal [21] and Tomlins et al. [22]
4.4.1 The SOFC terminal voltage The computer calculated and the experimentally obtained cell terminal voltages under different cell current densities are shown in Figure 9. The relative deviation of the model-predicted data from the experimental ones is no larger than 1.0% for the SOFC tested by Hirano et al. [19], 5.6% for that by Hagiwara et al. [20], and 6.0% for that by Tomlins et al. [22]. 1.0
0.5
Simulation Hirano et al. [19]
0.0
V
cell
1.0
0.5
Simulation Hagiwara et al.[20]
1.0
0.0
0.5
Simulation Tomlins et al.[22]
0.0 0
100
200
300
400
500
600
2
Current Density ( mA/cm )
Figure 9: Results of prediction and testing for cell voltage versus current density ( Operating pressure of the cell tested by Hagiwara et al. [20] and Hirano et al. [19] are 1.0 atm, and that by Tomlins et al. [22] is 5 atm ) It is interesting to observe from Fig.9 that, under the same cell current density, the cell voltage of the SOFC tested by Hagiwara et al.[20] is the highest and that by Hirano et al.[19] is the lowest. The mole fraction of hydrogen in the fuel for the SOFC tested by Hirano et al. [19] is low, which might be the major reason that this cell has the lowest cell voltage. Because the current must be collected circumferentially in a tubular type fuel cell, the large diameter of cell tube by Singhal [21] and Tomlins et al. [22] will lead to a longer current pathway; and thus the cell voltages of these cells are lower than that by Hagiwara et al.[20]; even though the former investigators tested the SOFCs at pressurized operation of 5 atm, which, in fact, helps to improve the cell voltage. Under a certain current density of 300 mA/cm2, the cell voltage and power increase with the increasing operating pressure, as seen in Fig.10. The agreement between our model-predicted results and the experimental one by Singhal [21] is quite good, showing a maximum deviation of 7.4% at a low operating pressure. When the operating pressure increases from 1 atm to 5 atm, the cell output power shows a significant improvement of 9%. However, raising the operating pressure becomes less effective for improving the output power, when the operating pressure is high; for example, the cell output power has only 6% increase when the operating pressure increases from 5atm to 15atm. The reason for this is that the operating pressure contributes to the cell voltage in logarithmic manner.
Nevertheless, pressurized operation of the fuel cell can improve the output power significantly, for example, increasing the operating pressure from 1atm to 15atm, the cell output power can have an increment of 15.8 percent. There is no doubt from the above investigation that investigators can predict the overall performance of an SOFC through numerical modeling and computation with satisfactory. On the basis of good agreement about the overall fuel cell performance, the internal details of flow, temperature, and concentration fields from numerical prediction is thus worth to refer. 1
250
0.8
200 Power (W)
V
cell
0.6 Current density =300 mA/cm2
0.4
Simulation Singhal [21]
0.2
150 Current density =300 mA/cm 2
100 Simulation Singhal [21]
50
0 0
2
4
6
8
10
12
14
0
16
0
2
Operating Pressure (atm)
4 6 8 10 12 14 Operating Pressure (atm)
(a)
16
(b)
Figure 10: Effect of operating pressure to the terminal voltage and power 4.4.2 Cell temperature distribution Because the measurement of temperature in SOFC is very difficult, only three experimental data, the temperature at two ends and in the middle of the cell tube, was available from the work by Hirano at al.[19]. Shown in Fig.11 is the simulated cell temperature distribution for the SOFC, for which Hirano et al.[19] provided three tested data. The agreement of the simulated data and the experimental results is good in the middle, where the hotspot locates; relatively larger deviation between the prediction and experiment appears at the two ends of the cell. Nevertheless, such a discrepancy is acceptable when guiding the designing of an SOFC concerning about the prevention of excessive hot of cell materials. o
Temperature of Cell Tube ( C )
1200 1000 800 600
Current density =185 mA/cm2
400
Simulation Hirano et al. [19]
200 0 0
0.1
0.2
0.3
x(m)
Figure 11: Longitudinal temperature distribution of fuel cell
150 mA/cm 100 mA/cm2 2
300 mA/cm
2
450 mA/cm
2
1200
300 mA/cm
2
400 mA/cm
2
500 mA/cm2
1100
o
1100
200 mA/cm
Temperature of Cell Tube ( C )
Temperature of Cell Tube ( oC )
1200
2
200 mA/cm2
1000 900 800 700 600 0
0.1
0.2
0.3
0.4
0.5
x(m)
1000 900 800 700 600 0
0.5
1
1.5
x (m)
(a) The cell Hagiwara et al. [20] reported
(b) The cell Tomlins et al. [22] reported
Figure 12: Predicted longitudinal temperature distribution for two SOFCs The predicted temperature distribution for the fuel cells tested by Hagiwara et al. [20] and Tomlins et al. [22] is given in Fig.12; unfortunately there was no experimental data on the cell temperature. Generally, the two ends of the cell tube have lower temperature than the middle of the cell tube. However, at low current densities, the hotspot locates closer to the closed-end of the cell. With the increase of current density, the hotspot shift to the open-end side and the hotspot temperature also decreases, which makes the uniformity of temperature distribution along the fuel cell to be better. It is worth to observe that the heat transfer between the cooling air and the cell tube at the cell closed-end region is dominated by a laminar jet impingement, because the exit velocity from the airinducing tube is quite low. However, the velocity of the exit air from the airinducing tube affects the heat transfer coefficient significantly. For the high current density case, the flow rate of air also becomes large accordingly; thus the heat transfer coefficient between the air and the fuel cell closed-end region is increased. This can suppress the temperature level of the closed-end region of the fuel cell significantly. Since the air receives much heat at the closed-end region, its cooling to the fuel cell at downstream region will becomes weak; thus the uniformity of cell temperature distribution becomes much better when the fuel cell operates at high current densities. 4.4.3 Flow, temperature and concentration fields Fig.13 is the flow and temperature fields for the SOFC tested by Hirano et al.[19] at a current density of 185 mA/cm2. The air speed in the air-inducing tube has slight acceleration because the air absorbs heat and expands in this flow passage. After leaving the air-inducing tube, the air impinges to the closed end of the fuel cell and then flows backward to outside. In this pathway, air obtains heat from the heat-generating fuel cell tube and transfers the heat to the cold air in the air-feeding tube. It is easy to understand that the electrochemical reaction at the closed-end of the fuel cell is strong because the concentrations of fuel and air are both high there. Therefore, the heat generation due to Joule-heating and entropy
change of the electrochemical reaction is high at the upstream of the fuel path. However, it is known from both experiment and computation that the fuel cell closed-end does not have the highest temperature; therefore, it is believed that the cooling of air to the closed-end of the fuel cell is responsible for this. After being heated at the closed-end region, air gets higher temperature and its cooling ability to the cell tube is low, when it is in the annulus between the air-inducing tube and the cell tube. At the cell open-end region, the air in the annulus can transfer heat to the incoming cold fresh air in the air-inducing tube and this will help it to cool the fuel cell tube. From this airflow arrangement, the hotspot temperature of the cell tube may mostly happen in the center region in the longitudinal direction of the cell tube. The airflow has two passes, incoming in air-inducing tube and outgoing in the annulus between air-inducing tube and cell tube; therefore, the heat exchange in between the two passes allows the air to mitigate its temperature fluctuation in the whole air path and thus the temperature field in the fuel cell might be maintained to be relatively uniform. Nevertheless, the heat generation, air and fuel temperature and air-cooling to the fuel cell will collectively affect the temperature field in the fuel cell. Therefore, the hot spot position in a cell tube might shift more or less away from the center region depending on the operating condition of a fuel cell.
r (m)
0.004 0.002 0
0
0.01
0.02
x (m)
0.03
0.04
0.05
r (m)
(a) The velocity vectors of air near the exit of the air-inducing tube 0.008
0.0075 0
0.1
x (m)
0.2
0.3
(b) The velocity vectors in fuel stream
0.006
848
0.004
861 935
0
960
911 87 3
92 3
97 3
8
0.002 0
998
948 94
r (m)
0.008
935
0.05
0.1
600
0.15 x (m)
0.2
0.25
0.3
(c) The temperature field (in oC) Figure 13: Predicted flow and temperature fields for the SOFC reported by Hirano et al. [19] at current density of 185 mA/cm2.
Figure 14 shows the gas species’ molar fraction contours in even difference for the same SOFC under the same operating condition as discussed in Fig.13. In the air path, oxygen consumption at the cell closed-end region is relatively large that leads to more densely distributed contour lines. The contour shape of oxygen also indicates a relatively larger difference of the molar fraction between the bulk flow and the wall of cathode/air interface. This also implies that the mass transport resistance on airside might be dominant in lowering the cell performance if the stoichiometry of oxygen is low. Feeding more air than needed is the way already well applied to in fuel cell technology.
r (m)
0.004
0.190
0.207
0.185
0.168
0.177
0.002 0.210
0
0
0.210
0.210
0.05
0.1
0.15 x (m)
0.2
0.25
0.3
r (m)
(a) Molar fraction of O2 in airflow
0.008
0.364 0.514 0.428 0.386 0.471
0.0075 0
0.05
0.321
0.1
0.300
0.278
0.15 x (m)
0.257
0.236
0.2
0.214 0.193
0.25
0.3
r (m)
(b) Molar fraction of H2 in fuel flow
0.008
0.409
0.321
0.463
0.498
0.533
0.569
0.586
0.622
0.640
0.0075 0
0.05
0.1
0.15 x (m)
0.2
0.25
0.3
r (m)
(c) Molar fraction of H2O in fuel flow
0.008
0.106
0.106
0.101 0.097
0.092
0.088
0.079
0.074
0.065
0.056 0.047
0.0075 0
0.05
0.1
0.15 x (m)
0.2
0.25
0.3
(d) Molar fraction of CO in fuel flow Figure 14: Predicted fields of molar fraction of species for the SOFC reported by Hirano et al. [19] at current density of 185 mA/cm2.
r (m)
0.008 0.060
0.0075 0
0.065 0.070
0.05
0.083 0.088 0.092 0.097
0.079
0.074
0.1
0.15 x (m)
0.2
0.106
0.25
0.115
0.3
(e) Molar fraction of CO2 in fuel flow Figure 14: continued
The consumption by electrochemical reaction and generation from water-shift of CO collectively determines the hydrogen budget. Since the consumption takes the upper hand, thus the hydrogen molar fraction decreases along the fuel stream. In correspond to this situation for hydrogen, the consumption due to water-shift and production due to electrochemical reaction make the water vapor increases gradually along fuel stream. The water-shift of CO proceeds gradually in the fuel path, and thus the mole fraction of CO decreases but the CO2 increases. The shape of contour lines of the species in fuel path is relatively flat from the cell wall to the bulk flow. This is the indication that mass diffusion in fuel channel is relatively strong than that in the airflow. H2 H2O
x 10
CO CO2
-5
Species Flow Rate (mol/s)
20
15 Current density =185 mA/cm 2
10
5 0
0
0.05
0.1
0.15
0.2
0.25
0.3
x (m)
Figure 15: Predicted streamwise molar flow rate variation for species in fuel channel for the SOFC reported by Hirano et al. [19] For further understanding of the variation of gas species, Fig.15 shows their molar flow rate variation along the fuel path. In one third of the length from fuel inlet, hydrogen flow rate shows faster decrease and water flow rate shows faster increase, indicating a stronger reaction at the upstream region. The flow rate of CO and CO2 vary roughly in a linear style and a small amount of CO still exists in the waste gas.
5 Concluding remarks Fuel cell technology is in rapid development currently. To improve the SOFC performance, for high power density and efficiency, people have to put effort to reduce the three over-potentials, activation polarization, ohmic loss, and concentration polarization. Better understanding of these three over-potentials is also very important for one to develop accurate computer model for predicting the overall performance and internal details in the SOFC. The activation polarization relates to the electrode porous structure and electrocatalyst materials. The state-of-the-art material and manufacturing process for electrodes and electrolyte has been reported by Singhal [21]. The reduction of ohmic loss also heavily relies on the reduction of electronic and ion resistance in electrodes and electrolyte; besides, a shorter current collection pathway also helps the reduction of the ohmic loss. A new design, referred to as high power density solid oxide fuel cell (HPD-SOFC), has been developed by Siemens Westinghouse Power Corporation [21, 28], which has a significant shorter current pathway; and thus the power density has been improved significantly. The planar structure is promising to have shorter current pathway and thus higher power density; measures for reducing the ohmic loss in a planar type SOFC has also been reported by Tanner and Virkar [51]. The reduction of mass transport resistance, or the concentration polarization, has not been given enough attention. In fact, the mass transfer enhancement has been reported effective for polymer electrolyte membrane fuel cell (PEMFC) to obtain higher cell current density [52] before the sharp drop of cell voltage, which is due to excessive concentration polarization. It might also be promising for SOFC to obtain higher current density by means of mass transfer enhancement. In a numerical model for an SOFC, the calculation of the over-potentials is very important to the accurate prediction of the overall current-voltage performance; also, the heat generation from over-potentials is decisive to the temperature, flow, and species concentration fields. About the activation polarization, studies elucidating the data and equation of the exchange current density are still needed. For the prediction of ohmic loss, reliable property data for electrodes are expected; also, method for analyzing a complex network circuit in an SOFC needs to be developed. The concentration polarization is considered in numerical computation by using the local molar fractions of species at the interface of electrode and fluid when calculating the electromotive force by the Nernst equation. Because the porous electrodes also serve as the reaction site, there is no well-described model for the mass transport resistance in the electrodes. Adopting a lower exchange current density, which induces larger over-potential of activation polarization, may be a way for considering the mass transport resistance in electrodes into the activation polarization. The method given by Hirano et al. [19] for the consideration of mass transport resistance in the electrodes gives convenience but is simple and may need more investigation. With the progress of computer modeling for SOFC, it is expected that cost for the research and development for SOFC will be significantly reduced by using computer simulation in the future.
References [1] Srinivasan, S., Mosdale, R., Stevens, P. and Yang, C., Fuel cells: reaching the era of clear and efficient power generation in the twenty-first century, Annual Review Energy Environment, Vol. 24, pp. 281-328, 1999. [2] Gardner, F.J., Thermodynamic processes in solid oxide and other fuel cells, Proc. Instn. Mech. Engrs., Vol. 211, Part A, pp. 367-380, 1997. [3] Suzuki, K., Iwai, H., Kim, J.H., Li, P.W. and Teshima, K., Solid oxide fuel cell and micro gas turbine hybrid cycle and related fluid flow and heat transfer, the 12th International Heat Transfer Conference, Grenoble, August 18-23, Vol. 1, pp. 403-414, 2002. [4] Harvey, S.P. and Richter, H.J., Gas turbine cycles with solid oxide fuel cells, Part I: Improved gas turbine power plant efficiency by use of recycled exhaust gases and fuel cell technology, ASME Journal of Energy Resources Technology, Vol. 116, pp. 305-311, 1994. [5] Layne, A., Samuelsen, S., Williams, M. and Hoffman, P., Developmental status of hybrids, International Gas Turbine & Aeroengine Congress & Exhibition, Indianapolis, Indiana, USA, No. 99-GT-400, 1999. [6] Ali, S.A., Moritz, R.R., The hybrid cycle: integration of turbomachinery with a fuel cell, International Gas Turbine & Aeroengine Congress & Exhibition, Indianapolis, Indiana, USA, No. 99-GT-361, 1999. [7] Bevc, F., Advances in solid oxide fuel cells and integrated power plants, Proc. Instn. Mech. Engrs., Vol. 211, part A, pp. 359-366, 1997. [8] Freni, S. and Maggio, G., Energy balance of different internal reforming MCFC configurations, International Journal of Energy Research, Vol. 21, pp. 253-264, 1997. [9] Wark, K., Thermodynamics, Second Edition, McGraw-Hill Book Company: New York, pp. 536-561, 1971. [10] Perry, R.H. and Green, D.W., Perry’s Chemical Engineer’s Handbook, Seventh Edition, McGraw-Hill: New York, 1986. [11] Onuma, S., Kaimai, A., Kawamura, K., Nigara, Y., Kawada, T., Mizusaki, J. and Tagawa, H., Influence of the coexisting gases on the electrochemical reaction rates between 873 and 1173 K in a CH4-H2O/Pt/YSZ system, Solid State Ionics, Vol. 132, pp. 309-331, 2000. [12] Li, P.W., Schaefer L. Wang, Q.M. and Chyu M.K., Computation of the conjugating heat transfer of fuel and oxidant separated by a heat-generating cell tube in a solid oxide fuel cell, Paper No. IMECE2002-32564, Proceedings of ASME International Mechanical Engineering Congress and Exposition, Nov. 17-22, New Orleans, Louisian, USA, 2002. [13] Hartvigsen, J., Khandkar, A. and Elangovan, S., Development of an SOFC stack performance map for natural gas operation, Proc. of the Sixth International Symposium on Solid Oxide Fuel Cells, Honolulu, Hawaii, USA, 1999; on Web site www.mtiresearch.com/pubs.html No. MTI00-07. [14] Massardo, A.F. and Lubelli, F., Internal reforming solid oxide fuel cell-gas turbine combined cycles (IRSOFC-GT), Part A: Cell model and cycle thermodynamic analysis, International Gas Turbine & Aeroengine Congress & Exhibition, Stockholm, Sweden, No. 98-GT-577, 1998. [15] Chan, S.H., Khor, K.A. and Xia, Z.T., A complete polarization model of a solid oxide fuel cell and its sensitivity to the change of cell component
thickness, Journal of Power Sources, Vol. 93, pp. 130-140, 2001. [16] Keegan, K., Khaleel, M., Chick, L., Recknagle, K., Simner, S. and Deibler, J., Analysis of a planar solid oxide fuel cell based automotive auxiliary power unit, SAE 2002 Would Congress, Detroit, Michigan, No. 2002-010413, 2002. [17] Li, P.W. and Chyu M.K., Simulation of the chemical/electrochemical reaction and heat/mass transfer for a tubular SOFC working in a stack, To be published in Journal of Power Sources (accepted on July 18, 2003). [18] Li, P.W., Schaefer L. and Chyu M.K., Investigation of the energy budget in an internal reforming tubular type solid oxide fuel cell through numerical computation, Paper No. IJPGC2003-40126, Proceedings of the International Joint Power Conference, June 16-19, Atlanta, Georgia, USA, 2003. [19] Hirano, A., Suzuki, M. and Ippommatsu, M., Evaluation of a new solid oxide fuel cell system by non-isothermal modeling, J. Electrochem. Soc., Vol. 139, No. 10, pp. 2744-2751, 1992. [20] Hagiwara, A., Michibata, H., Kimura, A., Jaszcar, M.P., Tomlins, G.W. and Veyo, S.E., Tubular solid oxide fuel cell life tests, Proceedings of The third International Fuel Cell Conference, Nagoya, Japan, pp. 365-368, 1999. [21] Singhal, S.C., Advances in solid oxide fuel cell technology, Solid State Ionics, Vol. 135, pp. 305-313, 2000. [22] Tomlins, G.W. and Jaszcar, M.P., Elevated pressure testing of the Simens Westinghouse tubular solid oxide fuel cell, Proceedings of The third International Fuel Cell Conference, Nagoya, Japan, pp. 369-372, 1999. [23] Bessette N.F., Modeling and simulation for solid oxide fuel cell power system, Georgia Institute of Technology, Ph.D Thesis, 1994. [24] Ahmed, S., McPheeters, C. and Kumar, R., Thermal-hydraulic model of a monolithic solid oxide fuel cell, J. Electrochem. Soc., Vol. 138, pp. 27122718, 1991. [25] Nagata, S., Momma, A., Kato, T. and Kasuga, Y., Numerical analysis of output characteristics of tubular SOFC with internal reformer, Journal of Power Resources, Vol. 101, pp. 60-71, 2001. [26] Ferguson, J.R., Fiard J.M., Herbin, R., Three-dimensional numerical simulation for various geometries of solid oxide fuel cells, Journal of Power Sources, Vol. 58, pp. 109-122, 1996. [27] Iwata, M., Hikosaka, T., Morita, M., Iwanari, T., Ito, K., Onda, K., Esaki, Y., Sakaki, Y. and Nagata, S., Performance analysis of planar-type unit SOFC considering current and temperature distributions, Solid State Ionics, Vol. 132, pp. 297-308, 2000. [28] Singhal, S.C., Progress in tubular solid oxide fuel cell technology, Electrochemical Society Proceedings, Vol. 99-19, pp. 40-50, 2001. [29] Sverdrup, E.F., Warde, C.J. and Eback, R. L., Design of high temperature solid-electrolyte fuel cell batteries for maximum power output per unit volume, Energy Conversion, Vol. 13, pp. 129-136, 1973. [30] Li, P.W., Suzuki, K., Komori, H., Kim, J.H., Numerical simulation of heat and mass transfer in a tubular solid oxide fuel cell, Thermal Science & Engineering, Vol. 9, No. 4, pp. 13-14, 2001. [31] Aguiar, P., Chadwick, D. and Kershenbaum, L., Modeling of an indirect internal reforming solid oxide fuel cell, Chemical Engineering Science, Vol.
57, pp. 1665-1677, 2002. [32] Bessette, N.F. and Wepfer W.J., A mathematical model of a tubular solid oxide fuel cell, Journal of Energy Resources Technology, Vol. 117, pp. 43-49, 1995. [33] Campanari, S., Thermodynamic model and parametric analysis of a tubular SOFC module, Journal of Power Sources, Vol. 92, pp. 26-34, 2001. [34] Haynes, C. and Wepfer, W.J., Characterizing heat transfer within a commercial-grade tubular solid oxide fuel cell for enhanced thermal management, Int. J. of Hydrogen Energy, Vol. 26, pp. 369-379, 2001. [35] Recknagle, K.P., Williford, R.E., Chick, L.A., Rector, D.R. and Khaleel, M.A., Three-dimensional thermo-fluid electrochemical modeling of planar SOFC stacks, Journal of Power Sources, Vol. 113, pp. 109-114, 2003. [36] Li, P.W., Schaefer L. and Chyu M.K., Interdigitated heat/mass transfer and chemical/electrochemical reactions in a planar type solid oxide fuel cell, Proceedings of ASME 2003 Summer Heat Transfer Conference, Las Vegas, Paper No. HT2003-47436, 2003. [37] Joon, K., Fuel cells – a 21st century power system, Journal of Power Sources, Vol. 71, pp. 12-18, 1998. [38] Watanabe, T., Fuel cell power system applications in Japan, Proc. Instn. Mech. Engrs, Vol. 211, Part A., pp. 113-119, 1997. [39] Veyo, S.E. and Lundberg, W.L., Solid oxide fuel cell power system cycles, Proceedings of the International Gas Turbine & Aeroengine Congress & Exhibition, Indianapolis, Indiana, USA, Paper No. 99-GT-365, 1999. [40] White, D.J., Hybrid gas turbine and fuel cell systems in perspective review, Proceedings of the International Gas Turbine & Aeroengine Congress & Exhibition, Indianapolis, Indiana, USA, Paper No. 99-GT-419, 1999. [41] Li, P.W., Schaefer L. Wang, Q.M. Zhang, T. and Chyu M.K., Multi-gas transportation and electrochemical performance of a polymer electrolyte fuel cell with complex flow channels, Journal of Power Sources, Vol. 115, pp. 90-100, 2003. [42] Cahoon, N.C. and Weise, G, ed., The Primary Battery, Vol. II, John Wiley & Sons, New York, 1976. [43] Bird, R.B., Stewart, W.E. and Lightfoot, E.N., Transport Phenomena, John Wiley & Sons, New York, 1960. [44] Williams, F.A., Combustion Theory, Benjamin/Cummings Publishing Co., 1985. [45] E.R.G. Eckert and R.M. Drake, Heat and Mass Transfer, 2nd Edition of Introduction to the Transfer of Heat and Mass, McGraw-Hill Book Company, Inc. 1966. [46] Turns, S.R., Introduction to Combustion: Concept and Application, 2nd Edition, McGraw-Hill Higher Education, New York, 1999. [47] F.P. Incropera and D.P. DeWitt, Introduction to Heat Transfer, Third edition, John Wiley & Sons, 1996. [48] S. V. Patankar, Numerical Heat Transfer and Fluid Flow, McGraw-Hill 1980. [49] Karki, K.C. and Patankar, S.H., Pressure based calculation procedure for viscous flows at all speed in arbitrary configurations, AIAA Journal, Vol. 27, No. 9, pp. 1167-1174, 1989.
[50] Beckermann, C., and Smith, T.F., Incorporation of internal surface radiation exchange in the finite-volume method, Numerical Heat Transfer, Part B, Vol. 23, pp. 127-133, 1993. [51] Tanner, C.W. and Virkar A.V., A simple model for interconnect design of planar solid oxide fuel cells, Journal of Power Sources, Vol. 113, pp. 44-56, 2003. [52] Nguyen T.V., A gas distributor design for proton-exchange-membrane fuel cells, Journal of Electrochemical Society, Vol.143, No.5, pp. L103-L105, 1996.
