Steady-State and Dynamic Modeling of Intermediate-Temperature Protonic Ceramic Fuel Cells

Protonic ceramic fuel cells employ material sets, such as yttrium doped barium zirconate (BZY) or BZY-doped with cerium and ytterbium (BCZYYb), which exhibit conduction of multiple charged defects and are, therefore, fundamentally different than oxygen-ion conducting SOFCs or MCFCs. The BZY perovskite material exhibits both non-stoichiometry and hydration where oxygen transport is possible via oxygen vacancies and proton conductivity occurs due to protonic defects. The material also exhibits a small amount of electronic conductivity (via small polarons) in which fuel and air can be consumed without electrons being conducted through the external circuit. The anode is constructed from a BZY-Ni cermet which promotes separate electronic and ionic conducting phases and charge transfer reactions at the three phase boundary (TPB) as shown in Fig. 1. The nickel surface can participate in the heterogeneous steam methane reforming and water-gas shift reactions such that operation on a humidified methane fuel source is possible. Cathode materials have recently received increased attention and compositions giving rise to 'triple conductivity' have been tested which substantially reduce cell polarization and enable high power densities at practical voltage and current density conditions.¹² This model considers a cobalt and iron-doped BZY cathode material, however, double perovskites¹⁶ have been investigated as well.

Zoom In Zoom Out Reset image size

The transport of multiple-charged defect species during both operation and at open circuit conditions complicates the electrochemical analysis of the device.¹⁷ Recently developed membrane models have tried to resolve the highly coupled transport processes, however, the models have yet to be coupled with electrodes to predict the open circuit voltage and current voltage characteristics. Typically, PCFCs exhibit an open circuit voltage (OCV) which is lower than the theoretical Nernst potential due to the dynamic nature of equilibrium where thermodynamic and transport properties influence the result. Additionally, protons as well as oxygen ions can contribute to the charge transfer which can result in the formation of water at both the anode and cathode. The reactions associated with proton transfer are given as Eqs. 1 and 2 below. PCFCs also exhibit a small amount of oxygen ion conduction, and the resulting electrochemical reactions for the mobile oxygen ion are given as Eqs. 3 and 4.

The overall global oxidation reaction () is permissable at either anode or cathode electrodes for mixed ion conductors. However, since the calculation of a Nernst potential can proceed using the water partial pressure in either the anode or cathode channels, a unique value for OCV cannot be established given the variation in these respective water concentrations. Thus, any modeling effort necessarily needs to provide a means for accurate estimation of the open circuit potential for cells involving appreciable mixed ionic conduction.

Even considering the complicated nature of PCFCs, significant efforts have already been undertaken in modeling devices of this nature. However, the majority of the models are restrictive due to one or more of the following limitations: only consideration of proton transport in the electrolyte; operation on only a hydrogen fuel source; assumption of either isothermal operation or a lumped analysis; and neglect of OCV reduction due to non-equilibrium conditions at open circuit. All of these limitations become important when considering an engineering model to be integrated into system-level analysis.

Historically, reduced-order modeling of SOFCs in either one or two dimensions has been thoroughly investigated and reported in the literature.^18–24 Both steady-state and dynamic models have been developed and verified against experimental results for both planar and tubular devices. In this study, the modeling methodologies commonly employed for SOFCs are extended to the physicochemical phenomena unique to PCFCs in order to develop a planar, interface charge-transfer channel model that can represent the expected performance of a single cell.

Several modeling efforts for PCFCs have been undertaken,^25–30 however, these models have only considered high temperature cells (>700°C) with exclusively proton flux in the electrolyte. Additionally, several of the studies only consider a lumped geometry while others only consider a hydrogen fuel source. Recently several BZY20 membrane modeling studies^17,31 have considered multicomponent charge transfer; however, electrodes have yet to be modeled where equilibrium boundary conditions are employed and a resulting full PCFC model based on first principles is posed. Overall, practical engineering models of PCFCs which enable cell operating parameter selection, system design, and performance estimation are needed.

Prior modeling efforts from the extant literature which provide context to the present work are highlighted here. Ni et al. has presented several modeling studies^26–29 on both 1-D and 2-D channel-level models for planar PCFC devices. The majority of this work focused on the modeling of concentration overpotentials with proton conduction only and comparisons to SOFC devices. More recent work has focused on modeling hydrocarbon fuel sources and the effect of internal reforming reactions on PCFCs, however, these studies have all been at higher temperatures with only proton transport. The more rigorous problems of charge transfer and mixed conduction have yet to be addressed.

