Scholarly article on topic 'Sulfur Poisoning of SOFCs: A Model Based Explanation of Polarization Dependent Extent of Poisoning'

Sulfur Poisoning of SOFCs: A Model Based Explanation of Polarization Dependent Extent of Poisoning Academic research paper on "Chemical engineering"

Share paper
OECD Field of science

Academic research paper on topic "Sulfur Poisoning of SOFCs: A Model Based Explanation of Polarization Dependent Extent of Poisoning"

Sulfur Poisoning of SOFCs: A Model Based Explanation of Polarization Dependent Extent of Poisoning

Vinod M. Janardhanana,z and Dayadeep S. Mondera,b,z

aDepartment of Chemical Engineering, Indian Institute of Technology Hyderabad, Telangana 502205, India b Department of Energy Science and Engineering, Indian Institute of Technology Bombay, Maharashtra 400076, India

Several experimental studies have shown that, 1) the extent of the poisoning effect due to trace amounts of sulfur compounds in the fuel is lower if a SOFC is operated at a higher current density, and 2) the performance drop due to sulfur poisoning is much lower for Ni-GDC or Ni-ScSZ anodes when compared to Ni-YSZ anodes. This work presents a first principles numerical model that simulates experimental studies of sulfur poisoning on SOFC button cells. The exchange current densities for the electrodes are determined using sulfur-free polarization data for cells fueled by humidified mixtures of H2 and N2. A detailed surface reaction model that predicts the fractional coverage of all adsorbed species at the three phase interface is coupled to the SOFC model and the sulfur coverage is used to alter the anode exchange current density. The resulting model predictions match experimental observations during both galvanostatic and potentiostatic operation. Our analysis shows that the observed lower performance drop at higher current density is due to the non-linear nature of the electrochemical rate equations, and that the lower impact of sulfur poisoning on Ni-GDC and Ni-ScSZ anodes (compared to Ni-YSZ anodes) is due to their higher electrochemical activity.

© The Author(s) 2014. Published by ECS. This is an open access article distributed under the terms of the Creative Commons Attribution 4.0 License (CC BY,, which permits unrestricted reuse of the work in any medium, provided the original work is properly cited. [DOI: 10.1149/2.0611414jes] All rights reserved.

Manuscript submitted September 29, 2014; revised manuscript received October 1, 2014. Published October 10, 2014.

Solid oxide fuel cell (SOFC) systems operating on natural gas or gasified coal can be significantly more energy efficient than today's combustion systems and can thus reduce CO2 emissions. However, a major problem associated with the use of hydrocarbon fuels is the presence of sulfur containing compounds that are converted to H2S under SOFC operating conditions.1 The H2S present in the feed gas can readily deactivate conventional Ni cermet anodes. One can certainly have a desulfurization unit before the fuel cell stack to deal with the sulfur poisoning. However, this leads to increased system complexity and decrease in overall efficiency. Nevertheless, it is important to understand and quantify the sulfur tolerance of SOFC anodes in the event of sulfur breakthrough from the desulfurizer. Sulfur poisoning of Ni-YSZ anodes in SOFC is well studied experimentally for model fuel compositions consisting of ppm levels of H2S in mixtures of H2, H2O and an inert2-4 as well as H2S in gas mixtures representative of natural gas or reformed natural gas.5-8 Most of the experimental literature deals with H2S concentrations well below 10 ppm and work that treats H2S concentrations higher than 20 ppm is rather limited.2,5,8 On the other hand it is well known that natural gas and biogas can have substantially higher sulfur concentrations and that catalytic steam reforming of methane can be carried out on Ni even with 100 ppm H2S without losing all catalytic activity9,10 e.g., a methane conversion (in model biogas) of ~35% can be maintained using a supported Ni catalyst with 100 ppm H2S at 800°C.10

For an electrolyte supported cell Zha et al.,2 studied the cell performance for a wide range of H2S concentrations using DC polarization and electrochemical impedance spectra. The performance drops were reported for 0.18 ppm to 50 ppm H2S in H2/H2S mixture and they found the cell performance to be recoverable, although not fully, even with 50 ppm H2S. A different paper from the same group reported the impedance response of the cell before and after poisoning for galvanostatic and potentiostatic conditions and the maximum H2S concentration examined was 10 ppm.3 For a similar cell configuration as well as operating conditions Sasaki et al.,11 observed complete cell failure with 20 ppm H2S on Ni-YSZ. They found scandia stabilized zirconia (ScSZ) based anodes to be more sulfur tolerant compared to yttria stabilized zirconia (YSZ). Hagen et al.,5 reported similar findings with anode-supported cells and observed a lower degree of poisoning for Ni-ScSZ based anodes (compared to Ni-YSZ anodes). Rasmussen and Hagen also studied the performance of anode-supported cells for pulse addition of H2S.12 They observed an increasing drop in voltage upto 20 ppm H2S in galvanostatic mode. However, as in most other studies


they found that the performance drop leveled out once the H2S concentration was increased beyond 20 ppm. In early work on S poisoning of SOFCs, Matsuzaki and Yasuda studied cell performance at different temperatures and H2S concentrations.4 They used a thin anode (25 |xm) for their study and reported the variation in lower threshold H2S concentrations at which poisoning was observed for different operating temperatures. Hagen investigated the effect of step-wise addition of H2S for an anode supported cell.8 H2S concentrations of 2, 4, 10, 30, and 92 ppm were introduced and maintained for 24 h and the drop in voltage was monitored at galvanostatic conditions. Again, they observed that the cell performance drops sharply at low H2S concentrations and then levels out at the higher H2S concentrations.

Hansen and Rostrup-Nielsen,13 and more recently Cheng et al., wrote excellent reviews on sulfur poisoning of SOFC anodes in which they summarized the current understanding of this process from an experimental perspective considering Ni based anodes as well as other more sulfur-tolerant oxide based anodes. It is generally accepted that the degree of sulfur poisoning is lower at higher temperatures, and lower pH2S /pH2. Both these effects can be explained by assuming that the primary mechanism behind sulfur poisoning is the chemisorption of H2S followed by the desorption of the H atoms as H2 leaving behind sulfur adsorbed on the catalytic surface. Other effects such as a lower degree of poisoning at higher current densities14 as well as for Ni cermet anodes where yttria stabilized zirconia (YSZ) is replaced by mixed oxides of higher conductivity are not yet fully explained. As pointed out by Cheng et al.,3,1 a proper definition of degree of poisoning is lacking in the present literature and this leads to a certain amount of confusion when comparing experimental work on sulfur poisoning by different groups. The authors show how some of the confusion is avoided if the relative increase in polarization resistance is used as a measure of poisoning.

There are very few articles on the modeling of sulfur poisoning of SOFC anodes in the research literature. Ideally we need a model that can predict the sulfur coverage on Ni and relate it to the observed performance degradation under different operating conditions. Gazzarri and Kesler15 developed a mechanistic model for calculating the impedance spectra for an operating SOFC and investigated various degradation modes including sulfur poisoning. Their focus was identifying and deconvoluting the signatures of different degradation modes. For sulfur poisoning they showed that the increase in polarization resistance does not scale linearly with the electrochemically "active area" lost due to poisoning. For example, a 50% decrease in anode area leads to ~ 20% increase in anode polarization resistance in their simulations. They attribute this nonlinear effect to a rearrangement in the current production across the thickness of the anode.

However, their starting point is a given decrease in the active area and they do not try to model the relationship between fuel composition/temperature and the sulfur coverage, or the corresponding drop in "area" available for electrocatalysis.

Attempts to describe the effect of fuel composition/temperature on sulfur poisoned SOFC performance have focused on the Temkin isotherm first reported by Alstrup etal.,16 which was fitted to a large set of experimental data for sulfur uptake on Ni based catalysts. Hansen17 used this isotherm and showed that the observed performance drop was clearly correlated to the predicted sulfur coverage at the given inlet fuel composition and temperature. Hansen was able to fit linear expressions in the predicted sulfur coverage to the performance drop. However, each set of experiments yields a different expression. This is not surprising considering that the details about actual gas phase composition in the reaction zone or the electrochemical reaction mechanism and the effect of sulfur on that mechanism cannot be considered in such a correlation model.

A micron-scale model was presented by Monder and Karan18 for studying sulfur poisoning near the three phase boundary (TPB) in a SOFC anode. This model was built for a patterned anode geometry with Ni stripes on a YSZ substrate. The model includes a set of elementary reactions for gas surface reactions as well as the surface transport of, and reactions between, adsorbed species on both the Ni and YSZ surfaces. Also included is a set of electrochemical charge transfer reactions to directly simulate the dependence of anode overpotential on current produced. The model results showed that although the adsorbed sulfur coverage increases slightly in the near TPB-region when the anode is polarized, the deviation from the coverage in regions further away from TPB is negligible. Monder and Karan were also able to simulate the relatively lower performance drop seen at higher current density in the experimental data but no explanation was offered for this effect.

More recently, Prasad and Janardhanan19 developed a Langmuir-Hinshelwood type rate equation for the adsorption and dissociation of H2S on Ni. This allows the calculation of sulfur coverage on the Ni, which is coupled to an equivalent circuit electrochemical impedance spectra (EIS) model for an SOFC to simulate the effect of sulfur poisoning on EIS. This model was validated against the low ppm H2S performance data in Zha et al.2 and was shown to capture the expected increase in poisoning with decreasing temperature and increasing Ph2s/Ph2.

None of the handful of models available for sulfur poisoning of SOFC anodes is able to clearly explain either the lower degree of degradation at higher current densities or when using Ni cermet anodes where YSZ is replaced with higher ionic conductivity oxides. One of the possible explanations that have been cited in the experimental literature for the above behavior is that during SOFC operation, H2O produced as a result of electrochemical oxidation of H2 can lead to the formation of surface adsorbed elemental oxygen by dissociative adsorption. This oxygen can further react with adsorbed sulfur to clean up the three phase boundary (TPB). Another hypothesis is that the adsorbed sulfur in the TPB region is electrochemically oxidized by O2- ions at high current densities.14,13 The objective of the present work is to provide a theoretically sound model that can adequately explain the experimentally observed dependence of degradation due to sulfur poisoning in H2 based fuels on operating conditions and material properties.

Modeling Methodology

The button cell model considered in this work simulates the membrane electrode assembly (MEA) in the time domain. Mass transport in the gas channels above the electrodes is not considered in this work. The species transport in the MEA is governed by Eq. 1 while the density of the mixture is given by Eq. 2.

^Mk) = JA + AswkSk,

