Bistability and Oscillations during Electrooxidation of H2−CO Mixtures

Aug 18, 2007 - Mechanistic aspects of oscillations during CO electrooxidation on Pt in the presence of anions: Experiments and simulations. S. Malkhan...
2 downloads 9 Views 344KB Size
Download PDF
J. Phys. Chem. C 2007, 111, 13481-13489

13481

Bistability and Oscillations during Electrooxidation of H2-CO Mixtures on Pt: Modeling and Bifurcation Analysis Jan Siegmeier, Nilu1 fer Baba, and Katharina Krischer* Physik Department E19, TU Mu¨nchen, James-Franck Strasse 1, 85748 Garching, Germany ReceiVed: April 5, 2007; In Final Form: July 6, 2007

A model describing the electrooxidation of H2-CO mixtures on rotating Pt electrodes is presented together with a bifurcation analysis for the main experimental parameters. We discuss dynamic instabilities partly caused by a chemical autocatalysis stemming from the CO subsystem and partly from an electrical autocatalysis with the electrode potential as the essential positive feedback variable. Hence, the system contains two different pairs of activator-inhibitor feedback loops. Depending on the parameters, one of the two interaction schemes is dominant, and we find either a so-called S-NDR or an HN-NDR system. However, there are also parameter regimes where the interplay of both loops gives rise to oscillations, greatly enlarging the dynamically unstable parameter range.

I. Introduction Carbon monoxide is a well-known poison for electrocatalytic reactions on Pt electrodes. One of the most prominent examples are low-temperature fuel cells, which exhibit a dramatic efficiency decrease when the hydrogen feed stream is contaminated with traces of CO. However, small CO impurities are inevitable in practical applications if hydrogen is generated from hydrocarbons. This prompted intense research activities on the electrooxidation mechanism of H2-CO mixtures on Pt and PtRu electrodes, mostly focusing on electrode kinetics and interfacial structures (see ref 1 and references therein). Only very few studies are concerned with dynamical instabilities, such as sustained voltage or current oscillations,2-9 although first reports of oscillatory responses under constant current conditions date back several decades.2 Yamazaki and Kodera4 made a first attempt to model the oscillations. Their model introduced attractive interactions between adsorbed CO molecules to induce the dynamic instability, an assumption which clearly contradicts surface electrochemical results.10 Zhang and Datta5-7 performed experiments in proton exchange membrane fuel cells (PEMFC) with PtRu as anode catalyst. Compared with the steady state, the PEMFC yields a higher output when operated under autonomous oscillatory conditions. Their model is based on surface steps and material balance equations for the concentrations of CO and H2 in the anode chamber and reproduces their experimental results. Kiss et al.9 demonstrate that two platinum electrodes may exhibit different reaction rates as a result of a symmetry breaking when coupled through a common external resistor. Here, we focus on the theoretical description of the electrooxidation of H2-CO mixtures at rotating Pt electrodes, with an emphasis on elucidating the essential feedback loops that generate the dynamic instabilities. Furthermore, we present a thorough bifurcation analysis. Our theoretical description is based on a lumped model for CO electrooxidation on Pt,11 extended by an appropriate description of the hydrogen oxidation current. * Corresponding author. E-mail: [email protected]. URL: http:// www.physik.tu-muenchen.de/lehrstuehle/E19/AG%20Krischer.html.

From a nonlinear dynamics point of view, this system is intriguing, as it consists of two subsystems that may undergo oscillatory instabilities: CO electrooxidation on rotating electrodes is known to possess an S-shaped current potential curve11,12 and thus belongs to the so-called S-NDR (negative differential resistance) class. S-NDR systems exhibit a hysteresis between a high and a low current state for low cell resistance, so they are bistable. The destabilizing autocatalytic loop is of a chemical nature, the electrode potential taking on the role of an inhibitor variable. The characteristic time changes of the electrode potential is essentially given by the time needed to charge or decharge the double layer. It is much shorter than characteristic times in which concentrations of chemical species change because of reaction and transport mechanisms. Thus, S-NDR systems can be characterized as activator-inhibitor type systems with slow activator and fast inhibitor dynamics.13 In general, such systems do not oscillate but tend to exhibit spatial instabilities.14-16 On the other hand, it was suggested that the H2-CO system behaves similar to the H2, Cu2+, Cl-|Pt system, which is a typical representative of the so-called HN-NDR class.17 Here, the adsorption of Cl- at high potentials causes a decrease of the hydrogen current with increasing electrode potential; that is, the stationary polarization curve possesses an N-shaped negative differential resistance (N-NDR). In the presence of Cu2+ ions, the number of free surface sites gained by an increased potential due to Cu stripping overcompensates the number of free surface sites lost because of Cl- adsorption. Thus, the N-NDR is hidden in the potential range in which Cu desorbs from the surface. In the H2-CO system, OH and CO seem to play the roles of Cl- and Cu2+ in the H2, Cu2+,Clsystem, respectively. The electrode potential is the activatory variable in HN-NDR systems, while the negative feedback originates from chemical reaction steps. Thus, we are facing the opposite situation compared with S-NDR systems: Here, the autocatalytic feedback loop (driven by electrical variables) is fast, while the negative inhibition occurs on the slower time scale of chemical reactions. Consequently, oscillatory instabilities occur easily. Below, we demonstrate that the interaction of