Patcharavorachot et al.³⁰ developed an electrochemical model for a cell operating on a hydrogen fuel source; the model was verified against high temperature button cell data. The model considered proton transport across the electrolyte and a global form of the Butler-Volmer equation. Concentration overpotentials were captured through a Fickian model, which was stated to be a better approach due to the explicit nature of the formulation in comparison to the dusty gas model. In their work the cell geometry was varied through the support structure of anode, cathode, and electrolyte, concluding that the anode supported design is the best. Dumortier (2012) et al. used commercial CFD tools to model a proton-conducting electrolysis cell which considered gas transport along the length of a channel and through the porous electrodes, however, thermal effects and carbonaceous byproducts were not considered. The model development focused on the charge and mass transport through the porous electrode materials and a parametric study was conducted to determine the effect of porous media design parameters on cell performance. Zheng et al.³² developed an electrochemical modeling approach which considered both proton and oxygen ion transport through the electrolyte which was termed 'co-ionic'. The model incorporated two Nernst potential calculations which were reduced by individual activation overpotentials, ohmic overpotential, and concentration overpotential. The contribution of the current density from protons and oxygen ions was determined by requiring that the potential of each electrochemical reaction was equal. The model was solved by a unique, iterative approach that varies proton current and oxide ion current until the net potential of each electrochemical reaction was equal to the desired cell operating voltage. A concern with this approach is defining and measuring exchange current densities for the electrodes that conduct multiple ions as well as electrons. Furthermore, Zheng et al. only used a single ionic conductivity for the co-ionic membrane when quantifying ohmic polarization in the fuel cell, even though multiple ions are mobile. Another concern is that the model does not take into account the loss of potential due to multiple charged defects being mobile under a zero current condition, and thus results in an over prediction of open-circuit potential.

Huang et al.³³ used fundamental charged defect chemistry to derive an expression for the open-circuit voltage (OCV) of a co-ionic SOFC. The mobile species considered were protons, oxide ions, and electrons. The model focused on the effect of multiple mobile species in the electrolyte and neglected losses that occur in the electrodes or the transport of reactants to the reaction zone. The resulting OCV becomes a function of the potentials of each electrochemical reaction (Eqs. 1–4) and the ionic conductivities of the two mobile ions in the membrane. The main difficulty with this approach is measuring or obtaining the individual ion conductivities in the electrolyte. There are multiple mobile charged species, so careful experiments must be performed to distinguish the conductivity of one ion from another in the electrolyte membrane.

More recently, Menon et al.³⁴ developed a 'quasi' 2-D channel-level model of a proton-conducting SOFC including internal reforming of methane to resolve temperature and species profiles over the length of the cell. However, the electrochemical model was based on pure proton conductivity and did not capture the deviation from experimentally measured OCV and theoretical Nernst potentials which led to overpredictions of the model voltage at low current densities. To correct for the over-predicted voltage at low current densities, a Butler-Volmer formulation for predicting activation overpotentials was implemented which led to concavity in the voltage-current model predictions not observed in the experimental data. Additionally, the simulations focused on higher temperature operation (>700°C) where internal reforming of methane is highly favored. Zhu and Kee³⁵ present a mixed ionic-electronic conductor model for a BZY10-type cell in a button cell configuration. Defect transport of protons, oxygen anions (vacancies), O-site polarons (~electrons), and lattice oxygen are captured. This model includes transport and reaction in the porous gas diffusion electrodes including on-anode reforming kinetics, mixed conduction in the membrane itself, and is a distributed charge transfer model that quantifies the electrostatic potential gradient across the membrane electrode assembly. However, while capturing oxygen anion transport, this model does not consider the possibility of electrochemical oxidation of hydrogen at the anode due to oxygen vacancy transport (anions) into the triple-phase boundary.

Unlike prior work, the present modeling effort is focused on an interface charge transfer model (i.e., no electrostatic potential gradients normal to gas flow) for intermediate-temperature ranges (∼450-600°C) that includes internal-reforming, multiple charged defect transport, and spatially resolves species and temperature profiles and their associated effect on ohmic and activation overpotentials in the gas streamwise direction. Electrochemistry within the model considers the electrochemical oxidation of hydrogen by oxygen anions in the anode-side triple-phase boundary as well. The electrochemical model is semi-empirical in nature, however, it has been formulated in such a way that allows for calibration from button-cell experimental data and rapid scale-up to cell-level model predictions. The resulting cell-level performance prediction accounts for the coupled effects of mixed proton and anion transport, reactant flow configuration, fuel depletion, non-isothermal behavior, and finite-rate reforming kinetics.

The modeling approach requires the solution of 1-D plug flow conservation equations for the anode and cathode gas channel mass, species, momentum, and energy balances. The equations have been written in their conservative form applying the product rule of differentiation to the temporal terms such that derivatives of state-space variables are isolated. The fuel and air gas channel mass balance equations are given by Eqs. 5 and 6, respectively. The equations capture the transient mass holdup, convective transport, and source term due to the electrochemical reactions.

where ρ is fuel or air density, is the gas velocity, M_i is the molecular weight of species i, and H is the height of the anode or cathode gas channel. The mass source term can be related to the local current density as,

where j_k is the current attributed to the k^th charge transfer reaction, F is Faraday's constant, and ν are the molar stoichiometric coefficients.