d t d y s k

k = 1,

9 (ip) sr-^djk ^ .

= -z. dy + AsWkSk ■

k=1 y k=1

In the above equations, e is the porosity, p is the density (kg m-3), Yk is the mass fraction of species k, t is the time (s), jk is the mass flux of species k (kg m-2 s-1), y is the spatial coordinate (m), As is the surface area available per unit volume (m-1), Wk is the molecular weight of species k (kg mol-1), sk is the molar production rate of gaseous species k due to surface reactions (mol m-2 s-1), and Ng is the total number of gas-phase species. The mass flux of species jk is calculated using the Dusty Gas model.20 The binary diffusion coefficients that are required for the evaluation of species mass flux are calculated using the kinetic theory of gases.21 The pressure p within the electrodes is calculated the from the ideal gas law (Eq. 3),

pM = p RT.

Here M is the average molecular weight (kg mol-1) of the gas mixture, T is the temperature (K) and R is the gas constant (J mol-1 K-1). The surface coverage 9k of species k adsorbed on the surface is calculated according to

d- = k = Ng + 1,..., Ng + Ns,

where ok is the number of sites occupied by the adsorbed by species k, sk is the molar production rate of surface adsorbed species k (mol m-2 s-1), Ns is the number of adsorbed surface species, and r is the total site density (mol m-2). The net molar production rate sk of a gaseous or surface adsorbed species due to heterogeneous reaction is calculated according to22-25

= J2 vkikfi H [fki

Here Rs is the number of surface reactions, [Xk] is the concentration of species k, kfi is the forward rate constant for reaction i, and vki is the difference in stoichiometric coefficient for species k between the products and reactants in reaction i. v'ki is the stoichiometric coefficient of species k in the forward reaction, which is equal to one for all species in the reaction mechanism considered here. The forward reaction rate kfi for reaction i is calculated by a modified Arrhenius

expression22 25

kfi = AiT V

Ei_ RT

kki exp

^ki Ok ~RT

Here Ai is the pre-exponential factor, P is the temperature exponent, E is the activation energy, 9k is the surface coverage of species k, |xki is a parameter that can be used to change the order dependency, and eki is the linear coverage dependent activation energy of reaction i on 9k (coverage of species k). The product in Eq. 6 runs only over those species which contribute to a surface coverage modification of the rate expression. For a reaction that is not modified by surface coverages the rate constant is simply

kfi = AiTV exp ^ - RT) ■

Electrochemical model.— Two different electrochemical models are implemented and compared in this work. The first model assumes that the charge transfer reactions occur only at the interface between the dense electrolyte and electrode layers, where as the second model allows charge transfer reactions to occur throughout the thickness of the composite electrode. The charge transfer reaction rate is expressed in the form of a Butler-Volmer like equation as derived in Zhu et al.20 The electrochemical rate equation for H2 oxidation takes the form given in (Eq. 8):