10.1021/jp072691s CCC: $37.00 © 2007 American Chemical Society Published on Web 08/18/2007

13482 J. Phys. Chem. C, Vol. 111, No. 36, 2007

Siegmeier et al.

both subsystems gives rise to much enlarged oscillatory parameter regions. The paper is organized as follows: The mathematical model is introduced in section II. In section IIIA, the dynamic behavior of the system in parameter space is characterized with various one- and two-parameter bifurcation diagrams, and in section IIIB, the mechanistic origins of bistability and different types of oscillations are discussed. A summary of the results and conclusions are given in section IV.

kccbCOcCO(0.99 - θCO - θOH) - krθOHθCO efφDL (7) des reac ) θ˙ OH ) νads OH - νOH - ν

ka efφDL(θmax OH (1 - θCO) - θOH) kdθOH e-fφDL - krθOHθCO efφDL (8)

CObulk T COsurf

(1)

with f ) (RF/RT), where R is the transfer coefficient, F, R, T have their usual meanings, and φDL is the double layer potential. θOH is limited to 0 < θOH < θmax OH < 1. The potential-dependent charge-transfer processes of CO2 formation and OH adsorption and desorption are assumed to be of the Butler-Volmer type, with rate constants kr, ka, and kd, respectively. The potential is measured with respect to the reversible hydrogen electrode, thus implicitly including the pH dependence of the reaction rates.11 The electrical properties are determined by the evolution equation for the double layer potential φDL:

COsurf + / f COads

(2)

Cφ˙ DL ) -ireac + imig )

II. Model The description of the H2-CO system used here is based on the Langmuir-Hinshelwood mechanism for CO electrooxidation on platinum18,19 and takes into account that CO electrooxidation is mass transport limited under certain conditions.11 Thus, we consider the following reaction steps:

des reac ) + iH2] + -[StotF(νads OH - νOH + ν

H2O + / T OHads + H+ + e-

(3)

COads + OHads f CO2 + H+ + e- + 2/

(4)

Equation 1 describes the transport of CO from the bulk electrolyte to the WE surface, where it may adsorb onto free sites denoted by an asterisk (eq 2). CO adsorption is assumed to be practically irreversible. Reversible oxidative dissociation of water (eq 3) yields an adsorbed oxygen-containing species OHads, which can react with adsorbed CO to form carbon dioxide (eq 4). Additionally, hydrogen is oxidized on free surface sites:

H2 + 2/ T 2H+ + 2e- + 2/

(5)

Transport of H2 to the electrode surface and hydrogen adsorption preceding its oxidation are not explicitly taken into account (see below and ref 20). Equations 1-5 can be cast into the following set of differential equations, the first three of which being identical to the ones introduced by Koper et al.11 to model CO electrooxidation on Pt. Assuming that a linear concentration profile across a diffusion layer of constant thickness δ adjusts instantaneously,21 the concentration of CO in front of the WE surface develops according to

c˘ CO ) -

reac θ˙ CO ) νads ) CO - ν

2Stot c 2D k cCO(0.99 - θCO - θOH) + 2 (cbCO - cCO) (6) δ δ

Here, Stot and cbCO are the total number of surface sites per square centimeter and the CO bulk concentration, respectively; θCO and θOH denote the fractions of the surface covered by adsorbed CO and OH, and kc is the CO adsorption rate constant. We neglect the very slow desorption of CO and set the maximum coverage to 0.99, so that there are always some sites available for OH and H2 adsorption. The evolutions of the surface coverages, θCO and θCO, are modeled by

1 (U - φDL) (9) wF

where C is the specific electrode capacitance, and the product of uncompensated cell resistance and electrode area is expressed through the specific resistance F and a cell geometry factor w. The last term of eq 9 assumes potentiostatic control conditions:

U ) φDL + i‚(wF)

(10)

The hydrogen current density iH2 is taken from Plenge et al.,20 where the hydrogen oxidation reaction (HOR) in the presence of Cu2+ and Cl- is examined. In particular, the authors showed that it is not necessary to treat the surface coverage of hydrogen and its concentration in front of the surface as variables. Instead, the hydrogen current density can be modeled as a function of φDL, fitted to simulations of the full model, and a term depending on the fraction of free surface sites, θ*:

(

iH2 ) jH 1 -

)

2 g(θ*) 1 + eγφDL

(11)