The species balance equations for the fuel and air gas channels are given by Eqs. 8 and 9, respectively. In addition to the terms captured in the mass balance equations, the mass species balance captures the consumption and production of mass species (y_i) via heterogeneous chemical reaction rates that are discussed in the fuel reforming section and in Eqs. 29−31.

The momentum balance equations for the fuel and air gas channels are not dynamic and are given by Eqs. 10 and 11, respectively. The momentum balance captures the pressure gradient (), wall shear stress (τ_w), and convective momentum transport. The shear stress experienced at the channel walls is determined from empirical correlations for rectangular ducts as given in Ref. 47.

Thermal energy equations for the fuel and air gas channels are necessary for resolving the temperature distribution. Terms associated with convective energy transport, heterogeneous reaction, electrochemical reaction, and gas-solid convection with the PEN and interconnect surfaces are included.

where λ_{pen, int} are the convective heat transfer coefficients between the respective gases and the solid surfaces of either the PEN or the interconnects of the cell, respectively. The terms associated with heterogeneous chemical reaction have been specified as heat fluxes which relate to the reaction rate and species enthalpy (h_{f, i}) as

The gas phase conservation equations must be coupled with the energy balances for the PEN and interconnect solid structures. Here a lumped analysis has been assumed where the PEN and interconnect temperatures are only resolved in the channel direction. The PEN energy balance captures terms associated with axial conduction, internal reforming, electrochemical oxidation and gas-solid convection via convective coefficients (λ_pen).

where t_pen is the thickness of the PEN repeating layer and t_pen = t_an + t_ca + t_el. The interconnect energy balance only captures the terms of axial conduction and gas-solid convection and is expressed as,

Performance of the electrochemical conversion device requires the channel-level conservation equations to be coupled to a submodel which is predictive of the current-voltage relationship of the device. Importantly, a fundamental difference between SOFC and PCFC modeling is that protonic ceramics exhibit continual charge defect transport even under zero electronic current (i.e., open circuit) conditions. The three mobile charged defects (protons, oxygen vacancies, and polarons) in the BZY electrolyte can exhibit a charge-balanced ionic current while still incurring no net electronic current³⁹. That is, oxygen and hydrogen are transported across the ceramic membrane even when the cell is not polarized. In an SOFC or PEMFC, there is only one mobile charge carrier in the electrolyte, so there can be no charged-balanced ionic current without there being a non-zero electronic current. Unlike SOFCs, PCFCs at open circuit are in a state of non-equilibrium, and as a result, the energy required to sustain defect transport effectively lowers the actual OCV below the theoretical value that is based on equilibrium thermodynamics.¹⁷ To properly capture this effect, an electrochemical model which has provisions for both reduced OCV and transport of multiple charge carriers is required. In particular, it is necessary to have transport properties (ion conductivities) to distinguish between the amount of charge transferred via oxygen ions and hydrogen protons. However, proton and anion conductivity in protonic ceramic materials is notoriously difficult to measure.

To address these challenges, the approach taken here leverages the work of Shimada et al.⁴⁰ who used the definition of the ionic transference number (t_ion) to reduce the cell OCV from the Nernst potential calculation of both global electrochemical reactions from Eqs. (1–4). The Nernst potential for each electrochemical oxidation reaction is calculated by Eqs. 17 and 18, where the subscript "O" means the term is in reference to Reactions 3 and 4 where the oxygen anion is transported through the membrane. The subscript "H" means the term is in reference to Reactions 1 and 2 where a proton is the charge carrier. The species partial pressures, p_{k, an} and p_{k, ca} are the partial pressures of the species at the fuel side triple phase boundary and air side triple phase boundary, respectively. The ionic transference number is the sum of the individual proton and oxide ion transport numbers (and ), and it can also be represented using Eq. 19. The true value of t_ion is inherently a function of proton and oxide ion conductivities in the electrolyte. However, open circuit button cell experiments can be used to empirically determine a t_ion for a given cell assembly. Eq. 20 adjusts the calculated Nernst potential associated with each global electrochemical reaction to a lower OCV value that takes into account the state of non-equilibrium at open circuit of PCFCs using t_ion. The Nernst potential of each oxidation reaction is weighted according to the amount of current contributed to the total current from each reaction. The parameter x_O/T represents the ratio of current coming from the transport of oxygen anions to the total current. It has been shown to be temperature dependent and ranges from 0.024 to 0.046 under reduced atmosphere at normal cell operating conditions.⁴¹ The fraction of current due to oxygen ion transport is calculated at any discrete position in the cell channel using Eq. 21.

