Bifurcation analysis of Gray-Scott reaction-diffusion cells

Aaron William Thornton

Bachelor of Mathematics Project



Nonlinear chemical kinetics have a long term interest across the globe. In particular, combustion theory involves factors such as temperature, space, time, interaction rates, growth rates, decay rates and many more. To understand how each factor affects the final products of a chemical reaction is very important for safe and efficient reactor operation. Especially when possible phenomena include oscillations, ignitions and chemical patterns. To maintain reactor performance the regions at which these phenomena take place needs to be determined in the parameter space. This gives the operator an ideal method of starting up and controlling the chemical reactor.

Many problems in chemical reactor engineering can be modelled by a system of coupled differential equations. These are established by considering molecular interaction and molecular diffusion using the law of mass action. Reaction rate coefficients are determined experimentally and can also be functions of temperature. Some famous chemical reactions include the Belousov-Zhabotinsky reaction, the H2+O2 reaction and the iodate-reductant reaction.

One of the more complicated chemical reaction schemes is the Belousov-Zhabotinsky (BZ) system, which involves the oxidation of an organic species such as malonic acid by an acidified bromate solution in the presence of a metal ion catalyst. When considering a closed system such as a reaction-diffusion cell where no chemicals are coming in or out, the reaction typically exhibits a short induction period followed by an oscillatory phase. For moderate flow rates in an open system like the continuously stirred tank reactor (CSTR), the reaction exhibits large amplitude oscillations very similar to that observed in a closed system. In an open system the oscillations are indefinitely sustained and each excursion having exactly the same amplitude and period as the previous one. The difference between the open and closed cases of the BZ system is that the build up of product concentrations and variations in background concentrations is now prevented by the continuous flow process. Also when considering low and high flow rates, complex oscillations emerge and in some cases the system can exhibit chaotic behaviour.

Once the governing equations are established, bifurcation theory is the standard technique used to separate regions in the parameter space where the dynamics of the system are qualitatively different. These different regions can include hysteresis and oscillatory behaviour. When considering systems like the CSTR, the governing equations will be ordinary differential equations (odes) where bifurcation theory can be directly applied. However for systems of partial differential equations (pdes) like the case considered in this thesis, semi-analytical methods are useful in transforming the pdes into an approximate model of odes. This gives rise to the opportunity of using local stability and singularity theory to determine the regions of parameter space in which the various bifurcation patterns and Hopf bifurcations occur. Semi-analytical methods are important in the division of these regions as numerical integration of the governing pdes would be computationally expensive, tedious and impractical.

An example of this analysis has been used on the porous catalytic pellet problem by Marchant & Nelson (2004). The model corresponds to a first-order exothermic reaction with Arrhenius kinetics and is described by two coupled reaction-diffusion equations for the temperature and the degree of reactant conversion. Here the Galerkin method is used to obtain the semi-analytical model that allows both steady-state temperature &l; conversion profiles and steady-state diagrams to be obtained. Stable limit cycles are discovered and analyzed while the hysteresis curve is created in the parameter space to distinguish where unique and breaking wave patterns occur.

An easily understood model scheme is the cubic autocatalysis reaction that is used in this thesis. Consider the concentrations of A, the reactant, and B, the autocatalyst, given by a & b respectively. There is only one step that can take place during the mixing process. That is,

A + 2B → 3B.

This means that a species (or molecule) A will interact with two species (or molecules) of B to form three molecules of B. The rate at which this will occur is given by,

Rate = kab2

where k is the reaction rate coefficient. Using this we know that concentration a decreases over time at the rate kab2 because a product A will disappear at this rate to form one more molecule of B. Also concentration b increases over time at the rate kab2 because an extra molecule of B is created during interaction. Therefore ignoring diffusion and decay of the products, the differential equation can be established.

da/dt = -kab2
db/dt = kab2

This cubic reaction will form the basis (together with diffusion and linear decay) of the model system studied here. A wide variety of behaviour like multiple steady-state and oscillations arises.

The Gray-Scott chemical reaction scheme, which represents cubic-autocatalysis with linear catalyst decay, has been much studied due to its multiple steady-state responses and oscillatory solutions. See Gray (1988) and Gray & Scott (1990) for reviews and descriptions of this work.

Steady-state stability and singularity analysis have been performed on the Gray-Scott model by Marchant (2002) where four qualitatively different bifurcation diagrams were discovered including oscillatory solutions. A fifth bifurcation diagram was also found in a extremely tiny parameter region with it's size being O(10-10). The two term semi-analytical model proved to have excellent accuracy in modelling the concentration profiles. Bifurcation predictions were also accurately calculated and local stability analysis found the location of Hopf bifurcations in the parameter space. Finally a stable limit cycle is considered in detail and demonstrates where the semi-analytical solutions are expected to be less accurate due to the concentration profiles inside the cell.

Finlayson & Merkin (1999) have completed an analysis on a special case of this model with an electric field that demonstrated the effects of changing the initial & boundary conditions. These results also proved the model to have multiple steady-state structures as well as periodic solutions. The study mainly focused on the effect the electric field had on travelling waves, with both pulse and front waves being considered. Extensive numerical simulations were used to test predictions from linear theory. Overall only travelling waves were found with periodic boundary conditions with speed and wavelengths increasing as the electric field strength increased.

The extension considered in this thesis will be of Marchant (2002) where an electric field is applied to the Gray-Scott system, where E is the magnitude of the electric field. This electric field causes a migration of ions that changes the concentration of the chemicals. By replicating the results for E=0 singularity diagrams will be compared with the extended E >0 results in chapter 4.3. The Galerkin method is used to reduce the governing pdes to a model consisting of odes in chapter 4.1. This is accomplished by approximating the spatial structure of the reactant and autocatalyst concentrations. Bifurcation analysis of the ode model is performed; singularity theory is used to analyse the static multiplicity while a stability analysis is used to determine the dynamic multiplicity. Comparisons with numerical solutions of the governing pdes confirms the accuracy and usefulness of the two term ode model while the one term ode model becomes less accurate as the applied electric field is increased. This will be explained by demonstrating that the concentration profiles flatten towards the center of the cell as the electric field is applied.

<< Move to the list of previous twleve credit-point projects.
<< Move to my start page.

Page Created: 17th October 2008.
Last Updated: 2nd December 2009.