In the absence of adsorbates, g(1) ) 1, and iH2 takes on the diffusion-limited current density jH for φDL larger than about 0.05 V. The expression in brackets on the right hand side of eq 11 accounts for a realistic description of the (partially) reaction controlled oxidation current at low potentials. For small and intermediate values of the total coverage, the dependence of the hydrogen current on the number of free surface sites is approximated well by a linear term.20 In our case, however, the CO coverage may attain values close to 1, and the fact that hydrogen adsorption requires two free surface sites cannot be neglected without qualitative errors. We thus chose the coverage dependent term

g(θ*) ) θ*2 ) (1 - θCO - θOH)2

(12)

We also checked that the results comply with a model explicitly taking into account the hydrogen concentration in front of the electrode and the hydrogen coverage, similar to the strategy used by Plenge et al.20 The contributions to the faradaic current from the adsorptiondesorption of OH and the CO2-formation reaction are negligible compared with the hydrogen current for physical parameter settings and not crucial for the dynamic phenomena examined here (see Figure 1). Thus, for the qualitative analysis of the

Electrooxidation of H2-CO Mixtures on Pt

J. Phys. Chem. C, Vol. 111, No. 36, 2007 13483 TABLE 1: Typical Parameter Values description

Figure 1. Location of saddle node (region I) and Hopf bifurcations (region II) in the F-U parameter plane for the full model (dashed line), the reduced model (solid line), and a two variable model in which cCO and φDL were adiabatically eliminated (dotted line). For the remaining parameters, see Table 1.

system, the respective terms in eq 9 can be omitted. In the following, this reduced system will be used unless otherwise noted. However, there are considerable deviations at higher electrolyte resistances; therefore, the full system has to be employed wherever quantitative statements are made. See the next section for some more comments about the limitations of the simplifications. The equations were analyzed using the bifurcation analysis software AUTO.22 Simulations were performed using the Livermore Solver for Ordinary Differential Equations (lsode).23 III. Results and Discussion A. Bifurcation Analysis. The dynamic behavior of the model (eqs 6-9) was characterized by one- and two-parameter bifurcation diagrams, using the external voltage U, the specific electrolyte resistance F, the diffusion layer thickness δ, the specific double layer capacitance C, and the percentage of CO in the H2-CO gas mixture as the key parameters. In a two-parameter bifurcation diagram for F and U with typical values for all other parameters, two disconnected regions of qualitatively different and nontrivial dynamic behavior can be distinguished (Figure 1): For high conductivity and low values of U, a bistable region bounded by two saddle-node (SN) bifurcations exists. Toward higher values of F, the region narrows and terminates in a cusp point (region I in Figure 1). If both U and F are higher, we find an oscillatory regime whose extension approximately coincides with the location of the Hopf bifurcations (HB) (region II in Figure 1). In contrast to the bistable region, the oscillatory region becomes wider when F is increased.In Figure 1, the saddle-node and Hopf bifurcation lines are shown for three different models: (i) the full system, in which all faradaic contributions to the reaction current in eq 9 are considered (long dashed line); (ii) the reduced system, that is, the standard system considered in the following, with the hydrogen current as the only contribution to the reaction current (solid line); and (iii) a further simplified system in which two variables, namely, cCO and φDL, were adiabatically eliminated (dotted line). In the bistable area (I), the different curves can hardly be distinguished; all three models are in quantitative agreement. For the oscillatory region (II), the location of the Hopf bifurcations in the F-U parameter plane match qualitatively; however, we observe quantitative differences: Coming from low voltages, the first Hopf bifurcation of the two reduced

transfer coefficient Faraday constant gas constant temperature diffusion coefficient of CO diffusion layer thicknessb rate constant of OH adsorption rate constant of OH desorption rate constant of CO adsorption rate constant of CO2 formation hydrogen current constantc maximum OH coverage fit parameterc total no. of surface sites CO bulk concentrationd spec. double layer capacitance spec. electrolyte resistance cell geometry parameter

parametera R F R T D δ ka kd kc kr jH θmax OH γ Stot cbCO C F w

value 0.5 96484 8.314 293 5 × 10-5 0.0067 10-4 105 108 10-5 1.226 × 10-2 0.333 118.6812 2.2 × 10-9 2 × 10-8 2 × 10-5 bif. param. 4

unit C mol-1 J mol-1 K-1 K cm2 s-1 cm s-1 s-1 cm3 mol-1 s-1 s-1 A cm-2 V-1 mol cm-2 mol cm-3 F cm-2 Ω cm cm

a Parameter values from 11. b Corresponding to 100 rpm. c From ref 20. d Hydrogen gas with 2% CO.