At this point, the open circuit potential of the cell has been established after accounting for the losses due to sustaining the transport of multiple charge carriers. Next, much like a typical SOFC, the calculated OCV is reduced by the activation polarization and ohmic overpotentials through the area specific resistance (ASR), which yields the cell operating voltage at a given current density (j) using Eq. 22. EIS measurements at OCV with various gas compositions and temperatures allow for the extraction of ohmic and electrode polarization resistance parameters that depend on cell temperature and key partial pressures at the anode and cathode. The derived functional form of the ASR equation allows for both Arrhenius temperature dependence in the form of an activation energy in the exponential term and gas composition dependence in the pre-exponential term. Button cell experiments show that ohmic resistance of the cell is not dependent on gas composition, but the polarization resistance is dependent on the hydrogen pressure in the anode and the water pressure in the cathode. These are presented as Eqs. 23 and 24. The use of a semi-empirical electrochemical model relies on the calibration of several modeling parameters to experimental data which can be extracted from button cell experiments; this is discussed further in the button cell model calibration section. The cell operating voltage, calculated in Eq. 22, accounts for four losses occurring in the PCFC; multiple ion transport losses at open circuit with t_ion, concentration losses that are incorporated into the Nernst potential via Dusty Gas transport to the TPB, ohmic resistance losses, and activation (kinetic) losses.

Although the fundamental theory to capture transport in multiple charged defect conductors (MCDCs) has been developed,¹⁷ it relies on transport and thermodynamic parameters which are material specific and must be determined from experimental measurements.⁴² Coupled experimental approaches are needed, such as thermogravimetric analysis, and conductivity and conductivity relaxation experiments. Since the determination of thermodynamic and transport parameters for a given MCDC is a prohibitively rigorous task, the development of a semi-empirical method for the scale-up of button cell results to commercial-scale cell platforms provides an immediately useful tool while the basic composition of the cell is still being investigated. For example, the use of various sintering aids in the electrolyte material have shown increased conductivity leading to improvements in power density; however, as noted, above the observed OCV is reduced due to energy requirements in supporting the increased charged defect transport under net zero conditions.¹² As cell-level performance is typically determined by the characteristic of the polarization curve between OCV and maximum power density, the implication of different sintering aids on cell performance is not obvious. Thus, a rapid approach which uses button cell experimental results to perform cell-level studies would provide a useful engineering tool.

To predict the gas phase composition at the TPB, necessary for resolving the OCV and ASR variables, a transport model is used to capture the concentration overpotential of gas phase mass transport through the porous anode and cathode from the bulk channel gas flow. Since the gas phase mass transport from the bulk channel to the TPB does not occur by equimolar counter diffusion in the anode as it does in an SOFC, the general form of the Dusty-Gas model⁴³ is used to provide the functional relationship between the species concentration profile and molar flux. The relationship is given as,

where D^e_{k, Kn} is the effective Knudsen diffusion coefficient for species k, D^e_kl is the effective binary diffusion coefficient, B_g is the permeability of the anode/cathode, μ is the vapor viscosity of the fuel/air gas stream, is the driving force for diffusion, [X_k] and [X_l] are the molar concentrations of the species in the gas stream, [X_T] is the total concentration, and J_k and J_l are the molar fluxes of species k and l.

While the dusty gas model (DGM) is used to describe gas diffusion in the porous electrodes, being a constitutive relation, the DGM itself does not include a transient term. It should be clear that gas diffusion dynamics within the porous electrodes are, therefore, not captured within the present model. It has been shown that time constants associated with gas transport dynamics within porous electrodes, being on the order of 10⁻¹ s, are simply much faster than thermal transport dynamics (order 10² -10³ s).^{17,23,44–46} As shown in subsequent sections, the transient simulation studies reported here were conducted with current density ramps over 1s intervals; the effects of electrode gas dynamics will thus not be resolved due to the faster timescales for gas diffusion. Nevertheless, not including the transient effects of mass transport through the porous electrodes does not alter the predicted transient response of cells and stacks in terms of voltage, temperature, and power output, which are driven by the thermal response of the solid cell layers. Given that this interface charge-transport, channel model was developed for engineering-level simulations that can support experimental cell/stack testing, process systems design and analysis, and control strategy development, inclusion of dynamic mass transport within gas diffusion electrodes is an appropriate model simplification and increases computational efficiency. Applications for which including dynamic mass diffusion terms are necessary include modeling electrochemical impedance spectra and electrode microstructural design.^44,45

In addition to the electrochemical submodel, closure relationships for the friction factors, convection coefficients, gas properties, and chemical reaction rates are necessary. The friction factor and convection coefficients are taken from Ref. 47 for laminar flow in rectangular ducts. The physical properties of the gas phase are based on the ideal gas equations of state and use MultiFlash⁴⁸ to determine thermodynamic and transport properties assuming ideal mixing. The internal reforming reaction kinetics are based on the work of Timmermann et al.⁴⁹

The solution of the partial differential and algebraic equations is obtained numerically through implementation in gPROMS modelbuilder.⁴⁸ The transient conservation equations, which take the form of partial differential equations, have been spatially discretized through a finite volume method and the temporal terms are resolved through an ordinary differential equation time stepping algorithm (i.e., method-of-lines). Coupling the differential equations to the electrochemical submodel results in a system of differential-algebraic equations, which can be solved through the backward differentiation formula (BDF) method available in gPROMS.

