Prévia do material em texto
<p>J. Chem. Phys. 131, 044105 (2009); https://doi.org/10.1063/1.3168400 131, 044105</p><p>© 2009 American Institute of Physics.</p><p>Relative abundance and structure of chaotic</p><p>behavior: The nonpolynomial Belousov–</p><p>Zhabotinsky reaction kinetics</p><p>Cite as: J. Chem. Phys. 131, 044105 (2009); https://doi.org/10.1063/1.3168400</p><p>Submitted: 08 April 2009 . Accepted: 12 June 2009 . Published Online: 22 July 2009</p><p>Joana G. Freire, Richard J. Field, and Jason A. C. Gallas</p><p>ARTICLES YOU MAY BE INTERESTED IN</p><p>Deterministic chaos in the Belousov–Zhabotinsky reaction: Experiments and simulations</p><p>Chaos: An Interdisciplinary Journal of Nonlinear Science 3, 723 (1993); https://</p><p>doi.org/10.1063/1.165933</p><p>Oscillations in chemical systems. IV. Limit cycle behavior in a model of a real chemical reaction</p><p>The Journal of Chemical Physics 60, 1877 (1974); https://doi.org/10.1063/1.1681288</p><p>“Black spots” in a surfactant-rich Belousov–Zhabotinsky reaction dispersed in a water-in-oil</p><p>microemulsion system</p><p>The Journal of Chemical Physics 122, 174706 (2005); https://doi.org/10.1063/1.1888386</p><p>https://images.scitation.org/redirect.spark?MID=176720&plid=1085727&setID=378408&channelID=0&CID=358608&banID=519827791&PID=0&textadID=0&tc=1&type=tclick&mt=1&hc=ec42be1e5157177578ec02e7b63d7170461d1664&location=</p><p>https://doi.org/10.1063/1.3168400</p><p>https://doi.org/10.1063/1.3168400</p><p>https://aip.scitation.org/author/Freire%2C+Joana+G</p><p>https://aip.scitation.org/author/Field%2C+Richard+J</p><p>https://aip.scitation.org/author/Gallas%2C+Jason+A+C</p><p>https://doi.org/10.1063/1.3168400</p><p>https://aip.scitation.org/action/showCitFormats?type=show&doi=10.1063/1.3168400</p><p>https://aip.scitation.org/doi/10.1063/1.165933</p><p>https://doi.org/10.1063/1.165933</p><p>https://doi.org/10.1063/1.165933</p><p>https://aip.scitation.org/doi/10.1063/1.1681288</p><p>https://doi.org/10.1063/1.1681288</p><p>https://aip.scitation.org/doi/10.1063/1.1888386</p><p>https://aip.scitation.org/doi/10.1063/1.1888386</p><p>https://doi.org/10.1063/1.1888386</p><p>Relative abundance and structure of chaotic behavior:</p><p>The nonpolynomial Belousov–Zhabotinsky reaction kinetics</p><p>Joana G. Freire,1,2 Richard J. Field,3 and Jason A. C. Gallas1,4,a�</p><p>1Instituto de Física, Universidade Federal do Rio Grande do Sul, 91501-970 Porto Alegre, Brazil</p><p>2Centro de Estruturas Lineares e Combinatórias, Universidade de Lisboa, 1649-003 Lisboa, Portugal</p><p>3Department of Chemistry, The University of Montana Missoula, Montana 59812, USA</p><p>4Rechnergestützte Physik der Werkstoffe IfB, ETH Hönggerberg, HIF E12, CH-8093 Zurich, Switzerland</p><p>�Received 8 April 2009; accepted 12 June 2009; published online 22 July 2009�</p><p>We report a detailed numerical investigation of the relative abundance of periodic and chaotic</p><p>oscillations in phase diagrams for the Belousov–Zhabotinsky �BZ� reaction as described by a</p><p>nonpolynomial, autonomous, three-variable model suggested by Györgyi and Field �Nature</p><p>�London� 355, 808 �1992��. The model contains 14 parameters that may be tuned to produce rich</p><p>dynamical scenarios. By computing the Lyapunov spectra, we find the structuring of periodic and</p><p>chaotic phases of the BZ reaction to display unusual global patterns, very distinct from those</p><p>recently found for gas and semiconductor lasers, for electric circuits, and for a few other familiar</p><p>nonlinear oscillators. The unusual patterns found for the BZ reaction are surprisingly robust and</p><p>independent of the parameter explored. © 2009 American Institute of Physics.</p><p>�DOI: 10.1063/1.3168400�</p><p>I. INTRODUCTION</p><p>It is well known that several chemical reactions are ca-</p><p>pable of displaying both periodic and chaotic oscillations in</p><p>the concentrations of reactive intermediate species. The</p><p>Belousov–Zhabotinsky reaction �BZR� is a paradigm of this</p><p>rich dynamical behavior, investigated in many papers and</p><p>featured in several books.1–4 However, the quantification of</p><p>the relative abundance and structural distribution of chaos</p><p>and periodicity in parameter space of this system has re-</p><p>mained poorly investigated. Rather than reflecting a lack of</p><p>interest, this spotty knowledge reflects the large computa-</p><p>tional effort required to construct phase diagrams for phe-</p><p>nomena represented by flows in phase space, i.e., by</p><p>continuous-time dynamical systems governed by sets of dif-</p><p>ferential equations. While it is easy to iterate discrete maps,</p><p>it is far harder and time consuming to integrate differential</p><p>equations. For this reason, most of the accumulated knowl-</p><p>edge about chaotic behavior in natural phenomena, phase</p><p>diagrams in particular, comes from investigations based on</p><p>discrete-time nonlinear maps.</p><p>Exploration of the parameter space of dynamical sys-</p><p>tems governed by differential equations has attracted some</p><p>attention recently, following a report5 that the phase diagram</p><p>of a loss-modulated CO2 laser, a flow, is surprisingly similar</p><p>to that of a textbook example of a discrete-time dynamical</p><p>system, the Hénon map. The immediate question is what</p><p>other sort of vector flows might produce isomorphically</p><p>similar phase diagrams and what new features they might</p><p>have. Surprising bifurcation phenomena have been reported</p><p>recently for systems across distinct disciplines and with vari-</p><p>ous motivations.5–15 For a survey, see Ref. 16.</p><p>Generically, the most obvious regularities observed so</p><p>far in phase diagrams of flows consist of sequences of self-</p><p>similar periodicity islands spread in chaotic phases,</p><p>“shrimps,”17,18 which based on studies of maps, are known to</p><p>underly period-doubling cascades. The rich bifurcation phe-</p><p>nomena found in flows so far contain novel features that</p><p>emerge organized in regular patterns not known in discrete-</p><p>time systems �maps�. This novelty is related to the distinct</p><p>manners that shrimps are “glued” together to form regular</p><p>patterns over extended regions in parameter space. However,</p><p>the basic organizational Leitmotiv is still based on shrimp</p><p>networks that accumulate systematically in one way or</p><p>another.6,12 Thus, is it possible to find macroscopic regulari-</p><p>ties of a different kind, i.e., other than glued shrimps?</p><p>A common feature of the models already investigated is</p><p>that they are mainly governed by polynomial equations of</p><p>motion. So, what can one expect from systems governed by</p><p>nonpolynomial equations of motion? We find novel and un-</p><p>expected features in a nonpolynomial system and describe</p><p>them in some detail.</p><p>Although nonpolynomial models of dynamical systems</p><p>exist abundantly in literature, their chaotic phases have been</p><p>explored mainly by plotting bifurcation diagrams along a</p><p>few specific one-parameter cuts. Moreover, phase diagrams</p><p>normally reported do not describe details of the chaotic</p><p>phases19 but, instead, focus mainly on bifurcation boundaries</p><p>between regions involving periodic oscillations of small pe-</p><p>riods. We investigate here the relative abundance and struc-</p><p>turing of the chaotic phases for a nonpolynomial model in-</p><p>volving three independent variables and 14 parameters and</p><p>describing BZR chaos.20 This model is complex in the sense</p><p>of Nazarea and Rice.21 Bifurcation diagrams for this model</p><p>display features resembling those recently found near certain</p><p>hubs in phase diagrams.12,16 Because not many hubs are</p><p>presently known, and there is no theoretical method to an-</p><p>ticipate the location of hubs, we perform a detailed numeri-a�Electronic mail: jason.gallas@gmail.com.</p><p>THE JOURNAL OF CHEMICAL PHYSICS 131, 044105 �2009�</p><p>0021-9606/2009/131�4�/044105/8/$25.00 © 2009 American Institute of Physics131, 044105-1</p><p>http://dx.doi.org/10.1063/1.3168400</p><p>http://dx.doi.org/10.1063/1.3168400</p><p>http://dx.doi.org/10.1063/1.3168400</p><p>cal investigation of the model. Although hubs were not</p><p>found, we do find a rather unusual structuring �described</p><p>below�: fountainlike global patterns consisting of alternate</p><p>eruptions of chaos and periodicities. Models of real chemical</p><p>reactions are particularly appealing in that they arise from</p><p>experiment and, therefore, dynamical behaviors predicted</p><p>from them should be amenable to experimental verification.</p><p>II. THE BZ REACTION</p><p>The classic BZR2,22,23 is the cerium-ion �Ce�IV�/Ce�III��</p><p>catalyzed</p><p>oxidation of malonic acid �CH2�COOH�2� by bro-</p><p>mate ion �BrO3</p><p>−� in aqueous sulfuric acid �H2SO4� media.24</p><p>In a well-stirred, closed reactor, damped temporal oscilla-</p><p>tions occur in the concentrations of various intermediate spe-</p><p>cies, e.g., Ce�IV�, Ce�III�, BrO2, HBrO2, HOBr, and Br−. The</p><p>concentration of bromalonic acid �BrMA� is an important</p><p>dynamic quantity that serves as a bifurcation parameter for</p><p>the onset and eventual disappearance of oscillation in a</p><p>closed reactor.</p><p>True chemical steady states may be achieved when the</p><p>BZR is run in a continuous flow, stirred tank reactor �CSTR�,</p><p>into which solutions containing the reactants are pumped</p><p>while reaction mixture overflows.2,3 This makes the critical</p><p>bifurcation species BrMA �vide infra� a dynamic, oscillatory</p><p>variable. Oscillatory or even chaotic4 stationary states often</p><p>appear.</p><p>Deterministic chaos was first observed in BZ-CSTR ex-</p><p>periments in 1977 by Schmitz et al.25 and better character-</p><p>ized by Hudson and Mankin.26 The appearance of periodic-</p><p>chaotic windows as the flow rate was monotonically</p><p>increased was soon observed by Turner et al.27 and by Vidal</p><p>et al.28 Several classic transitions from periodicity to chaos</p><p>were observed in 1983 by Roux.29 The observed chaotic</p><p>states generally run from mixed-mode systems at low flow</p><p>rates27,28 to more complex behaviors at higher flow rates.30</p><p>There is some uncertainty concerning the origin of BZ-</p><p>CSTR chaos. Low-flow-rate chaos can be well reproduced31</p><p>by a model based only on the homogeneous chemistry24 of</p><p>the BZR. However, there is considerable experimental evi-</p><p>dence that high-flow-rate BZ-CSTR chaos may be at least</p><p>strongly affected by imperfect mixing effects.32–35 The model</p><p>under consideration here is based only on BZR chemistry.</p><p>The basic chemistry of the BZR is referred to as the</p><p>Field, Körös, Noyes24 �FKN� mechanism. The FKN mecha-</p><p>nism may be reduced to a skeleton form referred to as the</p><p>Oregonator36,37 involving only the species HBrO2, Br−, and</p><p>Ce�IV�. However, while both the FKN mechanism and the</p><p>Oregonator at least qualitatively reproduce in simulations the</p><p>BZR oscillations, neither model has been found to generate</p><p>chaos. Györgyi and Field38 expanded the details of the reac-</p><p>tions of BrMA and the malonyl radical to produce a three-</p><p>variable model20 that reproduces well the periodic-chaotic</p><p>nature of the low-flow-rate chaos. This model is based upon</p><p>the nonstoichiometric chemical reactions �1�–�7� and leads to</p><p>a set of three nonpolynomial differential equations. The left</p><p>hand sides of Eqs. �1�–�7� are rate determining for the ap-</p><p>pearance of products on the right hand sides,</p><p>Br− + HBrO2 + H+→</p><p>k1</p><p>2BrCH�COOH�2, �1�</p><p>Br− + BrO3</p><p>− + 2H+→</p><p>k2</p><p>BrCH�COOH�2 + HBrO2, �2�</p><p>2HBrO2→</p><p>k3</p><p>BrCH�COOH�2, �3�</p><p>1</p><p>2HBrO2 + BrO3</p><p>− + H+→</p><p>k4</p><p>HBrO2 + Ce�IV� , �4�</p><p>HBrO2 + Ce�IV�→</p><p>k5</p><p>1</p><p>2HBrO2, �5�</p><p>BrCH�COOH�2 + Ce�IV�→</p><p>k6</p><p>Br−, �6�</p><p>Ce�IV� + CH2�COOH�2→</p><p>k7</p><p>inert products. �7�</p><p>The concentrations of the principal reactants �H+, BrO3</p><p>−,</p><p>CH2�COOH�2� in Eqs. �1�–�7� are held constant, leaving four</p><p>dynamic variables �Br−, HBrO2, Ce�IV�, and BrMA� de-</p><p>scribed by four differential equations. Thus this model ex-</p><p>plicitly takes into account the dynamics of �BrMA�. The dif-</p><p>ferential equation describing the behavior of �Br−� is</p><p>eliminated using the pseudo-steady-state approximation,39</p><p>leaving three dynamic equations in �HBrO2�, �Ce�IV��, and</p><p>�BrMA�, as well as an algebraic expression for �Br−� given</p><p>below, in Eq. �13�.</p><p>III. THE NONPOLYNOMIAL EQUATIONS</p><p>By introducing the following equivalences:</p><p>X � �HBrO2�, A � �BrO3</p><p>−� ,</p><p>Y � �Br−�, C � �Ce�III�� + �Ce�IV�� ,</p><p>Z � �Ce�IV��, H � �H+� ,</p><p>V � �BrMA�, M � �CH2�COOH�2� ,</p><p>as well as the rate constants k1−k7 specified in reactions</p><p>�1�–�7� one obtains the following set of three nonpolynomial</p><p>044105-2 Freire, Field, and Gallas J. Chem. Phys. 131, 044105 �2009�</p><p>differential equations:1,20</p><p>dx</p><p>d�</p><p>= T0�− k1HY0xỹ + k2AH2 Y0</p><p>X0</p><p>ỹ − 2k3X0x2</p><p>+ 0.5k4A0.5H1.5X0</p><p>−0.5�C − Z0z�x0.5</p><p>− 0.5k5Z0xz − kfx , �8�</p><p>dz</p><p>d�</p><p>= T0�k4A0.5H1.5X0</p><p>0.5</p><p>C</p><p>Z0</p><p>− z�x0.5 − k5X0xz</p><p>− �k6V0zv − �k7Mz − kfz , �9�</p><p>dv</p><p>d�</p><p>= T0�2k1HX0</p><p>Y0</p><p>V0</p><p>xỹ + k2AH2 Y0</p><p>V0</p><p>ỹ + k3</p><p>X0</p><p>2</p><p>V0</p><p>x2</p><p>− �k6Z0zv − kfv , �10�</p><p>where</p><p>x �</p><p>X</p><p>X0</p><p>, y �</p><p>Y</p><p>Y0</p><p>, z �</p><p>Z</p><p>Z0</p><p>, v �</p><p>V</p><p>V0</p><p>, � �</p><p>t</p><p>T0</p><p>, �11�</p><p>with the definitions</p><p>X0 =</p><p>k2</p><p>k5</p><p>AH2, �12a�</p><p>Y0 = 4X0, �12b�</p><p>Z0 =</p><p>CA</p><p>40M</p><p>, �12c�</p><p>V0 = 4</p><p>AHC</p><p>M2 , �12d�</p><p>T0 =</p><p>1</p><p>10k2AHC</p><p>. �12e�</p><p>Apart from terms with noninteger exponents, Eqs. �8� and</p><p>�10� contain an additional nonpolynomial dependence, Eq.</p><p>�13�, arising from elimination of the �Br−� equation using the</p><p>pseudosteady state approximation to yield ỹ, the scaled,</p><p>pseudo-steady-state39 value of �Br−�:</p><p>ỹ =</p><p>�k6Z0V0zv</p><p>�k1HX0x + k2AH2 + kf�Y0</p><p>. �13�</p><p>The quantities � and � are historical artifacts originally de-</p><p>fined to separate the reactions of Ce�IV� with BrMA and M;</p><p>they have no chemical significance. The �C−Z0z� terms in</p><p>Eqs. �8� and �9� are inserted intuitively to account for the</p><p>depletion of Ce�III� as Ce�IV� is produced. All parameters</p><p>appearing in the equations above are collected in Table I,</p><p>along with the basic numerical values from Ref. 20.</p><p>The most important experimental parameter in the BZ-</p><p>CSTR system is the inverse residence time, kf, �flow rate�/</p><p>�reactor volume�, which may be very precisely controlled</p><p>and to which the system is dramatically sensitive. Experi-</p><p>ments are typically monitored by measuring �Ce�IV���Z in</p><p>the CSTR. Thus the independent parameter in nearly all of</p><p>our calculations is kf, and z is sometimes used to demonstrate</p><p>the dynamic behavior of the model as kf is varied. The other</p><p>readily controlled experimental parameters are the concentra-</p><p>tions of reactants �A, C, H, and M� in the feed streams. We</p><p>assume the constant reactant concentrations in the CSTR are</p><p>identical to their concentrations in the feed streams.</p><p>Figure 1 shows bifurcation diagrams obtained by plot-</p><p>ting the local maximal values of z�t� as a function of kf. In</p><p>both diagrams, the z axis was divided into a grid of 600 and</p><p>the kf axis into 1200 equally spaced values. We integrated</p><p>Eqs. �8�–�10� using a standard fourth-order Runge–Kutta al-</p><p>gorithm with fixed time step h=2�10−6. The first 7�104</p><p>steps were discarded as transient. During the next 140</p><p>�104 steps we searched for the local maxima of z, which</p><p>were then plotted in the bifurcation diagram. This procedure</p><p>was repeated for each of the 1200 values of kf. Computations</p><p>TABLE I. Numerical values of rate constants and parameters fixed in our</p><p>simulations, taken from Ref. 20, in the same units.</p><p>k1=4�106 k2=2.0</p><p>k3=3�103 k4=55.2</p><p>k5=7�103 k6=0.09</p><p>k7=0.23 kf =3.9�10−4</p><p>A=0.1 C=8.33�10−4</p><p>H=0.26 M =0.25</p><p>�=666.7 �=0.3478</p><p>0.0003 0.00046k</p><p>0.6</p><p>2.7</p><p>z</p><p>f</p><p>(a)</p><p>0.00033 0.00037k</p><p>0.7</p><p>2.5</p><p>z</p><p>f</p><p>(b)</p><p>FIG. 1. Bifurcation diagrams illustrating the qualitative similarity of the cascading when either increasing or decreasing the bifurcation parameter. �a� Period-1</p><p>solutions exist on both extremes of the diagram. The box marks the region magnified on the right. �b� Period-3 solutions exists on both ends of the diagram.</p><p>044105-3 Quantifying chaos in the BZ reaction J. Chem. Phys. 131, 044105 �2009�</p><p>were started at the minimum value of kf from the initial</p><p>conditions x0=z0=v0=0.5 and continued by following the at-</p><p>tractor, namely, by using the values of x, z, and v obtained at</p><p>the end of one calculation at a given kf to start a new calcu-</p><p>lation after incrementing kf. This is a standard way of gen-</p><p>erating bifurcation diagrams, and the rationale behind it is</p><p>that generically, basins of attraction do not change signifi-</p><p>cantly upon small changes in parameters, thereby ensuring a</p><p>smooth unfolding of the bifurcation curve.</p><p>Comparing the bifurcation diagrams in Fig. 1 above with</p><p>Fig. 1�a� of Györgyi and Field20 one finds that the diagrams</p><p>are virtually identical, meaning that using Poincaré sections,</p><p>as done in Ref. 20, or using the maxima of the variable</p><p>produces identical sequences of bifurcations.</p><p>The procedure described above is used to compute all</p><p>phase diagrams presented</p><p>in Sec. IV, except that instead of</p><p>considering local maxima of z, all points after the transient</p><p>were used to calculate the Lyapunov spectrum, i.e.,</p><p>Lyapunov exponents for the three variables of the model.</p><p>The construction of Lyapunov phase diagrams is a demand-</p><p>ing computational task. For instance, to classify a mesh of</p><p>N�N parameter points demands the computation of N2 ba-</p><p>sins of attraction, needed to sort out all possible solutions for</p><p>each parameter set. The computation of each individual basin</p><p>of attraction requires investigating sets of initial conditions</p><p>over a M �M mesh in phase space. The quality of the final</p><p>diagrams depends sensibly on both N and M being as large</p><p>as possible. Typically, we used 400�400 grids, although</p><p>600�600 grids were also used.</p><p>IV. PHASE DIAGRAMS</p><p>In this section we present several two-parameter phase</p><p>diagrams discriminating the dynamic behavior �chaotic or</p><p>periodic� as kf and one of the reactant concentrations is var-</p><p>ied. All parameters are varied over experimentally accessible</p><p>ranges. Rate constant values are not readily variable experi-</p><p>mentally. However, we do present phase diagrams in which</p><p>the values of k2 and k4 are varied because these important</p><p>quantities control important properties of the oscillations and</p><p>chaos.</p><p>Figure 2�a� shows a bifurcation diagram similar to those</p><p>in Figs. 1�a� and 1�b�, but now considering the total cerium-</p><p>ion concentration, C, as the bifurcation parameter while</p><p>keeping A=0.1 constant. However what would happen for</p><p>other values of A close to this one? Figure 2�b� shows how</p><p>chaos and periodicity are distributed in the C�A parameter</p><p>section. It shows that for relatively large variations of A</p><p>around 0.1 the bifurcation diagram remains essentially the</p><p>same, apart from an overall stretching.</p><p>Figure 2�b� illustrates a typical Lyapunov phase diagram</p><p>obtained by computing the three Lyapunov exponents for</p><p>Eqs. �8�–�10�, ordering them such that �3��2��1, and plot-</p><p>ting the largest nonzero exponent. We plotted �2 whenever</p><p>��1��10. As is well known,40 Lyapunov exponents provide a</p><p>handy quantity allowing one to discriminate between chaos</p><p>�positive exponents� and periodic oscillation �negative expo-</p><p>nents�. Figure 2�b� contains a scale of colors that, after suit-</p><p>able renormalization to reflect the extrema of exponents, was</p><p>used to construct similar figures below where scales are</p><p>omitted. Noteworthy is the great spread of the magnitude of</p><p>the exponents, an indication of the stiffness41 of the model.</p><p>Figure 3 displays phase diagrams for representative pa-</p><p>rameter sections: kf �A, kf �H, and kf �M. The top row</p><p>displays large views of the parameter space, discriminating</p><p>chaos from periodicity. Two of the diagrams contain the let-</p><p>ter � to mark, as is done in Fig. 2�b�, domains where one</p><p>finds chaos overwhelmingly, not periodic oscillations as the</p><p>shadings seem to indicate. To characterize chaos properly,</p><p>one would need to integrate for much longer times. The</p><p>middle row shows magnifications of the boxes in the upper</p><p>row. The thin white spines represent the location of “super-</p><p>stable loci,”12 i.e., loci of local minima of negative Lyapunov</p><p>exponents. Noteworthy in the middle row are the two right-</p><p>most panels displaying fountainlike “eruptions of chaos,”</p><p>a sort of chaos geyser in specific parameter regions in the</p><p>0.0007 0.0017C</p><p>0.5</p><p>3.4</p><p>z</p><p>(a)</p><p>0.0007 0.0017C</p><p>0.05</p><p>0.125</p><p>A</p><p>-905.7 66.90</p><p>0.0007 0.0017C</p><p>0.05</p><p>0.125</p><p>A</p><p>χ</p><p>(b)</p><p>FIG. 2. �a� Bifurcation diagram along A=0.1, indicated by the horizontal line in the right panel. �b� Lyapunov phase diagram discriminating periodicity and</p><p>chaos in the C�A control space. Colors denote chaos �i.e., positive exponents� while the darker shadings mark periodicity. Chaos prevails in the accumulation</p><p>region around the letter � despite the coloration. The color scale is linear on both sides of zero but not uniform. Note the high spread of exponents, indicating</p><p>stiffness. The bifurcation diagram has a resolution of 1200�600 pixels, while the phase diagram displays 600�600 Lyapunov exponents.</p><p>044105-4 Freire, Field, and Gallas J. Chem. Phys. 131, 044105 �2009�</p><p>kf �M and kf �H planes. Eruptions are also present in other</p><p>parameter sections visible in Fig. 2�b�, in the figures below,</p><p>and in several additional cuts of the parameter space, e.g.,</p><p>M �A, M �H, M �C, and C�H �not shown here�.</p><p>The bottom row in Fig. 3 illustrates an inherent difficulty</p><p>of computing Lyapunov exponents: the accurate determina-</p><p>tion of zero exponents. Away from zero, approximation er-</p><p>rors intrinsic to the integration algorithm as well as compu-</p><p>tational errors occurring in the calculation of exponents</p><p>remain confined to less significant digits. However, the tran-</p><p>sition between periodic and chaotic behaviors is marked by</p><p>zero exponents, where error in the less significant digits is</p><p>precisely what survives. The bottom row in Fig. 3 is a replot</p><p>of the panels in the middle row, but using a sharp cut at zero,</p><p>not the tolerant cut that allows accommodation of numerical</p><p>uncertainties. In the middle row we plotted �2, the second</p><p>largest exponent, whenever �1, the largest exponent, could</p><p>not be plotted after being “declared” zero, i.e., when it</p><p>obeyed ��1��10 �the “tolerance limit”�. In the bottom row,</p><p>we plotted �2 whenever �1�0 �strict limit�. The granularity</p><p>exposes regions of nearly zero exponents and higher numeri-</p><p>cal inaccuracies. Even under these very strict plotting re-</p><p>quirements, the overall structure of the phase diagrams re-</p><p>mains discernible. Recall that Eqs. �8�–�10� are stiff, thus</p><p>exacerbating numerical difficulties near the line of zero ex-</p><p>ponents.</p><p>Figure 4 shows the relative abundance and organization</p><p>of chaos and periodicity in the parameter section kf �C, the</p><p>total cerium-ion concentration. As is clear from the panels in</p><p>the bottom row, macroscopically, the parameter cuts present</p><p>two main regimes: a region where chaos and periodicity al-</p><p>ternate regularly while experiencing a relatively moderate</p><p>compression toward each other �exemplified in the right</p><p>panel�, and regions where one finds domains to bend over</p><p>and strongly accumulate toward well defined limit curves, as</p><p>seen in the left panel. This structure is surprisingly indepen-</p><p>dent of the parameter explored and markedly distinct from</p><p>structures recently found for gas and semiconductor lasers,</p><p>for electric circuits, and for a few other familiar nonlinear</p><p>oscillators.10,11,15,16 In particular, the alternation of chaos and</p><p>periodicity here is rather different from that found around</p><p>focal hubs of periodicity.12,16</p><p>Figure 5 presents progressively more highly resolved</p><p>phase diagrams obtained when varying simultaneously the</p><p>rate constants k2 and k4. The overall organization is similar to</p><p>that seen above. However in this parameter plane the com-</p><p>pression is not as strong as in previous parameter sections,</p><p>allowing finer details to be recognized even at moderate</p><p>resolutions. These phase diagrams allow one to understand</p><p>how the regimes described in Fig. 4 interconnect with each</p><p>other: by conspicuous “necks,” clearly seen in the panels in</p><p>the upper row, that get less and less pronounced when com-</p><p>ing closer and closer to the limit curve, the thick accumula-</p><p>tion line dominated by chaos. In the center panel of the</p><p>middle row one sees period-adding cascades and accumula-</p><p>tions similar to those observed in lasers, electric circuits, and</p><p>8e-05 0.0006k</p><p>0.02</p><p>0.2</p><p>A</p><p>f 0.0003 0.0009k</p><p>0.2</p><p>0.6</p><p>H</p><p>f</p><p>χ</p><p>0.0003 0.0009k</p><p>0.1</p><p>0.3</p><p>M</p><p>f</p><p>χ</p><p>0.0003 0.00045k</p><p>0.07</p><p>0.14</p><p>A</p><p>f</p><p>-1081.8 60.10</p><p>0.0003 0.00045</p><p>0.07</p><p>0.14</p><p>f 0.0003 0.00045k</p><p>0.24</p><p>0.3</p><p>H</p><p>f</p><p>-1044.1 65.20</p><p>0.0003 0.00045</p><p>0.24</p><p>0.3</p><p>f 0.0003 0.00045k</p><p>0.2</p><p>0.27</p><p>M</p><p>f</p><p>-1015.1 68.50</p><p>0.0003 0.00045</p><p>0.2</p><p>0.27</p><p>f</p><p>0.0003 0.00045k</p><p>0.07</p><p>0.14</p><p>A</p><p>f</p><p>-1081.8 60.10</p><p>0.0003 0.00045</p><p>0.07</p><p>0.14</p><p>f 0.0003 0.00045k</p><p>0.24</p><p>0.3</p><p>H</p><p>f</p><p>-1044.1 65.20</p><p>0.0003 0.00045</p><p>0.24</p><p>0.3</p><p>f 0.0003 0.00045k</p><p>0.2</p><p>0.27</p><p>M</p><p>f</p><p>-1015.1 68.50</p><p>0.0003 0.00045</p><p>0.2</p><p>0.27</p><p>f</p><p>FIG. 3. Top row: phase diagrams illus-</p><p>trating structuring for three distinct pa-</p><p>rameter cuts</p><p>over wide parameter</p><p>ranges. Chaos prevails around the let-</p><p>ter �. Middle row: magnifications of</p><p>the boxes in the corresponding panel</p><p>above. Note the great similarity of the</p><p>two rightmost diagrams with that in</p><p>Fig. 2�b�. Bottom row: replot of the</p><p>panels in the middle row, but using a</p><p>strict cut for zero exponents, as indi-</p><p>cated by the color scales. To enhance</p><p>contrast, the intensity of red was</p><p>slightly increased. The granularity ex-</p><p>poses intrinsic difficulties of calculat-</p><p>ing exponents that are close to zero</p><p>�see text�.</p><p>044105-5 Quantifying chaos in the BZ reaction J. Chem. Phys. 131, 044105 �2009�</p><p>other nonlinear models.15 The bottom row in Fig. 5 illustrates</p><p>that at higher resolutions, the chaotic domains of the non-</p><p>polynomial model display the same recurrent accumulations</p><p>of shrimps found abundantly in several other models, au-</p><p>tonomous or not.6,14–18 Shrimps are also found in several</p><p>additional cuts of the parameter space, e.g., kf �M, kf �H,</p><p>and kf �C �not shown here�. It would be interesting to study</p><p>how periodicity evolves along such accumulations and to</p><p>quantify their metric properties.</p><p>V. CONCLUSIONS AND OUTLOOK</p><p>This paper reports a detailed numerical investigation of</p><p>the relative abundance of periodic and chaotic oscillations</p><p>for a three-variable, 14-parameter, nonpolynomial, autono-</p><p>mous model of the BZR. Although chaotic solutions for BZR</p><p>models have been known for many years,1–4 they were con-</p><p>fined to isolated parameter points or to specific one-</p><p>dimensional bifurcation diagrams. No global classification of</p><p>chaotic phases was attempted. The present work reports</p><p>phase diagrams discriminating chaos and periodicity along</p><p>several sections of parameter space. We also describe details</p><p>of the intertwined structuring of chaos and periodicity over</p><p>extended experimentally accessible parameter domains.</p><p>Globally, on a macroscopic scale, we have shown the</p><p>structuring of the BZR phase diagrams to display a recurring</p><p>unusual fountainlike global pattern consisting of eruptions of</p><p>chaos and periodicities. Such a pattern is very distinct from</p><p>any structuring based on shrimps and hubs reported recently</p><p>in literature for prototypical nonlinear oscillators such as gas</p><p>lasers, semiconductor lasers, electric circuits, or a low-order</p><p>atmospheric circulation model. The unusual patterns found in</p><p>the BZR are surprisingly robust and independent of the pa-</p><p>rameter explored.</p><p>Locally, on a microscopic scale, the BZR phase dia-</p><p>grams contain very peculiar hierarchies of parameter net-</p><p>works ending in distinctive and rich accumulation bound-</p><p>aries and structuring similar to that found recently for</p><p>semiconductor lasers with optical injection, CO2 lasers with</p><p>feedback, autonomous electric circuits, the Rössler oscillator,</p><p>Lorenz-84 low-order atmospheric circulation model, and</p><p>other systems. Since these accumulations of microscopic de-</p><p>tails were found in all these rather distinct physical models, it</p><p>seems plausible to expect them to be generic features of</p><p>flows of codimension two and higher.</p><p>In sharp contrast with discrete dynamical systems where</p><p>periodicities vary always in discrete unitary steps,42 chemical</p><p>reactions provide a natural experimental framework to study</p><p>how periodicities defined by continuous real numbers evolve</p><p>and get organized in phase diagrams when several param-</p><p>eters are tuned. Additionally, another enticing open question</p><p>is whether or not it is possible to infer the presence in phase</p><p>space of unstable and complicated mathematical phenomena,</p><p>e.g., homoclinic orbits, based on observation of regularities</p><p>computed/measured solely in parameter space. Are there</p><p>clear parameter space signatures of homoclinic orbits? Is it</p><p>possible to recognize in Lyapunov phase diagrams the loca-</p><p>3.0 8.0k</p><p>7.0</p><p>30.0</p><p>C</p><p>f</p><p>χ</p><p>3.0 7.0k</p><p>7.5</p><p>11.5</p><p>C</p><p>f</p><p>3.0 4.5k</p><p>7.8</p><p>9.0</p><p>C</p><p>f 3.5 5.5k</p><p>8.3</p><p>10.0</p><p>C</p><p>f</p><p>FIG. 4. Successive magnifications il-</p><p>lustrating the alternation of periodic</p><p>and chaotic solutions in the kf �C</p><p>space. Note the structural resemblance</p><p>to parameter cuts shown in Figs. 2�b�</p><p>and 3. Chaos prevails around the letter</p><p>�. Both axis were multiplied by 104 to</p><p>avoid unnecessarily long sequences of</p><p>zeros. Thus, the minimum values of kf</p><p>and C in the upper leftmost panel are</p><p>3�10−4 and 7�10−4, respectively.</p><p>044105-6 Freire, Field, and Gallas J. Chem. Phys. 131, 044105 �2009�</p><p>tion of parameter loci characterizing simple and multiple</p><p>�degenerate� homoclinic and heteroclinic phenomena?</p><p>Our choice of parameters is motivated and centered</p><p>around the set originally chosen by Györgyi and Field.20 The</p><p>phase diagrams presented here extend considerably the re-</p><p>gimes originally investigated. However, even if using vari-</p><p>ables reduced by scaling, not the original physical quantities</p><p>as done here, the effective volume of the parameter space</p><p>that needs to be explored is huge. Consequently, despite the</p><p>work reported here, the parameter space still remains mostly</p><p>unchartered and wanting much more investigation. A sure</p><p>bet, however, is that journeys through this vast space are</p><p>bound to reveal interesting dynamics and rich operational</p><p>points for sustaining individual chemical oscillators and os-</p><p>cillator networks. In particular, they might reveal the mecha-</p><p>nism responsible for the novel fountains of chaos which we</p><p>observed here so frequently and in so many phase diagrams.</p><p>ACKNOWLEDGMENTS</p><p>J.G.F. thanks Fundação para a Ciência e Tecnologia, Por-</p><p>tugal, for a Posdoctoral Fellowship, and the Instituto de</p><p>Física da Universidade Federal do Rio Grande do Sul, Porto</p><p>Alegre, for hospitality. R.J.F. thanks the Department of</p><p>Chemistry, The University of Montana, for their support.</p><p>J.A.C.G. thanks Hans J. Herrmann for a fruitful month spent</p><p>in his Institute in Zürich. He is supported by the Air Force</p><p>Office of Scientific Research under Grant No. FA9550-07-1-</p><p>0102, and by CNPq, Brazil.</p><p>1 R. J. Field and L. Györgyi, Chaos in Chemistry and Biochemistry �World</p><p>Scientific, Singapore, 1993�.</p><p>2 I. R. Epstein and J. A. Pojman, An Introduction to Nonlinear Chemical</p><p>Dynamics: Oscillations, Waves, Patterns, and Chaos �Oxford University</p><p>Press, New York, 1998�.</p><p>3 P. Gray and S. K. Scott, Chemical Oscillations and Instabilities �Oxford</p><p>Science, New York, 1990�.</p><p>4 S. K. Scott, Chemical Chaos �Oxford Science, New York, 1991�.</p><p>5 C. Bonatto, J. C. Garreau, and J. A. C. Gallas, Phys. Rev. Lett. 95,</p><p>0.9 2.5k</p><p>50</p><p>80</p><p>k</p><p>2</p><p>4</p><p>1 2k</p><p>54.8</p><p>60</p><p>k</p><p>2</p><p>4</p><p>1.5 1.8k</p><p>55.3</p><p>55.9</p><p>k</p><p>2</p><p>4</p><p>1.2 2.5k</p><p>50</p><p>70</p><p>k</p><p>2</p><p>4</p><p>1.5 1.8k</p><p>55.2</p><p>56.5</p><p>k</p><p>2</p><p>4</p><p>1.62 1.635k</p><p>55.45</p><p>55.6</p><p>k</p><p>2</p><p>4</p><p>1 2.5k</p><p>50</p><p>60</p><p>k</p><p>2</p><p>4</p><p>1.3 1.8k</p><p>55</p><p>55.5</p><p>k</p><p>2</p><p>4</p><p>1.64 1.672k</p><p>55.45</p><p>55.65</p><p>k</p><p>2</p><p>4</p><p>FIG. 5. Progressively more highly resolved phase diagrams illustrating the fine structure of the chaotic phases for selected regions of the k2�k4 parameter</p><p>space. In the center panel of the middle row one sees period-adding cascades and accumulations similar to those observed in lasers, electric circuits and other</p><p>models �Ref. 15�. The bottom row shows that under high resolution the chaotic phase of the BZR is also riddled with shrimps �Refs. 17 and 18�.</p><p>044105-7 Quantifying chaos in the BZ reaction J. Chem. Phys. 131, 044105 �2009�</p><p>http://dx.doi.org/10.1103/PhysRevLett.95.143905</p><p>143905 �2005�.</p><p>6 C. Bonatto and J.A.C. Gallas, Phys. Rev. E 75, 055204�R� �2007�.</p><p>7 Y. Zou, M. Thiel, M. C. Romano, J. Kurths, and Q. Bi, Int. J. Bifurcation</p><p>Chaos Appl. Sci. Eng. 16, 3567 �2006�.</p><p>8 V. Castro, M. Monti, W. B. Pardo, J. A. Walkenstein, and E. Rosa, Int. J.</p><p>Bifurcation Chaos Appl. Sci. Eng. 17, 965 �2007�.</p><p>9 H. A. Albuquerque, R. M. Rubinger, and P. C. Rech, Phys. Lett. A 372,</p><p>4793 �2008�.</p><p>10 G. M. Ramírez-Ávila and J. A. C. Gallas, Revista Boliviana de Fisica 14,</p><p>1 �2008�.</p><p>11 G.M. Ramírez-Ávila and J.A.C. Gallas, “Quantification of regular and</p><p>chaotic oscillations of Chua’s circuit,” IEEE Trans. Circuits-I �submit-</p><p>ted�.</p><p>12 C. Bonatto and J. A. C. Gallas, Phys. Rev. Lett. 101, 054101 �2008�.</p><p>13 J. G. Freire, C. Bonatto, C. DaCamara, and J. A. C. Gallas, Chaos 18,</p><p>033121 �2008�.</p><p>14 C. Bonatto,</p><p>J.A.C. Gallas, and Y. Ueda, Phys. Rev. E 77, 026217 �2008�.</p><p>15 C. Bonatto and J.A.C. Gallas, Philos. Trans. R. Soc. London, Ser. A 366,</p><p>505 �2008�.</p><p>16 J. A. C. Gallas, “Infinite networks of periodicity hubs in phase diagrams</p><p>of simple autonomous flows,” Int. J. Bifurcation Chaos Appl. Sci. Eng.</p><p>�to be published�.</p><p>17 J. A. C. Gallas, Phys. Rev. Lett. 70, 2714 �1993�; Physica A 202, 196</p><p>�1994�; Appl. Phys. B: Lasers Opt. 60, S203 �1995�; B. R. Hunt, J. A. C.</p><p>Gallas, C. Grebogi, J. A. Yorke, and H. Koçak, Physica D 129, 35</p><p>�1999�.</p><p>18 E. N. Lorenz, Physica D 237, 1689 �2008�.</p><p>19 A phase diagram including the inner structuring of chaotic phases for a</p><p>nonautonomous chemical flow was reported in Ref. 15.</p><p>20 L. Györgyi and R. J. Field, Nature �London� 355, 808 �1992�.</p><p>21 A. P. Nazarea and S. A. Rice, Proc. Natl. Acad. Sci. U.S.A. 68, 2502</p><p>�1971�.</p><p>22 B. P. Belousov, in Oscillations and Traveling Waves in Chemical Sys-</p><p>tems, edited by R. J. Field and M. Burger �Wiley, New York, 1985�.</p><p>23 A. M. Zhabotinsky, Chaos 1, 379 �1991�; Scholarpedia J. 2, 1435</p><p>�2007�.</p><p>24 R. J. Field, E. Kőrös, and R. M. Noyes, J. Am. Chem. Soc. 94, 8649</p><p>�1972�.</p><p>25 R. A. Schmitz, K. R. Graziani, and J. L. Hudson, J. Chem. Phys. 67,</p><p>3040 �1977�.</p><p>26 J. L. Hudson and J. C. Mankin, J. Chem. Phys. 74, 6171 �1981�.</p><p>27 J. S. Turner, J.-C. Roux, W. D. McCormick, and H. L. Swinney, Phys.</p><p>Lett. A 85, 9 �1981�.</p><p>28 C. Vidal, S. Bachelart, and A. Rossi, J. Phys. �Paris� 43, 7 �1982�.</p><p>29 J.-C. Roux, Physica D 7, 57 �1983�.</p><p>30 L. Györgyi, R. J. Field, Z. Noszticzius, W. D. McCormick, and H. L.</p><p>Swinney, J. Phys. Chem. 96, 1228 �1992�.</p><p>31 L. Györgyi, S. Rempe, and R. J. Field, J. Phys. Chem. 95, 3159 �1991�.</p><p>32 I. R. Epstein, Nature �London� 346, 16 �1990�.</p><p>33 I. R. Epstein, Nature �London� 374, 321 �1995�.</p><p>34 F. W. Schneider and A. F. Münster, J. Phys. Chem. 95, 2130 �1991�.</p><p>35 L. Györgyi and R. J. Field, J. Phys. Chem. 92, 7079 �1988�.</p><p>36 R. J. Field and R. M. Noyes, J. Chem. Phys. 60, 1877 �1974�.</p><p>37 R. J. Field, Scholarpedia J. 2, 1386 �2007�.</p><p>38 L. Györgyi and R. J. Field, J. Phys. Chem. 95, 6594 �1991�.</p><p>39 R. J. Field, Scholarpedia J. 3, 4051 �2008�.</p><p>40 T. Tél and M. Gruiz, Chaotic Dynamics: An Introduction Based on Clas-</p><p>sical Mechanics �Cambridge University Press, Cambridge, 2006�.</p><p>41 J. H. E. Cartwright, Phys. Lett. A 264, 298 �1999�.</p><p>42 J. A. C. Gallas, Physica A 283, 17 �2000�.</p><p>044105-8 Freire, Field, and Gallas J. Chem. Phys. 131, 044105 �2009�</p><p>http://dx.doi.org/10.1103/PhysRevE.75.055204</p><p>http://dx.doi.org/10.1142/S0218127406016987</p><p>http://dx.doi.org/10.1142/S0218127406016987</p><p>http://dx.doi.org/10.1142/S0218127407017689</p><p>http://dx.doi.org/10.1142/S0218127407017689</p><p>http://dx.doi.org/10.1016/j.physleta.2008.05.036</p><p>http://dx.doi.org/10.1103/PhysRevLett.101.054101</p><p>http://dx.doi.org/10.1063/1.2953589</p><p>http://dx.doi.org/10.1103/PhysRevE.77.026217</p><p>http://dx.doi.org/10.1098/rsta.2007.2107</p><p>http://dx.doi.org/10.1103/PhysRevLett.70.2714</p><p>http://dx.doi.org/10.1016/0378-4371(94)90174-0</p><p>http://dx.doi.org/10.1016/S0167-2789(98)00201-2</p><p>http://dx.doi.org/10.1016/j.physd.2007.11.014</p><p>http://dx.doi.org/10.1038/355808a0</p><p>http://dx.doi.org/10.1073/pnas.68.10.2502</p><p>http://dx.doi.org/10.1063/1.165848</p><p>http://dx.doi.org/10.1021/ja00780a001</p><p>http://dx.doi.org/10.1063/1.435267</p><p>http://dx.doi.org/10.1063/1.441007</p><p>http://dx.doi.org/10.1016/0375-9601(81)90625-3</p><p>http://dx.doi.org/10.1016/0375-9601(81)90625-3</p><p>http://dx.doi.org/10.1051/jphys:019820043010700</p><p>http://dx.doi.org/10.1016/0167-2789(83)90115-X</p><p>http://dx.doi.org/10.1021/j100182a038</p><p>http://dx.doi.org/10.1021/j100161a038</p><p>http://dx.doi.org/10.1038/346016a0</p><p>http://dx.doi.org/10.1038/374321a0</p><p>http://dx.doi.org/10.1021/j100159a012</p><p>http://dx.doi.org/10.1021/j100336a011</p><p>http://dx.doi.org/10.1063/1.1681288</p><p>http://dx.doi.org/10.1021/j100170a041</p><p>http://dx.doi.org/10.1016/S0375-9601(99)00793-8</p><p>http://dx.doi.org/10.1016/S0378-4371(00)00123-0</p>