systems occurs at lower values of U than in the full model, and the deviation between reduced and full models is only negligible at low electrolyte resistances F. First, let us consider the bistable region and its parameter dependence in some more detail. Figure 2a shows a typical oneparameter bifurcation diagram in which the stationary current is plotted versus U. Clearly, there is a voltage interval in which stable low- and high-current states coexist, corresponding to high and low CO coverage of the electrode surface, respectively. These two stable branches “sandwich” an unstable branch, corresponding to an intersection of the load line with the I-φDL characteristic in a region of negative differential resistance (NDR). The fixed points are created and destroyed in SN bifurcations. As explained in detail in ref 14, increasing the resistance narrows the bistable voltage interval such that above a critical value the system is monostable (cf. Figure 1). The location of the bistable region in the F-U parameter plane depends on all system parameters but the capacitance C, which only influences the time scale of the φDL equation and thus does not affect the location of stationary states in phase space (although it might change its stability). The influence of the CO bulk concentration on the location of the SN bifurcations in the F-U parameter plane is shown in Figure 2b. Roughly, they shift toward smaller values of U with decreasing cbCO, and correspondingly, increasing the fraction of CO in the H2-CO mixture shifts the bistability toward more positive voltages, up to the limit given by the pure CO subsystem (obtained by using the full model with eq 9 and setting iH2 to zero). Figure 2c demonstrates that the bistability depends crucially on δ: The thinner the diffusion layer, that is, the more efficient the mass transport, the narrower the bistable voltage interval and the lower the critical resistance at which the cusp occurs. Above a critical transport rate, the system has only a single steady state, irrespective of the other parameter values (as long as they remain within physically meaningful boundaries). Next, let us focus on the parameter dependence of the oscillatory region. The one parameter bifurcation diagram shown in Figure 3 depicts the current density versus U in a voltage interval covering the two Hopf bifurcations at F ) 40 Ω cm and reveals that both Hopf bifurcations are subcritical at these parameter values. The unstable limit cycles annihilate with stable partners in saddle-node bifurcations of periodic orbits (SNP), such that the actual parameter region in which stable oscillations exist is somewhat larger than region II of Figure 1.

13484 J. Phys. Chem. C, Vol. 111, No. 36, 2007

Siegmeier et al.

Figure 2. Typical bifurcation diagrams in the bistable parameter range: (a) Bistable current density vs external voltage U as parameter. Solid line, stable fixed points; dashed line, unstable fixed points. (F ) 2 Ω cm). (b) Location of SN bifurcations and cusp points in the two-parameter F-U plane for different CO concentrations. (c) Location of SN bifurcations and cusp points in the two-parameter F-U plane for different diffusion layer thicknesses δ. For the remaining parameters, see Table 1.

Figure 3. One parameter bifurcation diagram through the oscillatory region at F ) 40 Ω cm. Thick (thin) solid line: stable (unstable) stationary states. Dots: limit cycles (unstable and stable). For the remaining parameters, see Table 1.

Again, there are significant quantitative differences between the full and the reduced model: The voltage interval between SNP and HB is considerably wider in the full system than in the reduced system. Furthermore, in the full system, secondary bifurcations of the periodic orbit, such as period doublings, were found. The simplified system apparently does not reproduce these features, but it is sufficient to describe the bifurcations of the steady state and thus more suitable to uncover the mechanism producing the basic oscillations, which is our main concern here. The location of the Hopf bifurcation in various two-parameter planes is shown in Figure 4. The variation of the oscillatory region with the CO concentration is most important from the point of view of possible applications. As can be seen in Figure 4a, a decrease in the CO concentration shifts the oscillatory region toward lower voltages and resistances; that is, it enlarges the oscillatory parameter interval. Thus, oscillations are favored at experimentally more important parameter values. The location of the Hopf bifurcation set in the δ-U plane again reveals that a minimum diffusion layer thickness is necessary for oscillations to occur (Figure 4b). This suggests that either a slow mass transfer is a necessary ingredient of the oscillation mechanism or that the possibly constant CO concentration in front of the electrode must not exceed an upper threshold. As we will discuss below, this is an important hint to the mechanistic steps essential for the Hopf instability. In this respect, it is also useful to consider the dependence of the width of the interval where oscillations exist on the electrode capacitance, which determines the relative speed of the evolution equations of the electrode potential and the chemical species. From Figure 4 (c and d), two important points can be inferred: First, there is no bound