The solution also requires the definition of boundary conditions to specify the inlet fuel and air composition and flowrate as well as coupling gas channel conservation equations and electrochemical submodel to the electrode conservation equations. The inlet mass flow rate of the fuel and air are specified by the fuel utilization (FU) and air ratio (AR), respectively, given,

is the molar flow rate of fuel supplied to the cell; I_cell is the current of the cell; and , and x_{f, CO} are the mole fractions of those components supplied to the fuel inlet. In calculating the air ratio, is the molar flow rate of air supplied to the cell and is the mole fraction of oxygen at the air inlet.

Determining operating conditions which do not result in solid carbon formation can be obtained through thermodynamic equilibrium calculations. A Gibbs ternary diagram in the C-H-O compositional space is a useful graphic for the display of these calculations and is given for different temperatures in Fig. 2. Carbon deposition analysis has been previously presented for the operating regime of SOFC devices by Sasaki et al.;^50,51 however, as noted above, operation of PCFCs is fundamentally different where protons are removed from the fuel channel rather than the addition of oxygen. This difference results in moving the fuel composition of PCFCs toward the solid carbon deposition boundary with increasing fuel utilization. This result is evident from the C-H-O ternary diagram (Fig. 2) where the solid carbon deposition boundary is plotted in compositional space for operating cell temperatures ranging from 400 to 600°C. The inlet fuel composition for humidified methane lies on the CH₄-H₂O line and removing protons will result in the composition moving directly away from the pure hydrogen vertex, which is toward the deposition boundary. However, crossing the boundary will only occur if fuel utilization is sufficiently high. Beginning with the same inlet composition for SOFC operation, the movement is toward the pure oxygen vertex and away from the deposition boundary, which is the reason only inlet compositions for SOFC operation need to be investigated.

Zoom In Zoom Out Reset image size

The calculations for the solid carbon deposition boundary implements the LaGrange Method of Undetermined Multipliers to minimize the Gibbs free energy of the gas composition when solid carbon is formed.⁵² The method requires specification of the gas phase species in the mixture and solves a set of algebraic equations to determine the composition. The analysis was performed using Engineering Equation Solver⁵³ and the internal ideal gas and graphitic carbon thermodynamic properties. The solid carbon deposition boundaries for three temperatures in the operating range of intermediate temperature PCFCs are plotted in Fig. 2, in which decreasing temperature is shown to have restrictive effects on the compositional space available for operation. In addition to solving for the solid carbon deposition boundary, the inlet S/C ratio needed to avoid deposition for a given fuel utilization can be determined. The required S/C at 500°C for a fuel utilization between 0 and 100% are also plotted below the ternary diagram in Fig. 2. Interestingly, the S/C ratio is shown to peak at a value of 2.325 at nearly 80% utilization where operation above this point will result in avoiding solid carbon deposition at all fuel utilizations at 500^oC. As fuel utilization is increased from 0 to almost 80%, the required S/C increases from about 1.3 to 2.3. Further increases in operating fuel utilization lower the required S/C ratio as the red curve in the lower portion of the figure illustrates. This lowering of the required S/C ratio corresponds to the point of inflection in the 500°C isotherm plotted on the ternary diagram. As fuel utilization approaches 100%, hydrogen is nearly completely removed from the anode gas and the gas mixture approaches pure carbon dioxide as shown by the dashed red-lines in both upper and lower diagrams. However, the required S/C ratio to avoid solid carbon deposition will vary with temperature.

Figure 3 plots the required S/C ratio to avoid solid carbon deposition over the entire temperature and fuel utilization space. The desired operating regime of PCFCs is highlighted within the figure for the temperature range of 500 ± 25°C and utilizations of 80−85% which are desired to maintain the benefits of low temperature and high electric efficiency. Ratios of S/C greater than about 2.4 are necessary to avoid solid carbon deposition in the desired temperature range for any utilization. It is important that the device avoids solid carbon deposition over all utilizations since all compositions in between will be observed over the length of the fuel channel. Additionally, concentration polarizations could lead to apparent local utilizations above the cell utilization such that solid carbon deposition can occur unless solid carbon deposition is suppressed at all utilizations, including for those above the cell utilization.

Zoom In Zoom Out Reset image size