(1 + Va,an)F r|a RT

F na RT

whereas the rate equation for O2 reduction takes the form given in (Eq. 9):

ßc,an F RT

- exP -

ßc,caF nc

In the above equations, i is the current density (A cm-2), F is Faraday's constant (C mol-1), % and nc are the activation overpotentials (V) for H2 oxidation and O2 reduction reactions respectively, and pa>an and |a>ca are the symmetry factor for the anodic and cathodic directions of the anodic charge transfer reaction. The first subscript for the symmetry factors denotes the electrode (a for anode, c for cathode). Both sets of symmetry factors follow the relation |an + |ca = 1. The exchange current densities i0 and are lumped parameters that are functions of the three phase boundary (TPB) length, temperature, and the surface concentrations of species that take part in the electrochemical charge transfer reactions. The surface concentrations can be expressed in terms of the partial pressures of gas-phase species assuming that all elementary reactions except for the rate-limiting step are in quasi-equilibrium. The exchange current density for H2 oxidation takes the form

i0 = i

( PH / PH2fa'Ca'2( PH2O/p0)'

1 + ( PH2 / PH2)1/2


and for O2 reduction takes the form

ic0 = ic

( PO2 / PO2)ßc'W2 1 + ( PO2 / PO2 )1/2 "

Here i* and i* are fitting parameters that are functions of temperature, and TPB length. pH2, Ph2o, and pO2 are the partial pressures of H2, H2O, and O2 respectively at the reaction sites, and p0 is the standard pressure (1 bar). Expressions for pH , pO and pH O can be found in Zhu et al.20 i* and i* can be further expressed in the Arrhenius form according to Eqs. 12 and 13:

la = ia eXP

Ic = Ic eXP

In writing Eqs. 12 and 13 the initial TPB length at the anode and cathode is lumped into i*' and i*' respectively. The constants i*', EH2, i*', and EO2 are empirical parameters that are adjusted to reproduce experimental observations. The species flux at the electrode-electrolyte interface due to electrochemical charge transfer reactions is calculated according to

jk = ±

where jk is the species molar flux (moles cm-2 s-1) and ne>k is the number of electrons transferred in the charge transfer reaction involving species k.

Charge transfer at the interface between electrode and dense electrolyte (Model-I).—This model assumes that the charge transfer reactions occur only at the interface between the dense electrolyte and the electrodes.20 Several of the early SOFC electrochemistry models in literature are based on this assumption.26-32 In this approach, the current density i is related to the cell voltage Vcell according to Eq. 15.

Vcell = Erev - na — jîlc| - Rsi,

where Erev is the reversible potential given by the Nernst equation, na and nc are functions of i through Eqs. 8 and 9, and Rs is the series area-specific resistance (ASR in ^ cm2) that represents the sum of the resistance of the electrolyte and any contact resistances. In the current work, the ASR is determined from the given experimental EIS data.2

Charge transfer throughout the composite electrodes (Model-II).— In composite electrodes the electrochemical charge transfer occurs throughout the thickness of the electrode. In this approach the cell voltage is given by

Vcell = 9c - Va,

where and are respectively the cathode and anode current collector potentials. Charge conservation in each phase leads to the following equations.33,34

Electronic potential cathode:

d(Ve,c - Vi,c) ' dt

Ionic potential cathode:

d(Vi,c - Ve,c) ' dt

Electronic potential anode:

d(Ve,a - Vi,a)

Ionic potential anode:

d(Vi,a - Ve,a ) d t

,9 Ve,c ! d y

,9 Ve.c ' d y

,9 Vi,a

where Cdl is the effective double layer capacitance between the electron and ion conducting phases, aee and o^ are the effective conductivity of charged species in the electron and ion conducting phases respectively. Ve and ^ are the potentials of electron conducting and ion conducting phases and the second subscript denotes the electrode (a for anode, c for cathode). The volumetric current density i in the above equations is given by Eq. 8 in the anode and Eq. 9 in the cathode. It should be noted that in Model-II the parameters ia and ic are in A m-3. The local electrochemical overpotentials are then given by

nc = (Ve - Vi )c - EC1 and na = (Ve - Vi )a - Eeaq.

Here (Ve-ty ) is the difference in potential between the electron conducting phase and ion conducting phase, ECq is the equilibrium value for the above potential difference in the cathode, and E^q is the equilibrium value for the above potential difference in the anode. These equilibrium potentials are calculated using the Nernst equation and the local gas phase composition. Also,

Eeq - Eeq = Ere

where the open circuit potential or reversible potential (£rev) is given by the overall Nernst equation.

The electrolyte in the examples considered here (YSZ) is a purely ionic conductor and thus there is no faradaic current generation within the electrolyte membrane. Therefore, within the electrolyte layer that separates the anode and cathode, the governing equation for the ionic potential is Eq. 23.

( d Vi dy\ldy

Solution algorithm.— A flow chart showing the sequence of steps is shown in Figure 1. In the case of Model-I, Eqs. 8, 9 and 15 form a system of algebraic equations that can be solved to get the current density i, anode activation overpotential na and cathode activation overpotential nc by fixing the cell voltage. If the current density is fixed, the solution of the system of equations results in the cell voltage in addition to the activation overpotentials. The solution of these equations requires the knowledge of exchange current densities, which can be calculated once the species composition at the TPB is known (Eqs. 10 and 11). The species composition is evaluated by solving Eqs. 1 to 4. The above differential equations are numerically integrated using the ODE solver CVode.35 The solution is considered as converged when the residuals fall below 10-6.

Figure 1. Flowchart showing the sequence of calculations for the models used.

In the case of Model-II, Eqs. 17 to 23 are solved together with Eqs. 1 to 4 to get the distribution of current and potential within the electrodes. In the model results presented in this work, the steady state form of Eqs. 17 to 20 is used as the timescale of the double layer charging and discharging is of the order of milliseconds.32 The potential boundary condition required to solve the composite electrode problem is shown schematically in Figure 2.

Sulfur poisoning.— The exchange current density parameter i* appearing in Eq. 10 is a lumped parameter that includes the effect of available TPB, temperature, and rate constants for the electrochemical charge transfer reactions. Therefore, modeling sulfur poisoning requires a functional relation that relates the intrinsic value of i* with its value when the surface is poisoned. If l0 is the available TPB length before poisoning and 9S is the fractional coverage of sulfur at time t, we assume that the TPB length during poisoning is given by Eq. 24

l(t) = lo (1 - 9s). [24]

Although the expression in Eq. 24 is very similar to the one reported by Ryan et al.36 for the reduction in TPB length using a 'damage factor'

f (T, pH2S), unlike Ryan et al., we do not fit the data to estimate the drop in TPB length and instead use the sulfur coverage (9S) predicted by the detailed micro-kinetic model proposed by Appari et al.37 The model proposed by Appari et al., is validated against the experimental

Figure 2. Schematic representation of potential boundary conditions required for composite electrode problem.

Table I. Detailed kinetic model for reforming of biogas.37

R No Reaction A(cm,mol,s) ß Ea a

R1 H2 + (Ni) + (Ni) ^ H(Ni) + H(Ni) 0.01b 0 0

R2 O2 +(Ni) +(Ni) ^ O(Ni) +O(Ni) 0.01b 0 0

R3 H2S +(Ni) ^ H2 S(Ni) 0.6b 0 0

R4 H2O +(Ni) ^ H2O(Ni) 0.1b 0 0

R5 SO2 +(Ni) ^ SO2 (Ni) 0.02b 0 0

R6 H(Ni) +H(Ni) ^ (Ni) +(Ni) +H2 2.676x1019 0 81.40

R7 O(Ni) +O(Ni) ^ (Ni) +(Ni) +O2 4.143x 1023 0 474.93

R8 H2O(Ni) ^ (Ni) +H2O 3.823x 1012 0 60.78

R9 H2S(Ni) ^ H2S +(Ni) 1.108x 1010 -0.8 69.47

R10 SO2(Ni) ^ SO2 +(Ni) 2.709 x1009 0 102.50

R11 O(Ni) +H(Ni) ^ OH(Ni) +(Ni) 5x 1022 0 97.90

R12 OH(Ni) +(Ni) ^ O(Ni) +H(Ni) 1.793x 1021 0 36.14

R13 OH(Ni) +H(Ni) ^ H2O(Ni) +(Ni) 3 x 1020 0 42.70

R14 H2O(Ni) +(Ni) ^ OH(Ni) +H(Ni) 2.251 x1021 0 91.79

R15 OH(Ni) +OH(Ni) ^ O(Ni) +H2O(Ni) 3 x 1021 0 100.00

R16 O(Ni) +H2O(Ni) ^ OH(Ni) +OH(Ni) 6.276x 1023 0 210.85

R17 H2S(Ni) +(Ni) ^ SH(Ni) +H(Ni) 5.5 x104 1.2 29.31

R18 SH(Ni) +H(Ni) ^ H2S(Ni) +(Ni) 1.291 x1013 0 106.19

R19 SH(Ni) +(Ni) ^ S(Ni) +H(Ni) 7.9x 1011 0 25.79

R20 S(Ni) +H(Ni) ^ SH(Ni) +(Ni) 6.375x 1015 0 142.94

R21 SH(Ni) +OH(Ni) ^ H2S(Ni) +O(Ni) 1.053x 1013 0 29.72

R22 H2S(Ni) +O(Ni) ^ SH(Ni) +OH(Ni) 8x 1011 -0.5 27.84

R23 S(Ni) +O(Ni) ^ SO(Ni) +(Ni) 1 x 1018 1 296.82

R24 SO(Ni) +(Ni) ^ S(Ni) +O(Ni) 1.775x 1012 0 0.0

R25 SH(Ni) +O(Ni) ^ SO(Ni) +H(Ni) 1 x 1014 -1 206.05

R26 SO(Ni) +H(Ni) ^ SH(Ni) +O(Ni) 2.115 x105 0 0

R27 S(Ni) +OH(Ni) ^ SO(Ni) +H(Ni) 1 x 1021 1 229.02

R28 SO(Ni) +H(Ni) ^ S(Ni) +OH(Ni) 3.352x 1023 -2.0 0.00

R29 SO2(Ni) +(Ni) ^ SO(Ni) +O(Ni) 1 x 1018 -0.5 106.31

R30 SO(Ni) +O(Ni) ^ SO2 (Ni) +(Ni) 9.029x 1009 1.5 0.00

R31 S(Ni) +H2O(Ni) ^ SH(Ni) +OH(Ni) 1 x 1010 0 143.37

R32 SH(Ni) +OH(Ni) ^ S(Ni) +H2O(Ni) 1.652x 105 0 0.0

aArrhenius parameters for the rate constants written in the form: k = AT exp(-Ea/RT) The units of A are given in terms of moles, centimeters, and seconds. Ea is in kJ/mol

bSticking coefficient. Total available surface site density is T=2.66x 10-9 mol/cm2

measurements by scaling the specific surface area with (1 - 9S). To maintain consistency we follow the functional relation proposed by Appari et al., nevertheless, the general conclusions drawn in this study are independent of the functional dependency of TPB area on sulfur coverage. Assuming the above linear relation Eq. 12 can now be written as

i*a = (1 - ÖS)

ia eXP -

The term in the square brackets in Eq. 25 is the intrinsic value of the exchange current density. In order to use the above expression one needs to evaluate i*' and EH2. This is done by calibrating the exchange current density parameters to reproduce experimentally observed DC polarization data (i - V data) in the absence of H2S.

To simulate a gas mixture consisting of H2, H2O, and H2S, it is sufficient to consider the subset of reactions proposed by Appari et al.,37 given in Table I. This reaction set consists of 16 reversible reactions among 5 gas-phase species (H2, H2O, O2, H2S, SO2), 10 surface adsorbed species (H, O, S, OH, H2O, H2S, SO2, SH, SO, SH) on Ni, and the vacant sites on the catalytic surface of Ni. The total surface site density used is 2.66 x 10-9 mol cm-2. The micro-kinetic model enables the calculation of sulfur coverage throughout the thickness of anode, which alters the anode exchange current density according to Eq. 25. Since the detailed kinetic model taken from Appari et al., is validated only for H2S concentrations in the range of 20-100 ppm, the SOFC anode sulfur poisoning model validations and discussions in this work are limited to the high H2S concentration cases reported in the literature. The kinetic model proposed by Appari et al., is devel-

oped using experiments performed on Ni/y-Al2O3. Nevertheless, the mechanism can also be used for the Ni/YSZ system since the catalytic activity of zirconia is negligible compared to that of Ni as shown by both experimental as well as modeling studies.38,39 Hecht et al.38 used a detailed kinetic model developed based on steam reforming studies of CH4 on Ni/y-Al2O3 in monolith reactors to model experiments performed on Ni-YSZ and demonstrated excellent agreement between the model predictions and experimental observations. Shishkin and Ziegler39 used density functional theory to probe the catalytic activity of YSZ toward methane reforming and showed that methane and hydrogen would not adsorb and dissociate on YSZ unless the YSZ had excess oxygen.

Generally, one would expect the surface reactions to be dependent on the coverage of sulfur. As mentioned earlier, the experimental and empirical investigations of Alstrup et al.16 have shown that the enthalpy of adsorption for sulfur on nickel decreases with increasing sulfur coverage and they were able to fit a Temkin isotherm to experimental adsorption data. Hansen17 then used the above Temkin isotherm to correlate the predicted sulfur coverage and the performance drop seen in SOFC sulfur poisoning experiments. More recent atomistic modeling work by Monder and Karan41,42 has also proposed a strong coverage dependence for the adsorption of sulfur on Ni. However, in the interest of simplicity, the sulfur species reactions in the reaction mechanism presented here are assumed to be independent of the sulfur coverage. This is a valid assumption for the purposes of this work because: 1) this work focuses on explaining the effect of poisoning on SOFC anode side losses and uses experimental results from a particular range of H2S concentration2 and temperature; 2) we

