Small Parasitic Parameters and Chemical Oscillations BY €3. F. GRAY AND L. J. AARONS School of Chemistry University of Leeds Leeds LS2 9JT Received 23rd July 1974 In treating chemical oscillators it is common practice to use no more than two concentrations as dependent variables due to the mathematical difficulties with more than two. Other variables necessarily involved to make the system nontrivial-the "pool "chemicals-are treated as parameters which are independent of time. We investigate here under what conditions the "pool " chemicals can be treated as constants without qualitatively effecting the behaviour of the system. Mathematical methods developed by Tikhonov are used to study the effect of a small parameter on the roots of the characteristic equation.This small parameter may be chosen as the ratio of the initial concentrations of the reactants to the "pool " chemicals or for dilute solutions the reciprocal of the heat capacity. A number of interesting results are obtained in which the slow variation of the " pool " chemicals can either produce a limit cycle where there was none previously or place severe restrictions on the rate constants so as to exclude regions where interesting instabilities have been found in the two variable case. Multistability in the two variable system lends itself very well to the production of limit cycles of the relaxation type. Finally it is shown that it is possible to devise thermokinetic oscillators with very small temperature oscillations provided the energy equation is highly nonlinear.1. INTRODUCTION It has become common practise in designing models of chemical oscillators to limit the description to two dependent variables. The reasons for this are obvious being due to the formidable mathematical apparatus that exists for treating non-linear differential equations with only two variables. Several two dimensional kinetic schemes that exhibit undamped oscillations have been developed including the so-called "Brusselator ''.2* However the latter model employs a physically unrealistic termolecular step. In fact Hanusse has shown that it is impossible for any two variable kinetic scheme which has only unimolecular and bimolecular steps to exhibit undamped oscillations. The one exception to this the Lotka-Volterra ~cherne,~ is conservative and so is again physically unacceptable.It is implied in all these models that certain substances generally referred to as " pool " chemicals (such as fuel and end products) are held constant either by being in large excess or by flowing them in (out) as soon as they are used (produced). It is the purpose of this paper to investigate under what conditions the '' pool " chemicals can be treated as constant and what effect their variation may have. We will be concerned with a system of equations typically of the form where p is a small positive parameter. The solutions of these equations fall into two regions the region of slow motion and the region of rapid motion. In the first region a representative point moves comparatively slowly (ij bounded) within a and small neighbourhood (of order p) of the curve F(x;y) = 0.Outside this neighbour- hood the representative point moves in rapid jumps and in the limit as p+o the equations of rapid motion can be written 1 y = yo = constant X-= -F(x; yo). Lc s 9-5 129 CHEMICAL OSCILLATORS Furthermore for a point to stay within a small neighbourhood of F(x ;y) = 0 then this region must be one of stable equilibrium for the equations of rapid motion (2). This will be true if all s roots of the characteristic eqn (3) have negative real parts. The possibility of discontinuous oscillations then occurs if a representative point alternates between regions of slow and rapid motion. The points of transition between the regions of slow and rapid motion-the “jump ” points-are given by the inter- section of the curves F(x;y) = 0 and D(x;y) = 0 where D(x;y) is the Jacobian The remainder of this paper will be given over to four examples.First we consider the Lotka-Volterra scheme and show that if the fuel is allowed to vary the system can show non-conservative rather than conservative oscillations. In section 3 it is shown that when the “ pool ” chemicals are allowed to vary in the “ Brussel-ator” severe restrictions are placed on the rate constants so as to exclude certain regions where it was thought to oscillate. In section 4a kinetic scheme is devised in which slow variation of one of the “ pool ” chemicals gives rise to discontinuous oscillations of the type described above.Finally in section 5 a scheme originally proposed by Edelstein is modified and discontinuous oscillations are demonstrated when the temperature is allowed to vary slowly. 2. THE LOTKA-VOLTERRA SCHEME The set of differential equations originally proposed by Lotka in an ecological context can be reset in a chemical context via the following kinetic scheme kI A+X+2X k2 X+Y+2Y k3 A+Y+B. The rate equations for this scheme are dA -= -klAX-li,AY+k(A,-A) dt dX -= k,AX-k2XY dt where species A diffuses into the system from the outside where its value is assumed B. F. GRAY AND L. J. AARONS constant at Ao. It is convenient to change to new dimensionless variables A' = A/Ao X' = X/Xo Y' = Y/ Yowhere Ao Xo Yoare the initial values of A,X Y respectively and let Xo/Ao-Yo/Ao= 5 < 1.Eqn (6) become dA' __ = -k,XoA'X'-k,YoA'Y'+k(l -A') dt ; It can be seen that eqn (7) only reduce to the usual Lotka-Volterra scheme as 5+O provided k2Xo is of order 1/< while kl Xo and k3Yoare of order 1. Under these conditions eqn (7) possess a singularity at the point X' = Y' = A' = 3[-k/2+ (2k+k2/4)'] which is easily shown to be an unstable focus.8 These equations were programmed on a Solatron HS7 analogue computer and for a wide range of initial conditions the trajectories always wound onto a stable limit cycle. A typical run is shown in fig. 1. Thus by including the slow variation of the "pool " chemical A the system is made structurally stable. In the next section we will give an example of the opposite effect where the inclusion of a small parameter changes an oscillatory scheme to a nonoscillatory one.X' 0 1 2 Y' FIG.1.-Lotka-Volterra scheme. 3. THE BRUSSELATOR The following kinetic scheme was devised and has been exhaustively studied by the Brussels school 2* kl A+X k2 B+X+Y+D k3 2X+Y-+3X k4 X-+E CHEMICAL OSCILLATORS where in the original model A B D,E were assumed to be constant in time. Tyson later removed the necessity of the termolecular step by introducing another variable. The rate equations for the two dimensional model are dX -= k1A-k,BX+ksX2Y-kqX dt dY (9) -= k2BX-k,X2Y. dt The equations have a singularity for X = k,A/k4 Y = k2k4B/klk3A. It is easily deduced that the necessary condition for this point to be unstable is Lefever and Nicolis lo have shown that under certain conditions this stationary point is surrounded by a stable limit cycle.If A and B are allowed to vary subject to external fluxes we get the two additional equations Changing to dimensionless variables X' Y' A' B' as before letting Xo -Yo< Ao~Bo such that Xo/Ao = < and introducing a new time scale z = (l/c)t we get the four equations dA' -= <(-kl A' + J,/Ao) dz dY' -= k,XoB'X '-k,SXoYoX"Y'. dz It can be seen that eqn (12) only reduce to eqn (9) in the limit <+O provided k3 and k4are of order (115) and k and k are of order 1. Alternatively the flux terms have to become unrealistically large. Therefore under ordinary conditions this scheme will not oscillate in the region of parameters chosen by Lefever and Nicolis.lo This result shows clearly that one must be careful when maintaining various reactants at constant concentration.4. A CHEMICAL DISCONTINUOUS OSCILLATOR Consider the following kinetic scheme kl A+X-+2X k2 x+x-+x k3 X + B+products k4 B4products. B. F. GRAY AND L. J. AARONS The rate equations are dX --klAX-k2X2-ksXB dt dB = -k,XB-k,B+$B dt dA -7 -k,AX+$, dt where 4Aand 4Bare flow rates. Changing to dimensionless variables X' A' B' and letting Xo/Ao-Bo/Ao = 4 4 1 we obtain dA' -= -k,X,A'X'+&, dt where 4L = 4*/A0 and & = &/AO. If k2X0,k4,4Aand 4Bare all of order 1/5 and k,Xo is of order unity eqn (15) become dA' --A'X'+&k.dt They are then in the form of eqn (1). The equation of slow motion is given by F(x;y) = 0 which yields X' = 0,B' = 46 (17) and X'2+(1 -Af)Y+(4;-A') = 0 B' = $E;/(XI+f). (18) These solutions are sketched in fig. 2 in the X' A' plane. The branch X' = 0 is stable for A between 0 and P (A' = 4;) at which point it becomes unstable. The upper branch is stable down to the "jump " point R (X' = -1+ ,/#; ; A' = -1+ d&) and the middle branch between P and R is unstable. Thus for 41 < 1-2 J& + 46 a representative point will perform discontinuous oscillations as depicted in fig. 2. This limit cycle was verified for a wide range of initial conditions at the analogue computer. Discontinuous oscillators of this kind have been proposed by Lavenda Nicolis and Herschkowitz-Kaufman in connection with the Brusselator and by Rossler l2 for an autocatalytic scheme with feed back inhibition.However this scheme is perhaps simpler and is well suited as a model of a branched chain reaction of the type occurring in many combustion reactions,' I33 CHEMICAL OSCILLATORS 3 I 5 A’ FIG.2.-Equation ofslow motion and limit cycle for eqn (16) with 4’ = 5. 5. EDELSTEIN’S EQUATIONS A novel biochemical scheme that shows multiple steady states has been deviscd by Edelsteh6 The scheme is characterised by the following reactions ki A+X$2X k-1 kz X+E+C (1 9) k-2 k3 C+E+B k-3 where E is an enzyme and C is the enzyme-substrate complex.It is assumed that the total amount of enzyme is conserved i.e. E+C = constant (ET).It is trivial to produce a discontinuous oscillator of the type described in the last section merely by letting A vary slowly subject to a flux into the system. We have in fact succeeded in producing oscillations at the analogue computer this way. However a more inter- esting case arises if the first eqn in (19) is made exothermic and the reverse reaction has a strong temperature dependence. The rate equations for this system are (assum- ing we choose the rate constants such that A and B can be treated as constant) dT C -= C Ahi(RT -Ri)-L(T -To) dt i where Tis the reactant temperature Tois the ambient temperature Cthe heat capacity n.F. GRAY AND L. J. AARONS of the mixture Ahi is the enthalpy of the ith reaction R+ and R; are the rates of the forward and reverse ith reaction respectively and L is the heat loss proportionality constant.l4 Assume that k-l is the only significant temperature dependent rate constant and that it can be written as 2 exp( -E,/T). If the reactants are in dilute solution such that Cis large we can choose I /Cas our small parameter. The equation of slow motion is given by (with all rate constants other than k- put equal to unity) z exp( -E,/T)X3-[A-(2 +B)Z exp( -E,/T)]x~ + [ET-A(2+B)]X -BE = 0 (21) ZE, E = ______. X+B+2' For certain values of the parameters a hysteresis curve results such as the one shown in fig. 3. The variation of T is given by dT C __-= Ahl(AX-Z exp(-E,/T)X2)-L(T-T,).(22) dt Again for particular choices of the parameters a discontinuous limit cycle results typified by the one shown in fig. 3. 302 3 03 T FIG.3.-X T phase plane and limit cycle for Edelstein's equations (20). Parametric values chosen A = 10.2 B = 0.2 ET = 35 Ea = lo4 2 = 3~ lo" Ah1 = 1 1/C = To = 300,L = 4. One may be misled into thinking that an oscillator is basically isothermal merely on the grounds that the temperature oscillations are of small amplitude. As has been shown here that may not be the case and one suspects that oscillators of this type may be more common than suspected due to the highly nonlinear energy equation. 6. CONCLUSIONS By way of conclusion we wish to stress several points made earlier.Our main point is that one should be careful when setting as constant the concentrations of CHEMICAL OSCILLATORS various " pool "chemicals. It may be possible that the inclusion of a small parameter may change a stable equilibrium point into an unstable one (the reverse certainly cannot occur). More commonly if the rate constants are not treated correctly the scheme including the small parameter will not reduce to the degenerate scheme as the small parameter is reduced to zero. The results obtained are interesting and include the discovery of a number of counterexamples to Hanusse's conjecture (e.g. very small percentage variations in the pool chemicals are sufficient to allow sustained oscillations in the active chemicals).Multistability of the two variable system seems to lend itself very well to the pro- duction of limit cycles of the relaxation type and bistable systems are shown to oscil- late in this manner when the necessary variation of "pool " chemicals is taken into account. Biochemical models employing allosteric enzymes are particularly good candidates for discontinuous oscillators of this type. 5-Conditions are derived under which the " pool " chemicals can be treated as constants as a zero'th approximation and these conditions place severe restrictions on the rate constants and may exclude regions where interesting instabilities have been discussed in the two variable case ; e.g. the " Brusselator " appears to be unlikely to " Brusselate '' in the region where it is permissible to treat the two variable system.Finally it is shown that very small amplitude temperature oscillations in solution are not grounds for assuming that the oscillator is basically isothermal with secondary temperature effects. Since the energy equation is so highly nonlinear it should be far easier to devise thermokinetic oscillators than isothermal ones. Again if the two variable system shows multistability it is possibIe to obtain limit cycle behaviour by allowing small temperature variations. It is highly significant that we are able to treat systems of more than two variables exactly. This depends on being able to separate the system into regions of slow and rapid motion. It is surprising that little work has been done using this theory when the equations are analytically soluble.One of us (L. J. A) is grateful to I.C.I. for the support of a postdoctoral fellowship. A. A. Andronov A. A. Vitt and S. E. Khaikin Theory ofOscillators (Pergamon Press Oxford 1966). 'R. Lefever J. Chem. Phys. 1968,49,4977. P. Glansdorff and I. Prigogine Thermodynamic Theory of Structure Stability and Fluctuations (Wiley-Interscience London 1971). 4P. Hanusse Compt. Rend. 1972 274C 1245. R. Lefever G. Nicolis and I. Prigogine J. Chem. Phys. 1967 47 1045. B. B. Edelstein J. Theoret. Biol. 1970 29 57. A. J. Lotka J. Amer. Chem. Soc. 1920,42,1595. B. F. Gray and C. H. Yang Combustion and Flame 1969 13,20. J. J. Tyson J. Chem. Phys. 1973,58 3919. lo R. Lefever and G. Nicolis J. Theoret. Biol. 1971 30 267.B. Lavenda G. Nicolis and M. Herschkowitz-Kaufman,J. Theoret. Biol. 1971 32 283. l2 0.E. Rossler J. Theoret. Biol. 1972 36,413. l3 B. F. Gray Kinetics of Oscillatory Reactions (Specialist Periodical ed. P. G. Ashmore Reports The Chemical Society London 1974). l4 R. &is Efementary Chemical Reactor Analysis (Prentice-Hall New Jersey 1969). l5 H. R. Karfunkel and F. F. Seelig J. Theoret. Biol. 1972 36 237. l6 L. J. Aarons and B F. Gray to be published.