In addition to humidifying the methane fuel source, a small amount of oxygen can be added to move the inlet composition away from the solid carbon deposition boundary. The addition of oxygen can serve as a substitute for a portion of the steam required to avoid solid carbon deposition, however, a portion of the methane fuel source will be directly oxidized and not available for electrochemical conversion. The required oxygen-to-carbon (O/C) ratio needed to avoid solid carbon deposition in temperature fuel utilization space at a fixed S/C of 1.8 are plotted in Fig. 3. The required O/C to avoid solid carbon deposition at any PCFC fuel utilization in the temperature operating regime is 1.1. An O/C of 4.0 would result in complete oxidation of the methane fuel to carbon dioxide and water leaving none of the methane able to be converted to hydrogen and electrochemically oxidized. Thus, the required O/C ratio of 1.1 would result in greater than 25% of the methane fuel source being oxidized to suppress solid carbon formation, which would greatly reduce system electrical efficiency. It should be further noted that the results presented in Figs. 2–3 are for pure hydrogen proton conduction in PCFCs, and do not include the 1−3% of total cell current that may be expected to occur via oxygen anion conduction in the 500−600°C operating temperature range. As oxygen conduction will mitigate carbon deposition, the results discussed thus far represent an upper bound to the required S/C ratios determined from the parameter space explored herein.

The electrochemical performance of high temperature fuel cells is typically improved by higher hydrogen content in the fuel supply gas. When the system is operating on a humidified methane fuel source, this would be achieved through some portion of external steam methane reforming. However, the amount of external steam-methane reformation is fundamentally limited by reactor operating temperature. Tight thermal integration between fuel processing and fuel cell stack generally translates into closely matched operating temperatures between these unit operations. The global fuel reforming reactions, steam methane reforming (SMR) and water-gas shift (WGS) reactions, and rate expression are provided in Eqs. 29−31. The global reaction rate for methane reforming over a nickel catalyst is developed in Ref. 49, and is used over a detailed micro-kinetic model because the global reaction model gives sufficiently accurate results while still remaining computationally efficient. Furthermore, the global reforming model requires fewer cell parameters, such as active catalyst area, anode porosity, or anode pore size to give an accurate representation of PCFC system performance.

In Eq. 31, k is the temperature dependent rate constant, p_i are the partial pressures of the relevant species, and K_eq is the global equilibrium constant for the SMR reaction. The WGS reaction is fast under normal PCFC operating conditions, so the rate is assumed to be infinite, and the reaction progress is then limited by equilibrium. The rate of hydrogen production from the WGS reaction is determined by the reaction equilibrium constant under local cell conditions.

Sensitivity studies on the reforming within the PCFC anode suggest that the conversion of methane is not a strong function of the rate parameter at normal operating conditions (i.e., high fuel utilization). Figure 4 below depicts the model-predicted methane concentration at the anode channel exit under two different fuel utilization conditions. The steam-to-carbon ratio was fixed at 2.4, leading to an inlet methane concentration of 29.1 mol%. The average cell temperature was varied so that temperature effects could also be analyzed. At low fuel utilization (5% FU), the steam reforming is kinetically limited as shown by the deviation of predicted methane concentration at the cell outlet (i.e., a form of conversion) versus the equilibrium amount. Doubling the reaction rate constant leads to an increase in internal reforming of methane, and a halving of the rate constant dramatically reduces it. Under these conditions, the hydrogen produced is substantially greater than that electrochemically consumed from the anode gas, thereby reducing the shifting equilibrium driving force along the cell to produce more hydrogen. Thus, the fuel channel is kinetically limited and therefore, sensitive to the kinetic rate parameters at low fuel utilization. However, at high utilization (80% FU)–a more likely operating point–the rate constant has little-to-no effect on the methane conversion in the anode fuel channel. Under high fuel utilization situations, a low hydrogen concentration occurs, which improves the equilibrium conversion, consequently enhancing methane reforming by maintaining a condition of shifting equilibrium along the channel. Under conditions of either doubled or halved reaction rate constants, the methane concentrations leaving the anode are very close to equilibrium across all potential PCFC operating temperatures. The authors conclude that the PCFC model is not highly sensitive to the rate of methane reforming at normal operating conditions due to the more dominant effect of shifting equilibrium in the channel rather than the value of reaction rate parameters.

Zoom In Zoom Out Reset image size

The equilibrium distribution of species for a S/C of 2.4 are plotted in Fig. 5 as a function of operating reactor temperature where the PCFC and SOFC temperature ranges have been denoted. Methane conversion is shown to approach 100% in the SOFC temperature range, however, in the PCFC temperature range methane conversion is limited to 30−40%. This result suggests that external methane reforming to improve cell performance is limited by the operating temperature range and only partial external reforming is possible. Indeed, given the catalytic properties of nickel-based PCFC anodes, attempts to supply a fully reformed steam-methane mixture at temperatures significantly higher than that of PCFC operation could easily be undone once introduced to the cell via the opposite of steam reforming (i.e., methanation).

Zoom In Zoom Out Reset image size

The calibration of the parameters needed for the semi-empirical electrochemistry model can be obtained from experimental button cell data. The model has been specifically formulated in such a way that the unknown modeling parameters can be easily extracted from button cell experimental data. Polarization curves as well as electrochemical impedance spectroscopy (EIS) measurements can be used to fit the ionic transference number and ASR parameters. The button cell preparation and experimental procedure are documented elsewhere¹² and only the parameters resulting from the calibration of the data set are given here.