for the double layer capacitance C limiting the oscillatory region, which is very unusual for electrochemical oscillators. Second, there is a capacitance threshold, about 2 orders of magnitude larger than the capacitance of a Pt electrode, above which the oscillatory region extends to lower values of the external voltage U. B. Feedback Loops and Oscillation Mechanisms. Having discussed the parameter dependence of the bistability and the oscillations, we can use this information to deduce the essential feedback mechanisms that generate saddle node and Hopf instabilities. Chemical Autocatalysis and Bistability. Focusing first on the bistable region, we saw that the saddle node bifurcations continuously shift in the F-U parameter plane when the H2CO ratio changes (cf. Figure 2b). The same is true for constant CO bulk concentrations if the hydrogen bulk concentration is lowered (experimentally, H2 could be gradually replaced by an inert gas like Ar). Furthermore, a crucial dependence of the instability on the diffusion layer thickness is revealed (cf. Figure 2b). We thus conclude that the underlying feedback mechanism is the same as in the pure carbon monoxide system described in refs 11 and 12, namely, a chemical autocatalysis arising from the Langmuir-Hinshelwood mechanism coupled to mass transport limitations: Consider a CO poisoned state, such as the lowcurrent branch in Figure 2a. Owing to the negligible oxidation rates, the concentrations of the reactive species in the reaction plane (in particular of CO) are identical to the bulk concentrations. However, when increasing U, hydroxide adsorption eventually sets in and OHads reacts with COads. This yields free surface sites, for which CO and OH molecules compete. As the adsorption rate of CO grows due to the larger number of free surface sites, CO adsorption enters the transport limited regime, while the oxidative adsorption of water remains reaction controlled. Thus, with time, OH adsorption is favored more and more (compared with CO adsorption), which eventually leads to an avalanche-like stripping of COads. In this way, the system is driven to the reactive state. If we decrease U, the carbon monoxide coverage only increases again when the adsorption rate of OH is sufficiently low. This also reduces the CO oxidation and the CO adsorption rate with it. A lower adsorption rate, in turn, leads to a replenishment with CO in the reaction plane. The ratio between CO and OH adsorption rates increases further until the electrode is practically completely covered by CO and the poisoned state is attained again. The presence of hydrogen increases the total oxidation current on the reactive branch and thus also the parameter values at which the system is bistable in the F-U parameter plane, but it does not affect the positive feedback loop. So we see that the bistability is generated by an S-NDR mechanism.

Electrooxidation of H2-CO Mixtures on Pt

J. Phys. Chem. C, Vol. 111, No. 36, 2007 13485

Figure 4. Location of Hopf bifurcations in two-parameter spaces. (a) Specific electrolyte resistance F vs external voltage U for different CO bulk concentrations. (b) Diffusion layer thickness δ vs U for F ) 100 Ω cm. (c) Double layer capacitance C vs U. (d) Specific resistance F vs U for different double layer capacitances. The dashed lines enclose the parameter region in which the steady state lies on the unstable, autocatalytic branch of the (IR-corrected) polarization curve. For the remaining parameters, see Table 1.

Oscillations. In principle, the chemical mechanism described above can also lead to instabilities outside the bistable regime via Hopf or Turing bifurcations. As a precondition, the steady state has to be on the autocatalytic middle branch of the i-φDL curve. We numerically determine the current densities and potential values at the SN bifurcations for vanishing electrolyte resistance (i1, φDL,1 and i2, φDL,2) and combine this with the potentiostatic control condition. The location of the autocatalytic branch in the F-U space is then given by the two lines F ) (U - φDL,1/2)/i1/2w, shown in Figure 4d along with with the Hopf bifurcation line for three different values of the electrode capacitance. Three important features can be observed: First, regardless of the capacitance, the Hopf bifurcation at high external voltages does not lie on the autocatalytic branch but lies on the reactive branch. Close to the low resistance end of the oscillatory region and for low (that is, realistic) values of C, both Hopf bifurcations occur on the reactive branch. This is exemplified in Figure 5, where the one-parameter bifurcation diagram of Figure 3 is shown after IR correction (R denotes the total cell resistance, and I denotes the total current). Consequently, these oscillations must arise from feedback loops not related to the Langmuir-Hinshelwood mechanism. Second, at large values of C, the Hopf set at low voltages coincides with the low-voltage end of the autocatalytic branch. Since a large capacity shifts the time scale of φDL (which is an inhibitor in S-NDR systems) to larger values, the low voltage oscillations occur at parameter values typical for S-NDR oscillations. We thus expect the system to possess two different mechanisms that generate oscillations. Third, since for low capacitance the oscillations partially occur on the autocatalytic branch, while at high values of C they are also observed around the reactive branch, it is likely that there are parameter regions in which both oscillation mechanisms come into effect.

Figure 5. IR corrected one-parameter bifurcation diagram with U as main bifurcation parameter. Parameter values are identical to those of Figure 3. Solid (dotted) line: stable (unstable) stationary states.

To shed more light on the oscillation mechanisms, we first note that there are three types of oscillations with distinct phenotypes. Their approximate location in the C-U parameter plane is indicated in Figure 6. The vertical lines (a1 and a2) limit the voltage region in which the steady state lies on the autocatalytic branch. Region A: HN-NDR Oscillations. Typical time series and phase portraits of the low capacitance oscillations of region A are depicted in Figure 7a,b, respectively. Figure 7b shows θCO, θOH, hydrogen oxidation current and CO concentration versus φDL, and additionally, the polarization curve and the load line. The intersection of the latter is the (unstable) steady state, which

13486 J. Phys. Chem. C, Vol. 111, No. 36, 2007

