J. Phys. Chem. 1992,96,2676-2681
2676
Modeling Periodk and Chaotic Dynamics in Anodic Nickel Dissolution D. Haim: 0. Lev,* L. M. Pismen2 and M. Sheintuch**+ Department of Chemical Engineering, Technion-Israel Institute of Technology, Haifa 32000, Israel, and Division of Environmental Sciences, School of Applied Science and Technology, The Hebrew University, Jerusalem 91 904, Israel (Received: August 28, 1991)
A systematic approach for modeling and parameters estimation of the motions and singularities observed in anodic nickel dissolution in sulfuric acid is presented. The observed and modeled behaviors include bistability, periodic oscillations bounded by Hopf and saddle loop bifurcations, the disappearance of these transitions in a double-zero singularity, and the existence of chaotic solutions which are bounded by a period-doubling transition. The modeling and parameter estimation are based on comparison of the experimental and simulated bifurcation maps, showing domains of different characteristic behaviors in the plane of applied constant current vs acid concentration. Qualitative agreement between the observed and predicted maps has been obtained.
Introduction Periodicity and dynamic instabilities in electrochemical systems date back to the first observations of periodic deposition and dissolution in iron wires by Fechner in 1828.' Electrochemical instabilities may be easily monitored by simple current or voltage measurements, and parametric changes of a few orders of magnitude can be easily applied to electrochemical systems by current or potential control. Moreover, the overall reaction rate can be directly manipulated by current control. Nonlinear electrochemical systems were found to be most valuable as models for more complicated chemical or physiological nonlinear systems. In fact, the first neuronic excitation models" were based on the activity/passivity fronts which were observed in anodized iron wires. That subject has been reviewed by Heathcoat as early as 1907.5 This system still serves as a most powerful model for information transfer in mammals6 Recently, a plethora of nonlinear dynamics associated with passivating metal has been discovered.' The lumped current or potential oscillationsare often found to be associated with spatial structural changes of the electrode ~ u r f a c e . ~ >Despite ~ the fact that electrochemical instabilities serve as models for other more complex systems, it becomes clearer now that the underlying mechanism of the nonlinear electrochemical phenomena and the transitions between them are not fully understood. Academic interest in mathematical modeling of nonlinear electrachemical systems has been recently revived, and mechanistic models that explain periodicity and even chaos have appeared in the literature.iO*ii The approach taken here for studying electrochemical systems follows singularity theory. A rational approach for organizing experimental and theoretical analysis of dynamic behaviors is to draw bifurcation diagrams, showing domains of different character as one parameter is varied, or bifurcation maps, where the domains are presented in the parameter plane. Identifying singular (transition or bifurcation) points in the experimental diagram or in the map suggests a systematic approach for modeling the observed behavior. This approach is novel to electrochemical reactions; furthermore, the observed behavior is richer than that reported for most catalytic oscillators and better resolution has been achieved. Thus the present modeling problem, based on the singular points, is more challenging than that in other heterogeneous systems. The galvanostatic dissolution of nickel in sulfuric acid is known to exhibit a plethora of nonstationary states including steady-state multiplicity; periodic, multiperiodic, quasiperiodic, and chaotic potential oscillations; Hopf, saddle loop (infinite period), and semistable limit-cycle transitions to periodicity; period doubling and quasiperiodic transitions to chaos; and higher order singularities (codimension two). These states and singularities were identified in a systematic way and mapped in bifurcation diagrams
'
Technion-Israel Institute of Technology. *The Hebrew University.
and maps.12 Moreover, measurements of the local current distribution near the anodized nickel electrode by an array of microreference electrodesi3revealed that the nonstationary steady states are associated with spatial symmetry breaking, including antiphase oscillations and propagating wave fronts.I4 Switching between these types of behavior is possible by a simple shift from a constant current to a constant voltage control of the electrochemical cell. Accounting for such a rich dynamic behavior by a single model calls for a certain strategy. Here we adopt the simplest approach by which we build a two-variable lumped model that accounts for a simple periodicity and all its singularisies and then relax one of the assumptions to obtain a three-variable model that admits chaos. In a future publication we shall demonstrate that the same kinetic model accounts for the spatial patterns. Construding a model and finding a set of parameters that induce chaos or admits a global (in the phase plane) bifurcation is not a straightforward task except in the vicinity of certain singuiarities or when the time scales of the variables are widely different. It calls, therefore, for a trial and error search. Parameter estimation for the twovariable model is conducted in a systematic way and the methodology is outlined. We outline our considerations in building the three-variable model and in searching for the parameters that admit chaotic behavior. A more sophisticated approach would consider first whether chaotic solutions and certain singularities (like the double-zero, (1) Fechner, G. Th Schweigg J. Chem. Phys. 1828.53, 129. (2) Bonhoeffer, K. F. Gem. Physiol. 1948, 32, 69. Bonhoeffer, K. F. Corrosion 1955, 1 1 , 304. (3) Bonhoeffer, K. F. J . Vetter Z . Phys. Chim. 1950, 196, 127. (4) Lillie, R. S.J . Gen. Physiol. 1925, 7 , 473; Biol. Rev. 1936, 1 1 , 181; Science 1918, 48, 51. ( 5 ) Heathcoat, H. L. J . SOC.Chem. Ind. 1907, 26, 899. (6) Chizmadzhev, Yu. A.; Pastushenko, V. F. In Comprehensive Studies
of Electrochemistry; Bockris, J. OM., White, R. E., Eds.; Plenum Press: New York, 1985; Vol. 10, p 381. (7) Wojtowicz, J. In Modern Aspects of Electrochemistry; Bockris, J. O M . , Conway, B. E., Eds.; Plenum Press: New York, 1972; Vol. 8, p 47. Bassett, M. R.; Hudson, J. L. J. Phys. Chem. 1989,93, 2731; Physica D 1989, 35, 289. Lee, H. P.; Nobe, K.; Pearlstein, A. J. Electrochem. SOC.1985, 132, 1031.
(8) Tsitsopoulos, L. T.; Tsotsis, T. T.; Webster, I. A. Surf. Sci. 1987,191, 225. (9) Bassett, M. R.; Hudson, J. L. J . Electrochem. SOC.1990, 137, 992; 1990, 137, 1815. (IO) Franck, U. F.; FitzHugh, R. 2.Elektrochem. 1961,65, 156. (11) Wang, Y.; Hudson, J. L. J . Electrochem. SOC.1990, 137, 485. Pearlstein, A. J.; Lee, H. P.; Nobe, K. Ibid. 1985, 132, 2159. Pearlstein, A. J.; Johnson, J. Ibid. 1989, 136, 1290. Talbot, J. B.; Oriani, R. A.; Dicarlo, M. J. Ibid. 1985, 132, 1545. St-Pierre, J.; Piron, D. L. Ibid. 1987, 134, 1689. Russel, P.; Newman, J. Ibid. 1987, 134, 1051. (12) Lev, 0.;WolMberg, A,; Seintuch, M.; Pismen, L. M. Chem. Eng. Sci. 1988,43, 1339. Lev, 0.;Wolffberg, A,; Pismen, L. M.; Sheintuch, M. J . Phys. Chem. 1989, 93, 1661. (13) Lev, 0.;Sheintuch, M.; Yarnitzky, H.; Pismen, L. M. Chem. Eng. Sci. 1990, 45, 839. (14) Lev, 0.; Sheintuch, M.; Pismen, L. M.; Yarnitzky, H. Nature 1988, 336. 458.
OO22-3654/92/2096-2616$03.00 f 0 0 1992 American Chemical Society
The Journal of Physical Chemistry, Vol. 96, No. 6, 1992 2677
Anodic Nickel Dissolution
J (CURRENT)
BO L
40N
1
1.71
wUSN
40
500
700
a0 -
900
1100
1300
1500
1700
Ch (ACID CONC.)
0-
m I( ELECTRODE POTENTIAL
0
("d)
(4
0.
4.
'
UNIQUE STABLE STATE
-
U(cY1 Figure 1. Comparison of observed (a-c, after Lev et al. Chem. Eng. Sci. 1988,43, 1339) and simulated (d-f) characteristic bifurcation diagrams in anodic nickel dissolution at constant current: solid and broken thick lines denote stable and unstable states, respectively, and thin lines show end points of oscillations (an upper oscillatory domain exists in diagram (e) in similarity to (b); its lower boundary is marked SL).
see below) are due to spatial distribution. In that case we may account for the observed behavior with a kinetically simpler m.ode1. Dynamics of Anodic Ni Dissolution Potentiostatic Dissolution. The potentiostatic dissolution of nickel in acidic media has been a subject of extensive research.15 Five distinct operating regions can be distinguished in the potential-current curve: Active dissolution exhibiting a Tafel-like dependence occurs at low potentials, followed by a sharp decrease of the current at the primary passivity region. A transpassive branch starts well before oxygen evolution becomes thermodynamically feasible. A second descending potential range (secondary passivity) commences at higher potentials. This negatively sloped branch occurred only in concentrated sulfuric acid solution and was not found at nitric, hydrochloric, or hydrofluoric acids.16 Oxygen evolution occurs at higher potentials. The potentiostatic dissolution of nickel does not exhibit any multiple or nonstationary steady states. The present model restricts attention to the phenomena that take place after the passive layer is developed (i.e., transpassivity, secondary passivity, and oxygen evolution regions). Galvanostatic Dissolution. Figure 1a-c presents typical bifurcation diagrams observed during the galvanostatic (constant current) dissolution of nickel in sulfuric acid. The diagrams describe the dependence of the potential on the applied current at a fixed acid concentration. Such diagrams were repeated at different concentrations and the identified transitions (bifurcations) were noted in a current-potential map; their loci constitute the (15) Sato, N.; Okamoto, G. J. Electrochem. SOC.1964, 111, 897. Sato, N.; Okamoto, G. Zbid. 1963,110,605. MacDougal, B.; Cohen, M. Ibid. 1977, 124, 1185. Bockris, J. OM.; Reddy, A. K. N.; Rao, B. Ibid. 1966, 113, 1133. Delichere, P.; Goff, A. H.-L.; Yu, N. J. Electrochem. SOC.1986,133, 2106. (16) Jouanneau, A.; Keddam, M.; Petit, M. C . Electrochim. Acta 1976, 21, 287.
boundaries (lines)between domains of different dynamic behaviors. These boundaries terminate by coalescence with other borders. By verifying the continuity of the potential vs current diagrams and the continuity of the bifurcations in the current vs concentration plane, we were able to capture all the dynamic features of the experimental system. In dilute sulfuric acid ( 4 . 1 N), the potentiostatic current-potential curve is monotonic, secondary passivity still has not commenced, and the steady state is always stable. With larger concentrations, the potentiostatic currentvoltage curve becomes nonmonotonic. This results in steady-state multiplicity when operating in the galvanostaticmode of operation: The secondary passivity (and the primary passivity) branch becomes unattainable and sweeping the current up and down yields hysteresis between the oxygen evolution and the transpassive branches. The region of the unattainable (saddle) states is bounded by the saddle node boundaries which emanate from the cusp point (Figure 2a, 0.1 N, 4 mA/cm2). At higher concentrations an oscillatory region bounded by two soft Hopf bifurcations evolves (Figure la). Upon increasing the current, the amplitude of the potential oscillations gradually increases and then decreases again until a second Hopf bifurcation leads to a stable, stationary portion of the transpassive branch. Further increase of the current results in a hard transition to the oxygen evolution branch and a formation of a hysteresis loop. At still higher concentrations of sulfuric acid the oscillatory branch develops by a Hopf bifurcation scenario. Upon further increase of the current the period of oscillations increases dramatically, and a saddle loop transition to the oxygen evolution region is encountered. A second oscillatory branch bounded by a saddle loop bifurcation from below and a Hopf bifurcation from above exists at higher currents (Figure lb). Upon further increase of the sulfuric acid concentration the higher oscillatory branch disappears (Figure IC). The three boundaries, the Hopf, saddle node, and saddle loop bifurcations, coaleasce at a singularity with two vanishing eigenvalues (marked DZ (double zero) in Figure 2a). As this concentration is approached, the upper stationary stable branch gradually shrinks as the bounding transitions converge. Beyond the DZ concentration the upper branch is unstable and cannot
2678 The Journal of Physical Chemistry, Vol. 96, No.6,1992 be realized. The double zero singularity is a local bifurcation and in its vicinity the oscillations are of infinitisimal amplitude. The bifurcations, including the saddle loop, can be described analytically and the behavior can be unfolded (see Guckenheimer and Holmes''). At still higher sulfuric acid concentrations (6-8 M), the galvanostatic dissolution of nickel reveals an ordered transition to chaos via quasiperiodic oscillations12(two torus transition). The chaotic region is interrupted by windows of multiperiodic frequency-locked states. The oscillatory domain in still bound by Hopf and saddle loop bifurcations. Period doubling transitions to chaos were observed in shallow diluted (1 M) solutions.I2 The global potential oscillators are accompanied by antiphase oscillations of the local current; i.e., the current at one side of the electrode oscillates exactly out of phase with current fluctuations at the other side. The global current, which is the integral of the local current, remains constant throughout. When operating at a fmed potential across the cell and with an extemal small resistor in series with the cell, the standing antiphase oscillations are converted into a train ~ a v e . ' ~ J ~
The Model The model assumes the following mechanism for the current-potential behavior of nickel dissolution. Initially, the nickel electrode is covered by an oxide film formed in the primary passivity. The small current observed during the passivity and transpassivity are due to field-assistedmigration of nickelous ions through the oxide layer. However, nickel dissolution can take place only through the surface fraction which is not covered by adsorbed species. The phase crossing of the nickel ions from the oxide into the solution may take place only when nickel ions (hydrated or not) are directly exposed to the solution. At higher surfuric acid concentrationsthe bisulfate concentration of the solution becomes significant and since bisulfate has a greater tendency (than sulfate) to adsorb on electrode surfaces it competes for surface sites and blocks the nickel ion path to the solution. This is the reason for the secondary passivity and it is encountered only at high potentials where adsorption of the negatively charged bisulfate ions becomes favorable. This explanation accounts also for the existence of secondary passivity only in sulfuric acid solution and its absence when other ions can compete for the surface sites. At higher potentials oxygen evolution becomes possible by tunneling through some of the surface sites (exposed nickel atoms or adsorbed bisulfate). The model assumes that charge transfer is induced only by the dissolution of nickeous ions (Ni Ni2++ 2e-) or oxygen evolution. The adsorption of anions on the oxide surface changes the surface capacitor but does not contribute to the total charge transfer. The electrochemical transformations from adsorbed hydroxides to surface hydroxyls and the build up of the nickel oxides layer are assumed to contribute only a small fraction of the total current. Four different types of surface species are distinguished below: adsorbed bisulfate, surface oxide, surface hydroxide, and the exposed nickel atoms. This view of the system is very similar to the multiple adsorption mechanism p r o p e d by Jouanneue et a1.I6 Their model which assumed multiple adsorption of combinations of hydroxide, sulfate, and bisulfate accounted well for the impedance response of anodized nickel electrode. It was later replaced by a more simplistic picture containing fewer parameters.'* McDougal et al.19 describe the dissolution of nickel in the transpassive region as a result of leakage through defects in the oxide layer. The exact chemical nature of such defects is not specified, but they are likely to be initiated by localized surface species such as unstable adsorbed anions. The equations that are used in the model are equally suitable to describe the mechanism of defect formation and subsequent oxide dissolution.
-
Haim et al. NiHS0, (ads)
(e3
oxygen evolution
It
tl
F i p e 3. Schematic model for nickel dissolution.
Model Eqmtiom. The proposed model includes four variables: the electrode potential and three adsorbed species. F m 3 depicts the postulated reaction pathways of the adsorbed species. The three surface spwies (and their surface coverage) are adsorbed bisulfate (e,), NiO ( q ) , and NiOH (e - q ) (for reasons of convenience 0 is defined as the total coverage by the oxide and hydroxide). The surface fraction with accessible nickel atoms is (1 - 8 - e,). Since only three variables are required to account for chaos, bisulfate will be assumed to be in adsorption equilibrium, thus eliminating one variable from the model. Furthermore, since periodic behavior and all their singularities can be accounted for with only two dynamic variables, we later lump the surface oxide and hydroxide into one variable. Adsorbed bisulfate (e,) is formed by adsorption from the solution on the bare surface. Ni + HSO, $ NiHSO4-,* (la) Ni dissolves by a first-order, acid-catalyzed, Tafel-like reaction
-
Ni Ni2+ At high potentials oxygen evolves
H20
-
+ 2e-
(1b)
+
2H+ + Y2O2 2e-
(14
-
(Id)
Surface hydroxyl coverage is formed on the bare surface
+
Ni H 2 0 and is oxidized to NiO via
NiOH
+
+ H+ + e-
+
NiOH --c NiO H+ eNiO dissolves through acid-catalyzed dissolution
+
-
(le)
+
NiO 2H+ Ni2+ H 2 0 (If) Assuming that all electrochemicalreaction rates follow Tafel law, r exp(aFt$/RT) (see notation below), a mass balanoe over the adsorbed bisulfate gives
-
where F is Faraday constant, R is the gas constant, t$ represents the electrode potential, kafand kabdenote the coefficients of the adsorption and desorption of bisulfate (eq la), and at and ab are the transfer caefficients for the adsorption and dewrption. When no other ions are present in the solution then [H+] = [HS04-] + 0.S[S042-],and since K,[HS04-] = [H+][S042-] (where K, = 1.2 X 1W2 at 25 "C),[HSO,-] is approximated as being linearly dependent on sulfuric acid concentration or when only sulfuric acid is introduced into the cell, it is linearly dependent on the concentration of hydrogen ions. When fast equilibrium is established, then a Langmuir type isotherm describes the adsorption of bisulfate (c, denotes activity of the hydrogen ions). (3)
(17) Guckenheimer, J.; Holmes,P. Nonlinear Oscillations, Dynamical Systems and Bifurcation of Vector Fields; Springer-Vnlag; New York, 1983. (18) Keddam, M.; Takenouti, H.; Yu, N. J . Electrochem. Soc. 1985,132, 2561. (19) MacDougal, B. J . Electrochem. Soc. 1980, 127,789. MacDougal, B.; Cohen, M. Ibid. 1976, 123, 191.
The Journal of Physical Chemistry, Vol. 96, No. 6, 1992 2679
Anodic Nickel Dissolution The rate of Ni dissolution through the exposed nickel atoms follows
The other dominant reaction is oxygen evolution which is assumed to occur through surface sites that are not covered by NiOH or NiO
rc = k,u’(l - e)
(6)
Surface hydroxyl is formed at a slow rate of
where for mathematical simplicity the Tafel exponents 6 and d are assumed to be identical. The hydroxyl is oxidized to NiO following
re = ke(B - q ) exp(
$)
= ke(B - q)u4
(8)
and NiO surface sites are destroyed by
emerge at a double-zero singularity and turn around at a lower concentration. Note that every model that admits a DZ singularity will admit saddle node (hence multiplicity), Hopf, and saddle loop bifurcations. To achieve the observed features the following properties should also be obeyed: The D Z should bifurcate into stable limit cycles from the lower transpassive steady-state branch and the boundaries should show the orientation described in Figure 2a. The orientation condition is assured from the universal unfolding of the DZ if the previous two conditions are satisfied; Since at the DZ point the three bifurcations merge tangentially and since the branch emerging from the cusp point is a saddle stable node, the Hopf and saddle loop must bifurcate from the DZ point toward lower concentrations. At sufficiently small concentrations, however, any physically realistic model must predict only unique and stable steady states. In the absence of other DZ points, the Hopf and saddle loop bifurcations must turn around and continue for higher concentrations. Estimation of the parameters in the two-variable model can be based on the fitting of the observed boundaries and singularities in Figure 2. To simplify the algebra let us construct the bifurcation set for high concentrations (ch >> 1). Approximating f as
(9) Surface balances on the two oxide sites (NiOH and NiO) are then
the steady-state equation is
where CMis the surface molar capacity. The charge balance reads
Cl
d4
=I-rb-rc dt
Now, as
ch
unless u
-
(11)
where C, is the double layer capacitance per unit electrode area, I is the imposed current, and rb and r, are the chargetransfer rates by the two main reactions; the other reactions are assumed to be slow and do not contribute significantly to the charge balance. In dimensionless form, with the following variables
-
m,
+
J= auy 0 and t > v, where now 1 J= u6
(19)
(20)
C +bCh~‘
The saddle node bifurcation (dJ/du = 0) corresponding to these two asymptotes is
the balances take the form
These correspond to the lower and upper branches of the saddle node curve. The lower branch is independent of ch,in agreement with the experimental findings, while along the upper branch J varies like Ch6/(r+6). To obtain the cusp point, we have to solve simultaneously the equation To reduce the problem into a two-variable form, we note that if r2is sufficiently small, q can assumed to be in a pseudo steady state, 7 = u”e/(cch P ) and , the reduced model takes the form
+
-dJ-- -d 2-J - 0 du
du2
where J is defined in eq 18. This gives
+ a(6 + y + t ) * u * + ~ + ~ (24) - v)2u‘-Y (6 + €)2U*+C+ bCh2
t2uf
(t
Parameter Estimation. The model should exhibit the following features (see Figure 2a): (a) steady-state multiplicity and saddle node boundaries with the inclination shown in Figure 2a; and (b) oscillations bounded by Hopf and saddle loop bifurcations that
To derive this expression, note that if J =f,(u)/f2(u) then the cusp condition (eq 23) is equivalent to
2680 The Journal of Physical Chemistry, Vol. 96, No. 6,1992
where$ andf" are the first and second derivatives. The Hopf bifurcation occurs in a system of two differential equations if the trace of the Jacobian is zero. Contrary to the saddle node bifurcation case,the Hopf bifurcation depends on the capacities and time scale ratios rl and r2.The trace condition is satisfied when
Denoting C(u) E g(u)/f(u) and using the approximationfI(u) = chub reduces it to the form (1 + c y
+ r$h(ayU7+' - 6)c = o
We shall compute now the asymptotes for c h the Hopf bifurcation, eq 27 takes the form e2
- -
- (6qCh - 2)c + 1 = o
m.
(27)
If u
0 at
(28)
It follows that CI H (6rlch)-', and the lower Hopf boundary at high c h is Ut"'
At the upper branch e
-
= elC/bCh
(29)
= and u is finite at high concentrations
U " '
H
eC/bCh
(31)
and substitution in eq 27 yields
It follows that 42
- ch
J2 = ( 1 4 + ~ ~au7)-
e2
1 + c2
-
(33) Ct
(34)
The double zero point (DZ) is the point where the Hopf and SN bifurcations coalesce. Experimentally, it was found at a concentration that is 30-40 times higher than that of the cusp point. The approximation c h can be applied to the computation of this point, yielding
-
and
e=
(36)
and the corresponding c h and J can be determined now from eqs 22 and 32. Parameter estimation can be carried out by calculating the ratio of J and c h at one singular point (the cusp) to J and c h at another singular point (the DZ)and comparing these to the ratia of the current and concentration in the experimental bifurcation map (Figure 2a). A fine tuning of the parameters Can achieve a curve which fits the experimentaldata. For emhple, Figure 2b describes the bifurcation map of the full 2-v (two-variables) model with a = 0.3;
b=6X c = 0.001; rl = 0.01; y = 1; a = 1; 6 = 0.5; 6 = 0.5; c = 2; v = 1 (37)
(Note that redefining u = exp([(af - ah)F4]/Rq will not alter the equations, while a = 1 and defining 6 = b/(af - ab), E = e/(ar - ab), 'k E f/(a[- ah), etc.). The boundaries now are calculated exactly. The computed map show the boundaries and singularities orientation of the experi-
Haim et al.
vm(p)S
vm(w
vm(pk0
Figure 4. Motion in the vicinity of a saddle loop bifurcation. mental map. The computed DZ singularity occurs at c h Y 1400 and J 1.5 while the cusp point is at c h H 500 and oscillations exist for c h > 850. Experimentally, the oscillatory domain exists over a wider range of concentrations and currents. The three bifurcation diagrams in Figure 1, d, e, and f, are characteristic for the behavior at 850 < Ch < 1250,1250 < Ci, < 1400, and Ch > 1400, rwpactively. The qualitative agreement with the experimental observations (Figure 1a-c) is evident. The loci of the saddle node branches are found from the solution of the equation aJ/au = 0, where J is described in eq 18. The Hopf bifurcation line is found from eq 26. At the SL bifurcation (or the homoclinic orbit) a trajectory that arcapesfrom the suddie point will return to its original location thus formifig a loop with the saddle point us its corner. Analytical solution for this branch exists only for conservative systems and near singularities like the DZ point. The saddle loop bifurcation was computed by following the trajectory in the phase plane starting from a point on the uescaping" eigenvector (unstable manifold) of saddle point, until it reaches the closest point to the saddle point again. If u' describes the location of the trajectory (with the saddle point as the origin), then when the distance from it to the saddle point is minimal, the trajectory direction is perpendicular to u', Le., the scalar product Z-du changes from a negative to a positive value. By changing a parameter we search for a value where u' = 0 at the closest point. For that purpose we define the pseudoscalar function vm(p) as the vector product u' X du. As one can see from Figure 4, um(p) changes its sign when a saddle loop bifurcation occurs. Complex Behavior. The three-variable model was checked for its ability to exhibit chaotic solutions before its 2-v version was analyzed. One possible way of generating a 3-v model that admits complex behavior is to convert a bifurcation variable of a 2-v model that acquires a hysteresis loop bound by saddle loop and saddle node bifurcations (e.g., Figure lc,f) into a slow-changing parameter. We did not employ that or any other approach for generating chaos but rather searched for such solutions in the 3-v model. Several 3-51 models, for which complex behavior was not unraveled, were discarded. The 3-v model (eqs 13-15) with r2= 2 (other parameters are h range of defined in eq 37) admits complex behavior in the c 500-800. Several bifurcation diagrams showing the potential values of 60 upper peaks of the time trace (i.e., du/dt = 0, d2u/dt2 < 0, the initial 60 peaks were discarded) at increasing currents are shown if Figure 5. At c h = 500 the domain of complex behavior is bounded by perioddoubling transitions on both sides; within the complex domain the motion is mainly aperiodic but its structure is not coherent. At c h = 600 the complex domain is interrupted by a wide domain of period two oscillations and narrow domains of periods 4 and 5 . The hard transition from complex behavior to period two oscillations at J = 1.093 is portrayed in Figure 6. Wider frequency locking domains appear at higher amcentratiomi: domains of period 3,4, and 5 are evident at ch = 800. The complex domain is still bound by a perioddoubling cascade on one side but on the other a hard transition form aperiodic solution occurs. Concluding Remirks Modeling of oscillatory behaviors in catalytic systems has to rely on the physicochemical information on surface processes and has to account for the observed dynamic features. While there is a general agreement that observed oscillations in reaaion rate cannot reveal the underlying mechanism, the study of surface processes and the elementary steps responsible for the oscillations
The Journal of Physical Chemistry, Vol. 96, No. 6,1992 2681
Anodic Nickel Dissolution
1.55
I.5 2
r
8. b l
THEm
I.24
! 8.81
-8.12
- 1.34 1. b 5 j
ETR
,
0.32 .
.
..
8
I.
'I
1. 1 5
1. 0 7
1. 08
1. 0 9
1. 10
1. 11
1. 12
1.21
-
0.89
.
'
. . / : .
1. i 3
.
1. 14
'
J
Figure 5. Simulated bifurcation diagrams showing domains of aperiodic behavior at C, = 500 (upper diagram), 600 (middle), and 800 (lower). At each current 60 upper peaks are denoted after the initial 60 were discarded.
is a complex experimental task that has been conducted with partial success only for low-pressure CO oxidation on a single platinum crystal. In the absence of such information the modeling effort is limited to finding the simplest plausible mechanism that predicts qualitatively and quantitatively the observed behavior. The model presented here follows other steady-state studies of anodic Ni dissolution and accounts for three surface speciesadsorbed bisulfate, nickel oxide, and nickel hydroxide. Adsorption equilibrium is assumed for the adsorbed species and the autocatalytic step is the self-inhibitorybehavior of surface potential: at high potentials the rate of Ni dissolution declines due to surface blockage by the adsorbed bisulfate (eq 13,B < a). The oscillations are due to slow surface modificationsby the oxide and hydroxide. Once the dynamic variables are defined we present a systematic approach for modeling and parameter estimation of the observed motions and singularities. The observed and modeled behaviors include bistability, periodic oscillations bound by Hopf and saddle loop bifurcations, the disappearance of these transitions in the double-zero singularity, and the existence of chaotic solutions
-8.02
I.55 THETR
-8.12
CHAOTIC FLOW
Figure 6. Phase plane motions showing the hard bifurcation from chaotic (at J = 1.5, upper diagram, and J = 1.09) to periodic solutions ( J = 1.097, lower diagram) at C,, = 600.
which are bound by a period-doubling transition. The modeling and parameter estimation are based on comparison of experimental and simulated bifurcation maps, showing domains of different characteristic behaviors in the plane of applied constant current vs acid concentration. Qualitative agreement between the observed and predicted maps has been obtained. Better quantitative agreement can be achieved by fine-tuning of the parameters. Strategies for modeling or parameter estimation for predicting aperiodic solutions are still lacking. Due to the global nature of the motion in the phase plane one cannot expect a simple procedure for such tasks, but some classification of features and properties in aperiodic motions may be valuable.
Acknowledgment. This work was supported by the US-Israel Binational Science Foundation.