The ionic transference number (t_ion) can be determined by the ratio of the experimental OCV to the theoretical Nernst potential. Plotting the button cell data as a function of temperature is shown in Figure 7 below. In order to impose the physical limitations of the transference number being between zero and one, the data points were fit to a sigmoid function (Eq. 32). The resulting parameters (β, T_o) are listed in Table I.

Zoom In Zoom Out Reset image size

EIS measurements were performed under OCV conditions and allowed for the separate measurement of polarization and ohmic resistances which have been reported for various anode and cathode gas compositions in Fig. 8. Plotting the resistance terms in Arrhenius form allows for the calculation of the activation energies and pre-exponentials from linear fits. The known gas composition variations in the anode and cathode compartments of the button cell have been varied such that the partial pressure dependence can be established. The button cell experiments showed that ohmic resistance in the cell is not dependent on gas partial pressures, but the polarization resistance is dependent the hydrogen partial pressure in the anode and the water partial pressure in the cathode. The resistances are shown to be dependent on key gas phase compositions in the anode and cathode as well as temperature, which have been built into the semi-empirical electrochemical model as given in Eqs. 23−25. All of the parameters obtained from the button cell testing for the semi-empirical electrochemical submodel are reported in Table I below.

Zoom In Zoom Out Reset image size

Button cell experiments are performed under a wide variety of test conditions in order to properly calibrate the modeling parameters, and then a button cell model is compared to button cell data that is collected under conditions mimicking potential operating conditions for validation. In this case, a button cell was tested at three temperatures, a steam-to-carbon ratio of 2.5, and current densities ranging from 0.0 to 0.6 Acm^{− 2}. The tests were performed in a temperature controlled furnace under flooded fuel and air conditions, so fuel depletion and temperature variation effects were neglected in the model. The results of the button cell test and model prediction are displayed in Fig. 9. The model error is calculated as the relative error between the experimental voltage and the model predicted voltage at each current density, and is then plotted below the polarization curve. The model prediction maintained an error of less than 5% at low current densities, and trended toward higher error at higher current densities at all temperatures tested; with the error at 600°C being the highest. Studies have shown that PCFC system optimal operating conditions usually fall in the 0.75–0.85 V range because this area reduces balance-of-plant size and cost as well as increases efficiency.⁵⁴ The model had low error in this voltage/current range, thus indicating that the model is sufficiently calibrated and ready to be scaled up to simulate larger systems than just a button cell.

Zoom In Zoom Out Reset image size

Simulation results for the 10 cm channel are displayed below. The PEN consists of a 0.5 mm Ni-BZY20 anode, 20 µm BZY20 electrolyte, and 40 µm doped BZY cathode. Species, temperature, current density, and velocity profiles are given for simulation results and are reported for co- and counter-flow reactant gas configurations. The simulations are performed at nominal 500°C, 80% fuel utilization, and an air ratio of 6.0. The average cell temperature of 500°C is maintained by varying the gas inlet temperatures, which will be different in co- and counter-flow configurations. The fuel composition is taken to be steam and methane at a S/C ratio of 2.4 with a small amount of hydrogen and the air composition is taken to be nitrogen and oxygen in the ratio of ambient air with the a small amount of water vapor which enhances the numerical stability of the electrochemical submodel.

The co-flow configuration modeling results displaying anode and cathode channel composition, temperature, and current density distributions are given in Fig. 10. The anode gas composition distribution over the length of the cell provides insight into the coupled processes of internal reforming and electrochemical oxidation of the hydrogen. The inlet compositions are shown to consist of humidified methane followed by equilibrium-limited reforming over the inlet region of the cell. The lower PCFC operating temperature of 500°C tends to promote the water-gas shift reaction but limits methane reforming. The gradual conversion of methane over the length of the cell is different than that typically experienced in higher temperature SOFC operation where equilibrium conversion to carbon monoxide and hydrogen is rapid and nearly complete in the first 30% of fuel gas transit in the streamwise direction over the cell. As previously noted, appreciable methane conversion in PCFCs is accomplished by a shifting equilibrium process whereby the hydrogen product is electrochemically oxidized, promoting the forward direction of the steam reforming reaction (Eqs. 29 and 30).

Zoom In Zoom Out Reset image size

Additionally, the increasing cell temperature over the length promotes additional methane reforming such that the highest hydrogen mole fractions are observed near the outlet. At 80% utilization, the exit methane and hydrogen concentrations are less than about 4% and 17%, respectively. Due to the equilibrium limitations of the internal reforming reactions at the low temperature inlet region of the cell, the hydrogen mole fraction in the anode channel never achieves a value above 30%. The effect of low hydrogen partial pressures on cell performance is two-fold: it increases ASR and decreases OCV–both of which contribute to the operating voltage and power density of the cell at a given current. The co-flow cell temperature distribution is shown to increase approximately linearly with streamwise position; a characteristic which arises from the competing effects of the exothermic processes of electrochemical oxidation and irreversible heat generation from i²R losses (i.e., Joule heating), and endothermic reforming reactions. The equilibrium constrained reforming reactions have a beneficial effect on the cell by reducing the temperature gradients in the PEN structure (∼1°C/mm) to values lower than that estimated for internally reforming SOFC devices (e.g., typically 1.5 to 2.5°C/mm). A full-scale PCFC in a co-flow configuration is estimated to obtain a single-cell nominal operating voltage of 0.846 V for the operating conditions of Fig. 10.