Figure 6. Regions with oscillations of different phenotype in the C-U parameter plane. a1-a2: voltage interval, in which the steady state lies on the autocatalytic branch of the polarization curve (F ) 100 Ω cm).

Figure 7. Typical time series (a) and phase portraits (b) of HN-NDR type oscillations in region A of Figure 6. In b, also the (IR corrected) polarization curve and the load line are shown. C ) 20 µF cm-2, F ) 100 Ω cm, U ) 1 V, and for other parameters, see Table 1.

clearly lies on a branch with positive differential resistance. However, the hydrogen oxidation current decreases during the oscillation, while φDL increases and vice versa. This behavior is anticorrelated with the total coverage of the electrode and suggests that the oscillations are of the HN-NDR type. The double layer potential is the activator variable and the positive feedback loop stems from a hidden negative differential resistance of the polarization curve coupled to the potentiostatic constraint.

Siegmeier et al. We follow one oscillatory cycle (7a): Starting at low double layer potential (region 1), potential-independent CO adsorption increases the coverage slowly and linearly, while there is no hydroxide on the surface. The decreasing number of free surface sites prompts a decrease of the reaction current, which, in turn, entails an increase of φDL. Eventually, φDL is high enough for OH adsorption to set in (region 2). A growing θOH leads to a decrease in θCO due to CO2 formation. At first, the decline of θCO is overcompensated by the growing hydroxide coverage, leading to even higher φDL. However, this also enhances the oxidation rate of COads; therefore, at some point, the total coverage drops and initiates a reverse avalanche-like response (region 3). The shift of φDL toward smaller values not only reduces the OH adsorption rate but also ensures that OH eventually desorbs completely. The number of free surface sites again increases drastically, such that φDL is driven toward quite low values, from where the cycle starts anew (region 4). From Figure 7, it is apparent that cCO changes only by about 2% around its mean value during one oscillatory cycle. In fact, if we set the CO concentration to a constant value close to its mean value during the oscillations, the system remains oscillatory, though the S-shaped polarization curve (which stems from the chemical autocatalysis) turns into a sigmoid. This is additional evidence for the HN-NDR mechanism discussed above. Region B: S-NDR Type Oscillations Assisted by a PositiVe Feedback in φDL. Now, let us turn to the oscillations in region B in Figure 6. For the larger part of this region, the steady state lies on the middle branch of the S-shaped current potential curve. Furthermore, oscillations only occur if the capacitance is at least double that of region A. This is in accordance with classical S-NDR oscillations, where the autocatalysis is of a chemical nature, in our case stemming from the LH mechanism, with the double layer potential as the slow, inhibitory variable. However, this interpretation turns out to be too simple: We compare the time series and phase portraits of region B oscillations to those of the pure CO system, which is a classical S-NDR oscillator (Figure 8). In the latter system, φDL changes slowly, thereby driving the system back and forth between the reactive and the poisoned state, while the activatory subsystem, that is, θCO, θOH, and cCO, adjust quickly to the potential changes (Figure 8a). This is particularly apparent in the current-potential representation (Figure 8b), where the (IR corrected) polarization curve and load line are shown along with the hydrogen current density for the oscillatory state of Figure 8a. The latter follows the S-shaped polarization curve closely (except for the overshooting after the transition from the poisoned to the reactive state; see inset), which also implies that the amplitude of the electrode potential coincides approximately with the width of the hysteresis. For oscillations in region B (Figure 8c,d), the picture is different: Considering the time series, we note that the section where φDL increases is reminiscent of the corresponding one in the pure CO system, while the decrease occurs on a much faster time scale. The minimum potential falls well below the lowpotential limit of the reactive branch (Figure 8c,d). The origin of this behavior can be seen in the phase space representation of Figure 8d. The presence of hydrogen increases the current density of the reactive branch by 2 orders of magnitude. The changes in φDL are proportional to the difference between the faradaic current and the load line at the respective electrode potential; therefore, potential changes are about 100 times faster; their characteristic time scale is shorter than those of the chemical variables. As a consequence, the system does not relax

Electrooxidation of H2-CO Mixtures on Pt

J. Phys. Chem. C, Vol. 111, No. 36, 2007 13487

Figure 8. Time series and phase portraits of typical S-NDR oscillations (a,b) and region B oscillations (c,d; see Figure 6). (a,b: 2% CO and 98% inert gas.) Additionally in b,d, polarization curve (IR corrected) and load line. C ) 0.2 F cm-2, F ) 100 Ω cm, U ) 0.66, and for other parameters, see Table 1.