Table II. SOFC button cell membrane electrode assembly (MEA) parameters for Model-I.

Parameters Anode Cathode Reference

Thickness , l (Mm) 50 50 Liu40

Porosity, e 45% 45% Typical values

Tortuosity, T 3.8 3.8 Typical values

Pore radius, rp (Mm) 0.5 0.5 Typical values

Particle diameter, d p 1.0 1.0 Typical values

Exchange current 3.648X102 2.703 X102 Evaluated

density, i*' (A cm-2)

Exchange current 46.643 59.781 Evaluated

density, E (kJ mol-1

Symmetry factor, ßa 0.5 0.3 Evaluated

use a reaction mechanism that has been validated in the range of concentration and temperature used in Zha et al.2 Thus, the kinetics and thermodynamics of the surface reactions (independent of 9s) can be assumed to be valid within the range of operating conditions studied in this work.

Results and Discussion

We start by simulating the experiments before sulfur was introduced into the fuel stream in Zha et al.2 Both electrochemical models (Model-I and Model-II) are used for the simulations. The parameters used for Model-I are given in Table II and a comparison of the model predicted and experimentally observed DC polarization data is shown in Figure 3. The pre-exponential factors for the exchange current densities of the anode and cathode are tuned to fit our model to the data. The model predictions and experimental data are in very good agreement with each other. The series resistance in the model (resistance of electrolyte plus the contact resistance) is tuned to reproduce the Ohmic resistance reported by Zha et al.2 The physical properties of the electrodes are common between Model-I and Model-II, but to account for volumetric current generation, volumetric electrochemical parameters and solid-phase volume fractions are required to solve Model-II. These parameters for Model-II are given in Table III. Again good agreement is observed between the model predictions and experimental observation.

In order to emphasize the predictive capability of the reaction mechanism used, a comparison between methane reforming exper-

0.2 0.4 0.6

Current density / A cm

0.8 1 -2

Figure 3. Comparison between experimentally observed and Model-I predicted DC polarization data for 50% H2, 1.5% H2O and 48.5% N2. The model parameters used are given in Table II.

Table III. SOFC button cell MEA parameters for Model-II.

Parameters Anode Cathode Reference

Exchange current density, i *' (A cm -3) 17x 102 15x102 Evaluated Symmetry factor, ßa 0.5 0.62 Evaluated

Ionic phase volume fraction 22% 1% Zha et al2

0.2 0.1 0

<>«• ..................

,„.-;■ - " • * U HI •

... • • • • • • ■ —----

6 8 Time / h

Figure 4. Comparison of model predicted and experimentally observed mole fractions at the reactor exit as a function of time in a fixed bed reactor operating at 800°C with 12.5% CH4, 8.4% CO2, 25.2% H2O, 53.9% N2 and 50 ppm H2S. Reproduced from.43

iments performed by Appari43 in the presence of H2S, and model predictions using the mechanism borrowed for use in this work37 is shown in Figure 4. The figure shows the model predicted and experimentally observed gas-phase composition on a dry basis at the exit of a fixed bed reactor (Ni/y-Al2O3) operating at 800°C with 12.5% CH4, 8.4% CO2, 25.2% H2O, 53.9% N2 and 50 ppm H2S. Although the reaction mechanism does not consider the sulfur coverage dependence of individual elementary reactions, it is able to predict the exit gas composition reasonably well over the entire time span of the experiment.

Current distribution in Model-II.— Model-II allows charge transfer throughout the thickness of the anode, and Figure 5a presents the distribution of the ionic and electronic current flux across the thickness of the anode at an operating voltage Vcell = 0.65 V. As expected, the ionic current decreases from the anode-electrolyte interface to the current collector (left to right), while the electronic current increases in the same direction. For any steady-state operating point, the sum of

a) before poisoning