A PCFC in counter-flow configuration has a slight performance advantage over the co-flow operation as displayed in Fig. 10, obtaining an operating voltage of 0.849 V for the same operating conditions. Additionally, the counter-flow configuration has a significantly lower inlet gas temperatures at 706 K compared to 739 K for co-flow operation. This result is due to the recuperative effect of counter-flow cell operation which can maintain a nominal cell operating temperature of 500°C with lower gas inlet temperatures. The recuperative nature of the counter-flow cell operation also leads to a more uniform cell operating temperature and a lower peak temperature located near the cell inlet. The presence of the high temperature inlet region leads to increased internal reforming at the inlet and higher peak hydrogen concentration. The higher localized hydrogen mole fractions result in large variations in current density generating higher average local Nernst potentials and a decrease in cell ASR. At 80% utilization, the exit methane and hydrogen concentrations approach 7% and 11%, respectively. Although the counter-flow cell is shown to have better performance, the peak current density is 2.6% higher than the co-flow cell which could negatively impact cell degradation.

The 1-D channel model is also used to evaluate the local rate of the reforming reactions for a counter-flow configuration. As summarized in Fig. 11, methane reforming is quickly initiated at the channel inlet with a reaction rate nearly 4X that of the mildly exothermic WGS reaction. The rate of methane reforming rapidly reduces as the fuel gas proceeds along the channel. A similar trend occurs for the WGS reaction where there is a slight delay while carbon monoxide is produced in the steam-methane reforming reaction; but then after a rapid increase, the reaction rate decreases in the gas flow direction. After about 2.5 cm, the reaction rates for the SMR and WGS reactions are approximately equal and the fuel gas quickly approaches a condition of chemical equilibrium. This simulation was conducted at high fuel utilization, and the shifting equilibrium greatly contributes to favoring the fuel reforming reactions. The rate of electrochemical oxidation of hydrogen is also plotted and registers a nearly constant value with position at a rate about 3.5 to 4 times that of the reforming rates.

Zoom In Zoom Out Reset image size

Fig. 12a illustrates the magnitude of the rate of various heat source and sink terms in the PEN assembly energy balance as a function of position in the channel. All of the energy flows within the PEN assembly are shown in the control volume diagram of Fig. 12b. Positive values are associated with energy addition to the cell assembly and negative terms with energy rejected to either the bulk gas flows or leaving by conduction or power production. At the channel inlet, the endothermic reforming reactions and power production account for over 80% of the energy leaving the MEA, with electrochemical oxidation and some convective heat addition by the cathode gas balancing the energy flows. From Figs. 11 and 12, it is evident that as the rate of the reforming reactions decreases along the length of the channel, the local endotherm from the heterogeneous reactions decreases. Axial conduction is also significant at the channel inlet because larger temperature gradients in the PEN assembly are established from exposure to the cooler reactant gas flows. Other locations have minimal energy transfer via conduction because there is a minimal cell temperature gradient in the counter-flow configuration. This analysis illustrates several important aspects of PCFC operation, yet the choice of flow configuration in stack design will still be influenced by the overall system configuration and application (e.g., gas recycle and load-following). Nevertheless, the methodology developed here is capable of being implemented alongside balance-of-plant (BOP) components to simulate system-level performance.

Zoom In Zoom Out Reset image size

The dynamic response of the PCFC is of great importance when discussing its role in a distributed power generation application since high economic value is placed on dispatchable devices capable of large turndown and rapid load-following dynamics. The cell-only dynamic response for steps in current density under constant fuel utilization and air ratio control⁵⁵ is plotted in Fig. 13. The PCFC electrochemical response is observed to occur over a time scale of seconds, however, it is followed by a longer thermal transient due to the evolving temperature distribution in the cell. Overall, the transient response of the device is observed to be very rapid where low power density operation causes thermal transients of the device to occur over a much longer time scale. However, the temperature variations are shown to be slight (<7°C) for large load steps due to the low power density where the stack size will be approximately twice the size of a typical SOFC producing the same power. Additionally, the high heat capacity of the cell relative to the gas flow rates reduces the possibility of thermal shock due to the reduction of rapid temperature fluctuations and high temperature gradients. However, the dynamic response of a distributed generation device requires system-level simulation where the coupled interactions of the cell and BOP must be resolved. The model developed here is capable of being incorporated in such a system-level analysis, however, it is outside the scope of the present work.

Zoom In Zoom Out Reset image size