to the poisoned state when the potential has dropped below the value where the reactive branch vanishes. Instead, φDL is driven to even lower values, where the hydrogen current density is falling slowly because of an increasing overall coverage and slows down changes of φDL as well. Now, the chemical selfenhancing process described above can set in again: At low potential, OH desorbs and the CO adsorption outweighs the OH adsorption, which in turn causes a considerable decrease of the CO2 formation rate. The CO concentration in front of the electrode increases, promoting the buildup of the CO coverage. This chemical autocatalysis is reflected by the fast change in cCO. Obviously, it also occurs close to the potential maximum. On the other hand, the fast move into the negative region is accompanied by an autocatalytic, HN-NDR type feedback, showing up in the increasing current density with decreasing potential in the bistable region of the polarization curve. Thus, we observe a peculiar interplay between two positive feedback loops. The presence of the HN-NDR type feedback explains the existence of oscillations beyond the voltage range of the autocatalytic branch (cf. Figure 6). The interpretation above is supported by studies using a reduced model. On the one hand, the oscillations do not exist when cCO is fixed at a value within its oscillation interval. Thus, in contrast to region A oscillations, for region B oscillations, changes of cCO are essential. However, cCO, as well as OHads, can be eliminated adiabatically. Then, we are left with a two variable system, and the sign matrix of the Jacobian matrix (evaluated around the steady state) immediately yields information about the dominant feedback loops causing the primary instability; however, off the steady state, that is, on the limit

cycle, both of them are effective and interact in the manner described above. Close to the low potential border, the sign matrix is the same as what we would expect for an S-NDR oscillator: Using the evolution equation of φDL as the first equation and the one for θCO as the second, we obtain:

sign( J ) )

( )

- + r φ˙ DL - equation - + r θ˙ CO - equation

(13)

Whenever it is on the reactive branch, the steady state is of the HN-NDR type:

sign( J ) )

( )

+ + r φ˙ DL - equation - - r θ˙ CO - equation

(14)

A hybrid matrix with two positive diagonal elements is found close to the high-potential border of the chemically autocatalytic branch, that is, close to U ) a2:

sign( J ) )

( )

+ + r φ˙ DL - equation - + r θ˙ CO - equation

(15)

Region C: Hybrid Oscillations Resulting from Two PositiVe Feedback Loops. The oscillations in region C (Figure 6) are best described as hybrids between those in region A and those in region B. In the phase portrait representation (Figure 9b), the trajectories for high values of φDL resemble those of region B: The fast change of cCO close to the maximum of φDL, accompanied by large changes in θCO and θOH (Figure 9a), at

13488 J. Phys. Chem. C, Vol. 111, No. 36, 2007

Figure 9. Typical time series (a) and phase portraits (b) of hybrid type oscillations of region C of Figure 6. In b, also the (IRcorrected) polarization curve and the load line are shown. C ) 0.03 F cm-2, F ) 100 Ω cm, U ) 1.26V, and for other parameters, see Table 1.

first seem to hint at the chemical autocatalysis described above. However, since the capacitance is smaller than in region B, φDL decreases faster, enhancing OH desorption. During the time interval where θOH drops, the CO coverage has resumed linear (i.e., diffusion controlled) growth, while cCO remains at a low value. This shows that the OH desorption is driven only by the potential changes and not by a changing rate of CO and OH adsorption; the chemical autocatalysis has no effect here. Instead, the interplay of fast OH desorption (decreasing total coverage) and decreasing φDL goes on until the electrode surface is virtually free of OH. After another short interval, when φDL is at its minimum and CO coverage is sufficiently high, the hydrogen current cannot be maintained anymore. Now, as for the region A oscillations, φDL increases again. The shape of the oscillation is thus again determined by the interaction between chemical and electrical positive feedback loops. The difference to region B oscillations is that the chemical autocatalytic cycle is only in effect close to the high potential limit of the oscillations. As before, the oscillations survive under adiabatic elimination of cCO and θOH (and preserve their characteristic form) but vanish when cCO is set to a constant value within its oscillation interval. The sign matrix of the Jacobian evaluated at the steady state reveals that the dominant destabilization is due to the positive feedback in φDL. IV. Conclusions Our theoretical analysis of a homogeneous model for the oxidation of H2-CO mixtures on Pt electrodes revealed dynamical instabilities in wide parameter regions. The instabili-