ionic ---after poisoning


V b) ionic v Z \ \ /

\ electronic

// -- // // //

10 20 30 40

distance from electrolyte, microns

Figure 5. Distribution of ionic and electronic flux before and after poisoning across the thickness of the anode, V^n = 0.65 V.

Exp * Model-I ■ Model-II

& 0.3 -

0.2 0.25 0.3 0.35

Initial current density / A cm

0.4 0.45 -2

Figure 6. Comparison between experimentally observed and model predicted drop in current on exposure to 50 ppm H2S. 50% H2, 1.5% H2O and 48.5% N2, T = 800oC.

the ionic and electronic current at any point in the anode is identical. It is interesting to note the shape of the current distribution which shows that the electronic current rises more sharply before poisoning and thus the current profile through the electrode is somewhat flatter after poisoning. This is clearer in Figure 5b where the scaled current density profiles are plotted across the anode. The above observation is in agreement with the analysis presented by Gazzarri and Kesler15 who also showed that the current generation profile becomes flatter or more uniform after poisoning.

Effect of current density on degree of poisoning.— A comparison between the model predicted and experimentally observed drop in cell current at 800 oC on introducing 50 ppm H2S into the fuel (50% H2, 1.5% H2O and 48.5% N2) is shown in Figure 6. Excellent agreement is observed between the predictions of both the models with experimental measurements reported by Zha et al.2 The predicted surface coverages after poisoning are shown in Figure 7. The detailed kinetic model used here accounts for the partial pressures of all gas-phase species at the electrochemically active interfaces to calculate the coverages of the surface species resulting from the adsorption/dissociation of gas-phase species on the Ni surface at the TPB. At higher current densities, the partial pressure of H2 at the anode-electrolyte interface will be lower than that at the anode free surface while the H2O partial pressure will be higher. However, as seen in Figure 8, these variations

<u o> E


10 15 20 25 30 35 Position /(microns)

40 45 50

Figure 7. Major surface adsorbed species in the anode on exposure to 50 ppm H2S, 50% H2, 1.5% H2O and 48.5% N2 at 800°C, Vcell = 0.6 V.

H2, electrolyte supported

H2J anode supported

coverage, anode supported

- 0.76

h 0.74 D> G

- 0.72

coverage, electrolyte supported

0.1 H H20, anode supported <-

H30, electrolyte supported

0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 Normalized anode thickness, x/l.

Figure 8. ph2 , Ph2o, and 0 s at steady state across the thickness of the anode for Veen = 0.3V. 50% H2, 1.5% H2O and 48.5% N2, T = 800oC. The anode free surface is the left-hand edge and the anode-electrolyte interface is right-hand edge of the plot.

are insignificant for a thin anode in an electrolyte supported cell such as the one used by Zha et al.2 The horizontal coordinate gives the dimensionless length (x/ la) along the thickness of the anode. Interestingly the model predicts that the surface coverage of sulfur at the three phase interface is not influenced by cell voltage or current density either and all cases resulted in a final sulfur coverage of ~ 72% for the above conditions.

To test whether the sulfur coverage in the active TPB region is affected at the higher current densities seen in anode-supported cells, we ran the model for a cell with the following dimensions: la = 500 |xm, lc = 50 |xm, and le = 25 |xm, without changing the other parameters. The results for a highly polarized cell (Vcell = 0.3 V) in Figure 8 show that the S coverage is ~ 0.7% higher at the anode-electrolyte interface for this anode-supported cell. It should be noted that this effect only emerges because of the significantly larger change in the local gas phase composition, which in turn is due to the thicker anode (10x that for the cell used by Zha et al.) as well as the much higher current density (1.7x that for the cell used by Zha et al.).

Figure 9 shows the sensitivity of the calculated cell current to the various kinetic rate constants in the reaction mechanism. The normalized sensitivity (NS) is calculated for a given voltage as

(i - io)/ io

max{|(i - io)/io|}'

where io is the predicted current without any change in the kinetic rate constants, and i is the predicted current with ±10% change in

Figure 9. Normalized sensitivity of various reactions (see Table I) on cell power output. Red bars represent 10% increase in pre-exponential factors while blue bars represent a 10% decrease. 50 ppm H2S, 50% H2, 1.5% H2O and 48.5% N2 at 800°C, Vcell = 0.7 V.


0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1

Current density/ A cm"'

Galvanostatic •

Potentiostatic ■

J* __—-e---!

Current density/ A cm

Figure 10. Performance drop of the cell due to S poisoning during potentio-static and galvanostatic operation measured as relative drop in power. 50 ppm H2S, 50% H2, 1.5% H2O and 48.5% N2, T = 800°C.

Figure 11. Variation of anode activation overpotential na before (lower curve) and after poisoning (upper curve) for potentiostatic and galvanostatic operation. 50 ppm H2S, 50% H2, 1.5% H2O and 48.5% N2, T = 800°C.

the specified rate constant. The denominator in Eq. 26 is the absolute value of the largest change observed by a 10% perturbation in kinetic constants. In this work, a change in the H2 adsorption reaction rate constant (sticking factor) gives the largest change in current. A positive change of 10% in the sticking coefficient for H2 lead to an increase of 1.7% in current density (with S poisoning). The sensitivity calculations given here are performed at V'cell = 0.7 V. For any given cell voltage, the resulting current depends on the sulfur coverage. A lower 9s gives a higher current and vice versa. Therefore, it can be inferred that a positive change in the sticking coefficient of H2 (R1) leads to lower sulfur coverage. Similarly a positive change in the sticking coefficient of H2S (R3) leads to a decrease in current output from the cell. A positive change in the rate constant for OH dispro-portionation (R15) also leads to a reduction in current. The only other reaction rate constants (R18-R20) that had an effect on the current were associated with H2S dissociation and the current decreases as H2S and SH dissociation increases on the Ni surface. Only 8 out of 32 irreversible reactions were found to have significant influence on current output from the cell. It is worth noting that H2O adsorption (R4) does not have any significant influence on the cell current.

The model predictions above clearly show that the lower relative performance drop at higher current density (or lower cell voltage) is not due to a change in the sulfur coverage. The computed relative power drop (Eq. 27) of the cell for potentiostatic and galvanostatic operation upon introducing 50 ppm H2S into the fuel is shown in Figure 10. For galvanostatic operation, the percentage drop in performance is worse for increasing polarization or current. For potentiostatic operation however, the performance drop becomes worse with decreasing polarization or current. These contradictory trends are in very good agreement with the experimental literature. An explanation for these trends is provided by Cheng et al., using equivalent circuit analysis.3 As pointed out by Cheng at al., characterizing the degree of sulfur poisoning in terms of relative drop in cell power output is misleading. The apparent contradiction in trends disappears if the degree of poisoning is expressed in terms of the relative change in activation overpotential (or the polarization resistance) of the anode: see Eqs. 28-29.

A Pe, =

Arirel =

A Rp, rel =

n — n

Rp — R0p

Here, n0, R0p are the anode overpotential, and anode polarization resistance of the cell before poisoning while n, Rp represent the same measures after poisoning. The activation overpotential before and after poisoning for potentiostatic and galvanostatic operation is

shown in Figure 11. The increase in overpotential after poisoning at any given current density is due to the decrease in the exchange current density according to Eq. 25. Figure 11 clearly shows that the activation overpotential of the anode for both potentiostatic and galvanostatic operation falls on the same set of curves before and after poisoning. As expected, the overpotential (anode voltage loss) increases with current and is higher after poisoning.

On a closer examination of Figure 11, it is obvious that Ar|rel (see Eq. 28) changes as a function of current density because the polarization curves are non-linear. Clearly, the ratio (r|a-^10)/^ is greater at low current density than the same quantity calculated at a higher current density. Thus, the relative increase in anode overpotential due to sulfur poisoning is higher at lower current densities. It should be noted that this effect is solely due to the non-linearity of the electrochemical rate equation. This effect can also be explained in terms of polarization resistance of the anode before and after poisoning as shown in Figure 12. For the given conditions, at the open circuit voltage, the anode polarization resistance after poisoning is 3.4 times greater than the value before poisoning. However, for a current density i = 0 . 2 A cm-2, the Rp for the poisoned condition is only 1.7 times greater than the value before poisoning. In other words, the increase in Rp upon poisoning is twice as large at open circuit conditions than at i = 0 . 2A cm-2. This effect has also been explained from an experimental point of view by Cheng et al.,1 who suggested that the observed decrease in relative polarization resistance of the anode with increase in current density need not be due to reduced degree of sulfur poisoning and can be seen as a consequence of the Tafel equation.

A discussion of sulfur oxidation reactions.— Most experimental studies have explained the above reduction in anode Rp at higher current densities by hypothesizing electrochemical oxidation of the adsorbed sulfur near theTPB. , , These studies, among others, have suggested that the increased flux of O2- ions at the TPB due to the

4.5 4 3.5 3 2.5 2 1.5 1 0.5 0

Before poisoning After poisoning