Siegmeier et al. ties are due to a chemical autocatalysis, where the electrode potential acts like an inhibitor, and an electrical autocatalysis with the electrode potential as the activator variable, while chemical reactions provide a negative feedback. The system thus displays the characteristic features of both S-NDR and HNNDR systems. For some parameter settings, one of the two possible interaction schemes is active and thus determines the dynamics. For example, at low cell resistance and low voltages, the chemical autocatalysis renders the system bistable, whereas at higher resistance, higher voltage, and low capacitance, the system exhibits HN-NDR oscillations. In other parameter regions, however, the dynamic behavior is clearly determined by the interplay of both activator-inhibitor type interactions, a hitherto unique behavior in an electrochemical system. As for the reported experimental instabilities, the chemical autocatalysis most likely accounts for the hysteresis observed during cyclic voltammetry,18,19 while experimentally reported oscillations with bulk Pt, and thus a low capacitance electrode, should be due to an HN-NDR type instability.4 For other experimental conditions, the situation is not as clear, and careful analysis is required. So far, we only considered the lumped, spatially uniform system. Since S-NDR systems are very likely to exhibit spatial instabilities, it can be anticipated that wide parameter ranges with a variety of spatio-temporal patterns exist, for example, arising from the interaction of Hopf and Turing bifurcations. Examples of the latter are still rare in physical systems. The theoretical results also call for experimental validation. In this respect, it is most desirable to prove the existence of the three different types of oscillations depicted in Figures 6-9. Here, the electrode capacitance is a crucial parameter. Typical capacitances of metal electrodes are of the order of 10 µF/cm2. Large values can be realized by polymer coated electrodes or electrodes consisting of Pt particles on some inert support. An experimentally much simpler access to large values of the electrode (pseudo-)capacitance is to employ a tunable differential controller as proposed in ref 24. Such studies are currently being carried out in our laboratory. The mentioned experimental and theoretical extensions of the studies presented here form the basis for future applications, in which dynamic instabilities keep the poisoning of the electrode by adsorbed CO to a minimum and thus make anodes more CO tolerant. Acknowledgment. We thank P. Bauer for fruitful discussions and for proof reading the manuscript. This work was supported by the Deutsche Forschungsgemeinschaft (Project No. KR1189/10). N.B. acknowledges receipt of a fellowship through the HWPII program. References and Notes (1) Markovic, N. M. In Handbook of Fuel Cells Vielstich, W., Lamm, A., Gasteiger, H. A., Eds.; John Wiley & Sons Ltd.: Chichester, U.K., 2003; Volume 2, p 368. (2) Deibert, M. C.; Williams, D. L. J. Electrochem. Soc. 1969, 116, 1290. (3) Szpak, S. J. Electrochem. Soc. 1970, 117, 1056. (4) Yamazaki, T.; Kodera, T. Electrochim. Acta 1990, 36, 639. (5) Zhang, J.; Datta, R. J. Electrochem. Soc. 2002, 149, A1423. (6) Zhang, J.; Datta, R. Electrochem. Solid State Lett. 2004, 7, A37. (7) Zhang, J.; Datta, R. J. Electrochem. Soc. 2005, 152, A1180. (8) Zhang, J. X.; Fehribach, J. D.; Datta, R. J. Electrochem. Soc. 2004, 151, A689. (9) Kiss, I. Z.; Bracket, A. W.; Hudson, J. L. J. Phys. Chem. A 2004, 108, 14599. (10) Liu, P.; Logadottir, A.; Nørskov, J. K. Electrochim. Acta 2003, 48, 3731. (11) Koper, M. T. M.; Schmidt, T. J.; Markovic, N. M.; Ross, P. N. J. Phys. Chem. B 2001, 105, 8381.

Electrooxidation of H2-CO Mixtures on Pt (12) Strasser, P.; Eiswirth, M.; Ertl, G. J. Chem. Phys. 1997, 107, 991. (13) Krischer, K. J. Electroanal. Chem. 2001, 501, 1. (14) Mazouz, N.; Krischer, K. J. Phys. Chem. B 2000, 104, 6081. (15) Li, Y. J.; Oslonovitch, J.; Mazouz, N.; Plenge, F.; Krischer, K.; Ertl, G. Science 2001, 291, 2395. (16) Bonnefont, A.; Varela, H.; Krischer, K. J. Phys. Chem. B 2005, 109, 3408. (17) Krischer, K.; Varela, H. In Handbook of Fuel Cells Vielstich, W., Lamm, A., Gasteiger, H. A., Eds.; John Wiley & Sons Ltd.: Chichester, U.K., 2003; Volume 2, p 679-701. (18) Gasteiger, H. A.; Markovic, N. M.; Ross, P. N. J. Phys. Chem. B 1999, 103, 9616. (19) Gasteiger, H. A.; Markovic, N. M.; Ross, P. N. J. Phys. Chem. 1995, 99, 16757.

J. Phys. Chem. C, Vol. 111, No. 36, 2007 13489 (20) Plenge, F.; Varela, H.; Lu¨bke, M.; Krischer, K. Z. Phys. Chem. 2003, 217, 365. (21) Koper, M. T. M.; Sluyters, J. H. J. Electroanal. Chem. 1991, 308, 151. (22) Doedel, E. J.; Kevernez, J. P. AUTO: Software for continuation and bifurcation problems in ordinary differential equations applied mathematics report; California Institute of Technology: Pasadena, CA, 1986. (23) Hindmarsh, A. C. In Scientific Computing; Stepleman, R. S., Ed.; North-Holland: Amsterdam, 1983; Vol. 1, p 55. (24) Kiss, I. Z.; Kazsu, Z.; Gaspar, V. J. Phys. Chem. A 2005, 109, 9521.