Current density/ A cm

Figure 12. Anode polarization resistance before and after poisoning as a function of cell current density. 50 ppm H2S, 50% H2, 1.5% H2O and 48.5% N2, T = 800° C.

C ument density / A cm"2

Figure 13. Relative change in anode polarization resistance due to sulfur poisoning at different operating temperatures. 50 ppm H2S, 50% H2, 1.5% H2O and 48.5% N2.

higher current could clean the Ni surface near the TPB by oxidizing the sulfur at the TPB. However, the model simulations here show that this hypothesis is unnecessary and that the inherent non-linearity of the hydrogen electrochemistry is enough to explain the observed experimental trends at higher current densities. Our model does NOT include electrochemical pathways for the transfer of oxygen from the YSZ to the Ni or for the oxidation of sulfur on YSZ: reactions Eqs. 30-35. The phases in the bracket for the species in these reaction equations represent the surface on which each species is adsorbed. These reactions would augment the reaction mechanism given in Table I and lead to the electrochemical oxidation of sulfur.

O2-(oxide) + V(Ni) ^ V(oxide) + O(Ni) + 2e-(Ni) [30]

O2-(oxide) + S(Ni) ^ V(oxide) + SO(Ni) + 2e-(Ni) [31]

SO(Ni) + O2-(oxide) ^ SO2(Ni) + V(oxide) + 2e-(Ni) [32]

O2- (oxide) + S(Ni) ^ SO(oxide) + V(Ni) + 2e- (Ni) [33]

SO(oxide) + O2-(oxide) ^ SO2(oxide) + V(oxide) + 2e-(Ni)

SO2(oxide) ^ SO2 + V(oxide) [35]

Although the model used in this work allows for the catalytic oxidation of sulfur by adsorbed oxygen, the only source allowed for the adsorbed oxygen on the Ni is the dissociative adsorption of water vapor. In summary, although the electrochemical model used in this work does not account for electrochemical charge transfer between adsorbed sulfur and O2-, it is still able to predict the observed experimental trend of a reduction in performance losses due to poisoning at higher current densities with sufficient accuracy.

Effect ofT on sulfur poisoning.— Our model also reproduces the well-known trend of a lower extent of poisoning with increasing temperature as seen in Figure 13. This is easily explained by the decrease in sulfur coverage with increasing temperature and the resulting effect on the exchange current density through Eq. 25. Another reason for the lower extent of poisoning is that the exchange current density increases exponentially with temperature (Eq. 12) and thus the anode polarization resistance is lower at higher temperatures. The lower initial Rp leads to a lower relative drop in power density on poisoning. This effect is explained in more detail in the next section.

Sulfur poisoning in Ni cermet anodes with oxides other than YSZ.— A number of studies of SOFC anode sulfur poisoning have shown that Ni cermet anodes with oxides other than YSZ, especially those with a higher oxygen ion conductivity (e.g., doped ceria), exhibit a lower degree of poisoning than Ni-YSZ anodes.11,5,44,45 The principal

Table IV. Effect of increasing exchange current density on relative

impact of sulfur poisoning. T =800 oC, i = 0.2 A cm-2, all other

parameters from II.

i*a' (A/cm2 ) Vcell (V) Vell (V) ARp,rel A Prel

45.86 0.755 0.683 1.02 0.095

458.6 0.814 0.789 2.33 0.031

hypotheses to explain this have been similar to the ones proposed to explain the lower degree of poisoning at higher current density: a) possible electrochemical spillover of O2- ions from the oxide surface to the Ni surface at the TPB (reaction Eq. 30); b) the electrochemical reaction of adsorbed sulfur species with O2- ions at the TPB (reactions Eq. 31-35).

We do not consider reactions (Eqs. 30-35) in our model and are nevertheless able to show that the experimental observations can be explained without the electrochemical oxidation of sulfur. Table IV presents simulation results for two cases where the only difference is that the anode exchange current density is 10 times higher in the second case compared to the base case (the anode used by Zha et al.2). Ni-SSZ and Ni-ceria (Ni-GDC or Ni-SDC) anodes have been demonstrated to have much higher activities than Ni-YSZ anodes and thus can be represented by the more active anode in Table IV. Whether a particular Ni-oxide anode is exactly 10 times more active than a Ni-YSZ is not important; what is important is the consequence of having an anode that has a lower anode Rp and thus a lower anode overpotential than Ni-YSZ. We would like to emphasize that the electrochemical reaction mechanism used in this work (originally proposed by Zhu et al.20) assumes that the only charge transfer steps are elemental hydrogen transfer reactions from the Ni to the oxide at the TPB.

In Table IV, V^ is the cell voltage at the given conditions (T=800 oC, i = 0.2 A cm-2) before poisoning while Vcell is the cell voltage after poisoning at the same temperature and current density. Vc0ell is understandably higher for the cell with the more active anode because of the lower anode polarization resistance. Interestingly, although the relative increase in anode polarization resistance is more than twice as high for the more active anode, the relative performance drop is three times lower. This is because the initial polarization resistance ( R0p ) is substantially lower for the more active anode and therefore the increase in Rp due to poisoning has a much lower impact on the cell voltage and power density.


Sulfur poisoning and the resulting performance degradation for galvanostatic as well as potentiostatic operation of solid oxide fuel cells is studied using a 1-D SOFC model coupled to a detailed micro-kinetic model and the results are compared with experimental data as well as trends reported in the literature. The micro-kinetic model accounts for the adsorption, dissociation and desorption of all gas phase species as well as surface reactions among adsorbed species and is thus used to calculate the sulfur coverage at the TPB. The SOFC model assumes that the sulfur absorbed in the TPB region results in a reduction in TPB length which can be expressed in terms of reduction in the anode exchange current density. The model predicted performance drop is in very good agreement with experimental data.

This model is capable of predicting seemingly contradictory trends when comparing the relative performance drop of an SOFC as a function of cell polarization in galvanostatic or potentiostatic operation. The much discussed contradictions vanish if the performance drop is expressed in terms of the relative change in anode overpotential or polarization resistance. The model simulations also show that, under typical SOFC operating conditions, the sulfur coverage at the anode TPB region does not change appreciably with current density. This study then presents an analysis demonstrating that the observed behavior of a lower degree of poisoning at high current density can be

completely explained by the non-linearity of the electrochemical rate equations without recourse to the usual hypothesis of the electrochemical oxidation of adsorbed sulfur in the TPB region.

Finally the model simulations are used to show that the experimentally observed behavior of a significantly lower degree of poisoning for Ni cermet anodes such as Ni-GDC and Ni-ScSZ compared to Ni-YSZ can also be explained by the much lower anode polarization resistance (both before and after poisoning) of such anodes. It should be noted that the model parameters are not adjusted to reproduce the experimental observations during poisoning. Rather our model predictions are a result of the multi-scale modeling of coupled interactions of various physical, chemical, and electrochemical processes in an operating SOFC.


We acknowledge the funding received from the Department of Science and Technology, Government of India under the project SR/RC-UK/Fuel-Cell-03/2011/IITH (G).


1. Z. Cheng, J. H. Wang, Y. Choi, L. Yang, M. Lin, and M. Liu, Energy Environ. Sci, 4, 4380 (2011).

2. S. Zha, Z. Cheng, and M. Liu, J. Electrochem. Soc., 154, B201 (2007).

3. Z. Cheng, S. Zha, and M. Liu, J. Power Sources, 172, 688 (2007).

4. Y. Matsuzaki and I. Yasuda, Solid State Ionics, 132, 261 (2000).

5. A. Hagen, J. F. Rasmussen, and K. Thyden, J. Power Sources, 196, 7271 (2011).

6. A. Kromp, S. Dierickx, A. Leonide, A. Weber, and E. Ivers-Tiffee, J. Electrochem. Soc., 159, B597 (2012).

7. A. Weber, S. Dierickx, A. Kromp, and E. Ivers-Tiffee, Fuel Cells, 13, 487-493, 10th European SOFC Forum, Lucerne, Switzerland, JUN 26-29, 2012.

8. A. Hagen, J. Electrochem. Soc., 160, F111 (2013).

9. M. Ashrafi, C. Pfeifer, T. Pro, and H. Hofbauer, Energy and Fuels, 22, 4190 (2008).

10. S. Appari, V. M. Janardhanan, R. Bauri, and S. Jayanti, Int. J. Hydrogen Energy, 39, 297 (2014).

11. K. Sasaki, K. Susuki, A. Iyoshi, M. Uchimura, andN. Imamura, J. Electrochem. Soc., 153, A2023 (2006).

12. J. F. B. Rasmussen and A. Hagen, J. Power Sources, 191, 534 (2009).

13. J. B. Hansen and J. Rostrup-Nielsen, In Handbook of Fuel Cells: Fundamentals, Technology, Applications; Wolf Vielstich, H. A. G. , and Arnold Lamm, Ed.; Wiley, 2009; Vol. 6; Chapter 65 Sulfur poisoning on Ni catalyst and anodes.

14. E. Brightman, D. G. Ivey, D. J. L. Brett, andN. P. Brandon, Journalof Power Sources, 196,7182 (2011).

15. J. Gazzarri and O. Kesler, J. Power Sources, 167, 100 (2007).

16. I. Alstrup, J. R. Rostrup-Nielsen, and S. Roen, Applied Catalysis, 1, 303 (1981).

17. J. B. Hansen, Electrochem. Solid-State Lett., 11, B178 (2008).

18. D. S. Monder and K. Karan, ECS Transactions, 35, 977 (2011).

19. B. V.R. S.N. Prasad and V.M. Janardhanan, J. Electrochem. Soc., 161, 208 (2014).

20. H. Zhu, R. J. Kee, V. M. Janardhanan, O. Deutschmann, and D. G. Goodwin, J. Electrochem. Soc., 152, A2427 (2005).

21. R. J. Kee, M. E. Coltrin, and P. Glarborg, Chemically reacting flow: Theory and Practice; John Wiley & Sons, 2003.

22. S. Appari, V. M. Janardhanan, S. Jayanti, L. Maier, S. Tischer, and O. Deutschmann, Chem. Eng. Sci., 66, 5184 (2011).

23. V. M. Janardhanan and O. Deutschmann, Z. Phys. Chem, 221, 443 (2007).

24. V. M. Janardhanan, A detailed approach to model transport, heterogeneous chemistry, and electrochemistry in solid-oxide fuel cells. Ph.D. thesis, KIT, 2007.

25. R. J. Kee, F. M. Rupley, J. A. Miller, M. Coltrin, J. Grcar, E. Meeks, and H. Moffat, SURFACE CHEMKIN-III; 1999.

26. H. Zhu and R. J. Kee, J. Power Sources, 117, 61 (2003).

27. P. Costamagna, A. Selimovic, M. Del, and G. Agnew, Chem. Eng. J., 102, 61 (2004).

28. P. Aguiar, C. S. Adjiman, and N. P. Brandon, J. Power Sources, 138, 120 (2004).

29. P. Iora, P. Aguiar, C. S. Adjiman, and N. P. Brandon, Chem. Eng. Sci., 60, 2963 (2005).

30. W. G. Bessler and S. Gewies, J. Electrochem. Soc., 154, B548 (2007).

31. M. Ni, D. Y. C. Leung, and M. K. H. Leung, J. Power Sources, 185, 233 (2008).

32. K. Nanaeda, F. Mueller, J. Brouwer, and S. Samuelsen, J. Power Sources, 195, 3176 (2010).

33. Y. Shi, N. Cai, C. Li, C. Bao, E. Croiset, J. Qian, Q. Hu, and S. Wang, J. Electrochem. Soc., 155, B270 (2008).

34. H. Zhu and R. J. Kee, J. Electrochem. Soc., 155, B715 (2008).

35. A. Hindmarsh, P. Brown, K. Grant, S. Lee, R. Serban, D. Shumaker, and C. Woodward, ACM Trans. Math. Softw., 31, 363 (2005).

36. E. M. Ryan, W. Xu, X. Sun, and M. A. Khaleel, J. Power Sources, 210, 233 (2012).

37. S. Appari, V. M. Janardhanan, R. Bauri, S. Jayanti, and O. Deutschmann, Appl. Catal. A, 471,118 (2014).

38. E. S. Hecht, G. K. Gupta, H. Zhu, A. M. Dean, R. J. Kee, L. Maier, and O. Deutschmann, Appl. Catal. A Gen., 295, 40 (2005).

39. M. Shishkin and T. Ziegler, J. Phys. Chem. C, 112, 19662 (2008).

40. M. L. Liu, Personal communication (2013).

41. D. S. Monder and K. Karan, J. Phys. Chem. C, 114, 22597 (2010).

42. D. S. Monder and K. Karan, ECS Trans, 57, 2449 (2013).

43. S. Appari, Experimental and Theoretical Investigation of Catalyst Poisoning and Regeneration During Biogas Steam Reforming on Nickel. Ph.D. thesis, Indian Institute of Technology Hyderabad, 2014.

44. L. Zhang, S. P. Jiang, H. Q. He, X. Chen, J. Ma, and X. C. Song, Int. J. Hydrogen Energy, 35, 12359 (2010).

45. S. K. Schubert, M. Kusnezoff, A. Michaelis, and S. I. Bredikhin, J. Power Sources, 217, 364 (2012).