Chemical Physics Letters 377 (2003) 269–278

www.elsevier.com/locate/cplett

Modelling pattern formation in CO+O2 on Pt{1 0 0}

I.M. Irurzun a,b,*, R.B. Hoyle c, M.R.E. Proctor a, D.A. King b

aDepartment of Applied Mathematics and Theoretical Physics, Centre of Mathematical Sciences,

University of Cambridge, Wilberforce Road, Cambridge CB3 0WA, UKbChemistry Department, University of Cambridge, Lensfield Road, Cambridge CB2-1EW, UK

cDepartment of Mathematics and Statistics, University of Surrey, Guildford, Surrey GU2 7XH, UK

Received 28 April 2003

Published online: 25 July 2003

Abstract

We extend a detailed kinetic model for CO+O2 on Pt{1 0 0} to describe pattern formation. The model includes: (i) a

non-linear power law to describe the phase transition, (ii) trapping and untrapping processes explicitly considered, and

(iii) experimentally determined coverage-dependent sticking probabilities and rate constants. This model is extended to

include diffusion and gas global coupling. Diffusion is included through a mass-balance equation which couples the

migration of CO with the phase transition. Gas global coupling is introduced considering realistic values of the

pumping flow, the reactor volume and the size of the crystal.

Crown Copyright � 2003 Published by Elsevier B.V. All rights reserved.

1. Introduction

Spatiotemporal pattern formation in non-equi-

librium systems, such as surface chemical reac-

tions, is a well-documented phenomenon today[1–8]. For example, CO oxidation on single-crystal

platinum surfaces under ultrahigh vacuum condi-

tions exhibits a wide variety of spatial patterns and

waves, which have been observed experimentally

at high spatial resolution using photoemission

electron microscopy. In this system spatial cou-

pling is realized via two basic mechanisms: locally

by surface diffusion of CO and globally by pres-

* Corresponding author. Fax: +44-1233-765-900.

E-mail address: i.m.irurzun@damtp.cam.ac.uk (I.M. Irur-

zun).

0009-2614/$ - see front matter. Crown Copyright � 2003 Published

doi:10.1016/S0009-2614(03)01079-0

sure changes in the gas phase. Even before these

spatial features could be resolved, the spatially

averaged temporal features had already been un-

der scrutiny for some time. Work function and

mass spectrometric measurements revealed regularas well as chaotic oscillations in the average sur-

face coverage of the reacting species and in the rate

of CO2 production. The mechanism underlying the

oscillations on the {1 1 0} and {1 0 0} surfaces of Pt

has been well-established, and can be explained by

the existence of a reversible adsorbate-induced

phase transition in the crystalline structure of the

surface, which is attributed either to a criticalvalue of adsorbate coverage or to a strongly non-

linear dependence of the rate of phase transfor-

mation on adsorbate coverage. Rate oscillations

occur under conditions where oxygen adsorption

by Elsevier B.V. All rights reserved.

270 I.M. Irurzun et al. / Chemical Physics Letters 377 (2003) 269–278

is rate-limiting and since the oxygen sticking co-

efficient is very sensitive to Pt surface structure, the

phase transition causes a periodic switching be-

tween two states of different catalytic activity.

Based upon this mechanism, kinetic models were

developed in a so-called mean field approximation.Extensions to describe pattern formation typically

introduce Fickean terms to consider diffusion of

the mobile species (CO), and gas global coupling.

Major theoretical investigations were carried out

on the Krischer–Eiswirth–Ertl (KEE) model [9],

initially developed to describe oscillations in CO

oxidation on Pt{1 1 0}, in which the phase transi-

tion occurs between a 1� 1-structure and a 1� 2(missing row) structure. The model consists of a

set of three non-linear, coupled differential equa-

tions, describing variations in CO and oxygen

coverage and in the fraction of the surface in 1� 1structure. The inclusion of diffusion leads to the

reaction–diffusion equation (RDE) formulation,

which can also account for anisotropic diffusion

[10–16]. While this model describes nicely a num-ber of experimentally observed patterns, it is also a

simplistic representation of the real kinetic mech-

anism governing the chemical reaction. Other de-

velopments based in Monte-Carlo simulations

assume unrealistically low diffusion coefficients in

order to simulate patterns on computationally

accessible length and timescales [17]. To make

significant progress towards an adequate descrip-tion of the dynamics of these systems it is neces-

sary to address more detailed models that correctly

describe the reaction kinetics and account for all

the available experimental information. The de-

scription of pattern formation is also an additional

challenge to the development of detailed kinetic

mechanisms. In this Letter we study Pt{1 0 0}, in

which oscillations are observed associated with aphase transition between 1� 1 and hex surface

structures. The first kinetic model developed by

Ertl and coworkers [18] was later substantially

improved by King and coworkers [19–25]. The

main features included were: (i) a non-linear law

governing the phase transition as a function of the

local coverage of CO on hex phase, (ii) experi-

mentally determined, coverage-dependent stickingcoefficients and rate constants, and (iii) both

trapping and untrapping processes considered ex-

plicitly, independently of the phase transition.

Preliminary attempts to extend this model to de-

scribe pattern formation consider trapping and

untrapping processes in the development of a

general mass-balance equation that, by contrast

with RDE approaches, explicitly couples the phasetransition with diffusion [26]. The aim of this work

is to explore this new approximation. In order to

ensure mathematical stability, corrections to the

previous system of ordinary differential equations

(ODEs) are presented as we explain below. These

modifications show some general aspects that must

be considered in the formulation of a kinetic

model. The major achievement of the current workis the development of a spatially extended model

that takes into account the detailed experimental

evidence regarding the reaction kinetics. This

Letter mainly concerns the presentation of the

model. The behaviour of the system is studied

numerically in one spatial dimension. It was de-

sirable to focus initially on this simple case, be-

cause the complexity of the model presents certainnumerical challenges. The extension to two di-

mensions, though computationally intensive, is

straightforward, and together with a detailed in-

vestigation of the parameter space will form the

subject of a future work. We also study how gas

global coupling may affect pattern formation. The

present work is organized as follows. In Section 2,

the kinetic mechanism is described, ODEs arepresented and both local and global spatial cou-

pling mechanisms are introduced. Section 3 deals

with the details of the simulations. In Section 4

numerical results are presented and discussed.

Conclusions are summarized in Section 5.

2. The model

The kinetic model developed by King and co-

workers [19–25] to describe CO oxidation on

Pt{1 0 0} employs a scheme based on a Langmuir–

Hinshelwood mechanism with nine elementary

steps:

COðgÞ þIhex ()k1;k2

COhex ð1:1Þ

COðgÞ þI1�1 !k1CO1�1 ð1:2Þ

I.M. Irurzun et al. / Chemical Physics Letters 377 (2003) 269–278 271

CO1�1 !k3COðgÞ þI1�1;f ð1:3Þ

COhex þI1�1 ()k4;k5

CO1�1 þIhex ð1:4Þ

1

2O2ðgÞ þI1�1 !

k6O1�1 ð1:5Þ

1

2O2ðgÞ þI1�1;f !

k6O1�1 ð1:6Þ

CO1�1 þO1�1 !k7CO2ðgÞ þ 2I1�1;f ð1:7Þ

nCOhex þ mIhex !k8 nCO1�1 þ mI1�1 ð1:8Þ

I1�1 !k9

Ihex ð1:9Þwhere I symbolizes a free adsorption site and the

indices hex and 1� 1 refer to the hexagonal and1� 1 phases, respectively. Steps (1.1)–(1.3) repre-sent CO adsorption onto and desorption from hexand 1� 1 areas. The subscript f indicates freed

sites to distinguish them from free sites. Freed sites

are created by desorption of CO and reaction.

Their existence was incorporated into the kineticmechanism from experimental measurements of

the sticking probabilities of O2 and CO on COad

and Oad prepared surfaces. These experiments

showed that the sticking probabilities on the

CO-freed and oxygen-freed 1� 1 phases, called

S1�1;fO2, are almost identical. But it turned out

that S1�1;fO2is slightly higher than S1�1O2

, the sticking

probability on a successively oxygen precovered(1� 1) surface. A conversion process from freed to

free sites was neglected in our model assuming that

its timescale is long compared with the oscillatory

period. The number of freed sites increases with an

increasing fraction of the surface in the 1� 1phase, and both freed and free sites are destroyed

during the phase transition. We refer the reader to

[22] and references therein for a detailed descrip-tion of the mechanism and relevant experimental

information. Migration of CO from the hex phaseonto the 1� 1 phase (trapping) and the reverse

process (untrapping) are included in step (1.4).

CO reaction with oxygen to form CO2 appears in

step (1.7). The adsorbate induced Pt{1 0 0}

hex ! 1� 1 surface transition is represented by

steps (1.8) and (1.9). The experimental fact that the

1� 1-CO island growth rate follows a strongly

non-linear power law with an apparent reaction

order n � 4, was particularly taken into account.

The formation of subsurface oxygen [27,28] was

not included in our model. On Pt{1 0 0} subsurfaceoxygen formation is important at temperatures

above 540 K [29], which is far from the tempera-

ture range described by our model.

Based on this mechanism the following set of

ODEs is proposed:

d hhexCOhhex

� �dt

¼ k1pCOShexCOhhexð1 hhex

COÞ k2hhexCOhhex

k4BhhexCOð1 h1�1O h1�1CO Þ

þ k5Bh1�1CO ð1 hhexCOÞ; ð2:1Þ

d h1�1CO h1�1� �

dt¼ k1pCOS1�1CO h1�1ð1 h1�1O h1�1CO Þ

k3h1�1CO h1�1 k7h

1�1CO h1�1O h1�1

þ k4BhhexCOð1 h1�1O h1�1CO Þ

k5Bh1�1CO ð1 hhexCOÞ; ð2:2Þ

d h1�1O h1�1� �

dt¼ k6pO2fS1�1O2

ð1 h1�1f h1�1O h1�1CO Þ

þ S1�1;fO2h1�1f gh1�1 k7h

1�1CO h1�1O h1�1;

ð2:3Þ

d h1�1f h1�1� �

dt¼ k1S1�1CO pCO

hþ k6S

1�1;fO2

pO2ih1�1f h1�1

þ 2k7h1�1CO h1�1O h1�1 k4Bhhex

COh1�1f

þ k3h1�1CO h1�1; ð2:4Þ

d h1�1ð Þdt

¼ k8 hhexCO

� �nhhex if cP 1;

k9 1 cð Þh1�1 if c6 1;

�ð2:5Þ

c ¼ h1�1CO

hcritCO

þ h1�1O

hcritO

; ð2:6Þ

where hhexCO, h1�1CO , and h1�1O are the adsorbate cov-

erages on the hex and the (1� 1) phases, hhex and

h1�1 are the fraction of the surface in hex phase and(1� 1) phase, respectively, and h1�1f is the fraction

of the (1� 1) surface that is temporarily freed

by desorption or reaction. pCO and pO2 are the

272 I.M. Irurzun et al. / Chemical Physics Letters 377 (2003) 269–278

partial pressures of CO and O2, respectively.

Since the total surface is either (1� 1) or hex,h1�1 þ hhex ¼ 1. The adsorbate coverages refer to

local coverages on the hex and the (1� 1) phases.The sum h1�1 þ h1�1f h1�1 is the whole ensemble ofvacant sites on that fraction of the surface in the

(1� 1) phase h1�1. This yields the further expres-sion for the fraction of the surface in the (1� 1)phase h1�1 ¼ h1�1 þ h1�1f h1�1 þ h1�1O h1�1 þ h1�1CO h1�1.The fraction of the surface in the hex phase is

similarly given by hhex ¼ hhex þ hhexCOhhex, where hhex

is the fraction of vacant sites in the hex phase.Experimentally determined coverage dependencefor sticking probabilities (Sphasegas ) and rate constants

ki were used in our calculations. Their values arequoted in Table 1. Expressions for Sphasegas were ad-

justed from experimental data in [22]. A parameter

B is included in Eqs. (2.1), (2.2) and (2.4), as a

spatially averaged quantity, to take into ac-

count that the efficiency of CO migration from

the hex phase onto the 1� 1 phase (trapping)and the reverse process (untrapping) depends

on the boundary between the hex and (1� 1)areas. Since the local coverage of both phases

Table 1

Parameters used in the mathematical model (Eqs. (2))

Description Parameter

CO desorption hex k2CO desorption 1� 1 k3CO trapping k4CO untrapping k5Reaction k7hex ! 1� 1 k81� 1! hex k9Description Parameter

CO impingement rate k1O2 impingement rate k6

Sticking probabilities:

CO on hex ShexCO

CO on 1� 1 S1�1CO

O2 on 1� 1 S1�1O2

S1�1;fO2

(1� 1) boundary length B

Critical coverages:

for (1� 1)! hex hCOcrithOcrit

hex ! 1� 1, reaction order n

remains virtually constant throughout island

growth, B must be constant. A detailed discussion

of the determination of the parameter B has been

given [19].

This set of equations is very similar to thatinitially proposed by King and coworkers [22], but

we explicitly consider the fraction of sites available

to adsorption, trapping and untrapping processes.

It is possible to show that if such terms are ne-

glected the temporal behaviour of the system is not

significantly modified, and then they can be effec-

tively dropped out in a mean-field approach. In a

spatially extended system, however, their exclusionintroduces instabilities due to the violation of mass

balance in the equations. We note that the fact

that our model and that of King and coworkers

are very similar in the mean field limit means that

the experimental evidence supporting the model-

ling of the reaction kinetics in the earlier model is

also in agreement with the current work. An fur-

ther improvement of these equations should alsotake into account CO transformation from hex to1� 1 phase according to step (1.8) in Eqs. (2.1)

and (2.2), and consumption of freed sites during

Ea (kJ mol1) m (s1)

E2 ¼ 105 m2 ¼ 3:7� 1012

E3 ¼ 154ðh ¼ 0Þ m3 ¼ 1:0� 1015

0 m5m2S1�1CO =m3ShexCO

E3 E2 ¼ 49 m5 ¼ 1:0� 104

E7 ¼ 58:6 m7 ¼ 2:0� 109

0 m8 ¼ 4:9� 104

108 m9 ¼ 2:5� 1011

Value

2:22� 105 ML mbar1 s1

2:08� 105 ML mbar1 s1

0.78

0:91ðh ¼ 0Þ0:28ðh ¼ 0Þ0:31ðh ¼ 0Þ1

0.25 ML

0.4 ML

4.17

I.M. Irurzun et al. / Chemical Physics Letters 377 (2003) 269–278 273

the phase transition in Eq. (2.4). The inclusion of

these terms implies however, additional mathe-

matical difficulties related to the sharpness of the

phase transition. We expect that their inclusion

will not produce qualitative modifications in the

examples shown below, and will be a subsequentimprovement of the present model.

The above model has solutions that typically

take the form of temporal oscillations of the

variables as shown in Fig. 1. We have included a

comparison with previous model of King and co-

workers. The period and amplitude of the oscilla-

tions, which are dependent on partial pressures

and temperature, are in good agreement with ex-periments, as were those found in the model of

King and coworkers [22].

In order to include spatial coupling, we used

mass balance equations derived from trapping and

untrapping processes, which are written as:

Fig. 1. (a) Temporal oscillations obtained by numerical inte-

gration of Eqs. (2). The amplitude and the period of oscillations

depend on temperature and partial pressures. The dotted line

shows temporal oscillations obtained with the original model

developed by King and coworkers [21] for comparison. (b)

Homogeneous oscillations on a one-dimensional extended

system without GGC (L ¼ 200). Parameters: T ¼ 500 K,

pCO ¼ 105 mbar, pO2 ¼ 14:03� 105 mbar.

d hhexCOhhex

� �dt

¼ Bdr:ðk4rðhhexCO½1 h1�1O h1�1CO �ÞÞ

k5BdDðh1�1CO ð1 hhexCOÞÞ; ð3:1Þ

d h1�1CO h1�1� �

dt¼ k5BdDðh1�1CO ð1 hhex

COÞÞ

Bdr:ðk4rðhhexCOð1 h1�1O h1�1CO ÞÞÞ;

ð3:2Þ

d h1�1f h1�1� �

dt¼ Bdr:ðk4rðhhex

COh1�1f ÞÞ: ð3:3Þ

The fact that k4 depends on coverages via stickingprobabilities (see Table 1) is considered explicitlyin these equations. The parameter Bd is related to

B in Eqs. (2) via Bd � Bdl2, where dl is a typicallengthscale for the migration processes. We take

Bd ¼ 1 in our simulations. It must be pointed out

that these equations do not represent diffusion of a

single chemical adsorbate, but they represent more

general coupled mass balance equations involving

different chemical species. We note that by addingEqs. (3.1) and (3.2) the absolute coverage of CO

remains constant in accord with the kinetic equi-

librium between trapping and untrapping pro-

cesses. This zero mode can explain the sharpness

of the diffusion profiles shown below, and may

also be the origin of the numerical complexities of

this model. Also, these equations couple the phase

transition with the migration processes throughh1�1. A detailed explanation of a similar approach

has been given [26]. Gas global coupling (GGC) is

introduced taking into account variations in the

partial pressures of both CO and O2. Therefore,

pCO and pO2 are treated as additional variables.

The partial pressure balances consist of inflow and

outflow as well as terms due to reaction. They are

written as:

dpCOdt

¼ K½p0CO pCO�

a0

ZA

½k1pCOShexCOð1

hhex

COÞ k2hhexCO�hhex

a0

ZA

k1pCOS1�1CO ð1

h1�1O h1�1CO Þh1�1

þ a0

Zk3h

1�1CO h1�1

; ð4:1Þ

A

274 I.M. Irurzun et al. / Chemical Physics Letters 377 (2003) 269–278

dpO2dt

¼K½p0O2 pO2 �

a0

ZA

k6pO2S1�1O2

ð1n

h1�1f h1�1O h1�1CO Þh1�1o

a0

ZA

k6pO2S1�1;fO2

h1�1f h1�1n o

; ð4:2Þ

where p0CO and p0O2 are the partial pressures of COand O2 in the stationary state. The integral terms

model variations of pCO and pO2 due to adsorptionand desorption processes from both hex and

(1� 1) areas, according to Eqs. (1.1)–(1.3), (1.5)

and (1.6). The integrals extend over the whole

crystal surface and embody the possibility of glo-

bal synchronization. The parameter K is given by

K ¼ 1=s, where s is the residence time, which isdependent on the experimental setup (the pumping

flow and the volume of the chamber), and specify

the mean time needed to reach stationary state. In

our simulations, a pumping flow J ¼ 200 l s1 and

a volume of the chamber V ¼ 60 l yield K ¼ 10=3s1, corresponding to a realistic value of s ¼ 0:3 s[10]. The parameter a0 is given by:

a0 ¼ 1:379� 1019N0T=V ðmbar� cm2Þ; ð5Þ

where T is the temperature and N0 is the number ofmolecules per square centimetre of a complete

monolayer of adsorbed molecules.

We note at this point that while in Eqs. (4)

the integration area A corresponds to the whole

surface of the crystal, in the simulations the in-

tegrals are evaluated on an area which is actually

very small compared with a real crystal. For

example, while simulations represents areas of

� 100 lm2, a crystal has an area of � 100 mm2.

For this reason a0 must be corrected by multi-

plying by the coefficient c which is roughly

written as:

c ¼ AcrystalAsimulations

: ð6Þ

The value of a ¼ ca0 (and hence the effect of

GGC) is frequently underestimated because c is

not considered. For example for T ¼ 500 K,

V ¼ 60 l, N0 ¼ 2� 1015, Acrystal � 250 mm2,

Asimulation ¼ 50 lm2, we obtain:

a ¼ ca0 ¼ 11:5 mbar� cm2: ð7Þ

In our simulations, a is used as a control param-eter and varied around this value to modify the

strength of the gas global coupling.

3. Simulations

Simulations were made on one-dimensional sys-

tems with periodic boundary conditions. This

geometry represents a ribbon with an area of

A ¼ Dx� ðLDxÞ ’ 50 lm2. The differential equa-

tions were explicitly integrated using a second-order

finite difference scheme, and Eq. (2.5) was approx-imated using a continuous function of c to inter-polate between the two expressions for d h1�1ð Þ=dton either side of c ¼ 1. The spatial integration in

Eqs. (4) was performed with a procedure of order

O(1/N 4) (where N is the number of grid points in

one dimension). The minimum number of grid

points (L) used was 200. Results were checked

with larger lattice sizes to test convergence. Twodifferent initial conditions were used in this

work. To describe them and our results we use the

phase picture of oscillatory systems. The evolution

of systems describing temporal oscillations can be

represented in the phase space, in which time is

considered as an implicit variable. In this phase

space the system moves on a limit cycle and its state

can be completely described by a phase variable /and an amplitude variable R. At any instant thevalue of / gives the position on the limit cycle, and

one whole rotation corresponds to / varying be-

tween some initial value /0 and / ¼ /0 þ 2p. Inspatially extended systems / depends on the space

coordinate r: / ¼ /ðr; tÞ. The initial conditions inour simulations always correspond to phase inho-

mogeneities. This means that all variables haveinitial values lying on the limit cycle. Similarly Rdescribes the amplitude of the limit cycle. For a

system always moving on a limit cycle R is constant,but amplitude defects can be detected through a

spatial dependence of R : R ¼ Rðr; tÞ.

4. Results

To evaluate our model, we performed simula-

tions with different initial conditions and different

Fig. 3. Effect of GGC on travelling waves. From top to bot-

tom: (a) a ¼ 0 mbar� cm2, (b) a ¼ 3:33 mbar� cm2, (c)

a ¼ 33:3 mbar� cm2, (d) a ¼ 1:6� 102 mbar� cm2 (e)

a ¼ 3:3� 102 mbar� cm2. Other parameters are: T ¼ 500 K,

pCO ¼ 105 mbar, pO2 ¼ 14:03� 105 mbar. The initial condi-

tion is a phase perturbation in the first cell of the lattice. The

time interval shown is 1000 s, and the lattice size is L ¼ 200.

I.M. Irurzun et al. / Chemical Physics Letters 377 (2003) 269–278 275

intensities of gas global coupling (GGC). The

strength of gas global coupling was modified using

a as a control parameter. Even though this pro-

cedure does not represent the real situation in

which only temperature and partial pressures can

be changed for a given experimental setup, itprovides a simple way to understand general as-

pects of GGC. Our purpose here is to show that

general aspects of pattern formation in dissipative

systems can be reproduced with this model, while a

systematic exploration of the behaviour as func-

tion of temperature and partial pressures will be

the subject of a further work.

Starting with an almost uniform surface trav-elling waves are formed. The initial perturbation,

which in our case occupies the first cell in the lat-

tice, then propagates across the surface. Patterns

stabilize after decaying of transients and a mini-

mum integration time of about 500 s. Fig. 2 shows

the space–time diagram of a collision of two

travelling waves in the absence of gas global cou-

pling at two different pO2 . In two dimensions, theequivalent patterns to travelling waves are targets.

While gas global coupling obviously modifies

patterns on the surface, the existence of patterns

also influences the gas global coupling. The effect

of gas global coupling on travelling waves is shown

in Fig. 3. We also show the temporal evolution of

the partial pressure of CO. At sufficiently high

values of a homogeneous oscillations of the sur-face are recovered. The amplitude and period of

the oscillations however, decrease as a is increased.This phenomenon is accompanied by similar

changes in the evolution of pCO as is shown in

Fig. 4. Decreasing a, GGC breaks up and patterns

are stable structures on the surface. It was shown

Fig. 2. Space–time diagrams of travelling waves obtained

without GGC. (a) pO2 ¼ 14:03� 105 mbar and (b)

pO2 ¼ 15:0� 105 mbar, respectively. Other parameters are:

T ¼ 500 K, pCO ¼ 105 mbar. The time interval shown is 1000 s,

and the lattice size is L ¼ 200.

that the breakdown of gas global coupling is a

threshold phenomenon [12] dominated by some

critical value of acrit such that for a > acrit the as-ymptotic state of the system is a spatially uniform

oscillating state. Our findings are qualitatively in

accord with this result, and for T ¼ 500 K,

pCO ¼ 105 mbar, pO2 ¼ 14:03� 105 mbar, acritwas estimated to be acrit ¼ 12� 2 mbar� cm2.

The critical value, acrit, should be dependent on

temperature and partial pressures and an explo-

ration as a function of these parameters is inprogress. In our simulations, by increasing a, allpatterns disappear through the formation of an

amplitude instability to give homogeneous oscil-

lations, as has been also reported [14]. To test the

stability of the homogeneous oscillation solution,

simulations starting with stabilized travelling

waves were performed. In all cases, homogeneous

oscillations were obtained again, demonstratingthat this is the stable solution.

We also made one-dimensional simulations

starting with a phase gradient of 2p across the

surface and as a function of a. Our results aresummarized in Fig. 5. With this initial condition

and in the same range of parameters in which

travelling waves have been observed, the existence

of phase flips has previously been reported. Theyhave been systematically investigated for complex

Ginzburg–Landau equation (CGLE) in the

Fig. 4. Variations of pCO accompanying the simulations shown in Fig. 3. From top to bottom: (a) a ¼ 0 mbar� cm2, (b) a ¼ 3:33

mbar� cm2, (c) a ¼ 33:3 mbar� cm2, (d) a ¼ 1:6� 102 mbar� cm2 (e) a ¼ 3:3� 102 mbar� cm2. Other parameters are: T ¼ 500

K, pO2 ¼ 14:03� 105 mbar.

276 I.M. Irurzun et al. / Chemical Physics Letters 377 (2003) 269–278

presence of a global feedback, and on the KEE

model with GGC [13,14].

Phase flips are phase fronts with a height of

2p and a space–time dependence of the form/ðx; tÞ ¼ /ðx stÞ. Here s denotes the velocity ofthe front. Phase flips could be formed only in the

presence of a moderate GGC. They also become

unstable if the strength of GGC is sufficiently high

and disappear through the formation of an am-

plitude defect to give homogeneous oscillations. In

our model phase flips are not obtained for the

values of parameters studied so far. WithoutGGC, travelling waves are the stable pattern in

accord with previous results [12]. We note that for

the simulations studied so far, the frequency and

wavelength of travelling waves are independent of

initial conditions and only depend on the systemparameters. For the intermediate values of ahowever, phase flips could not be stabilized and

irregular patterns appear. This is ascribed to a

limitation of the model as we discuss below. At

sufficiently high values of a, homogeneous oscil-lations appear and patterns disappear through the

formation of an amplitude defect.

An important additional level of complexity inour model is the introduction of experimentally

Fig. 5. Effect of GGC. The initial condition is a phase gradient

of 2p across the surface. From top to bottom: (a) a ¼ 0

mbar� cm2, (b) a ¼ 33:3 mbar� cm2, (c) a ¼ 100

mbar� cm2. Other parameters are: T ¼ 500 K, pCO ¼ 105

mbar, pO2 ¼ 14:03� 105 mbar. The time interval shown is

1000 s, and the lattice size is L ¼ 200.

I.M. Irurzun et al. / Chemical Physics Letters 377 (2003) 269–278 277

determined, coverage-dependent sticking proba-

bilities and rate constants. Because of experimental

uncertainty these quantities should be considered

as control parameters. In this sense an explorationof the behaviour of the system as a function of

these parameters may lead to additional im-

provements of the model. Starting with a gradient

of 2p across the surface, and without GGC, we

found that travelling waves can be stabilized only

by maintaining non-vanishing sticking probabilities

in the entire coverage range 06 h < 1. We em-

phasize that this requirement only becomes evidentwhen considering spatially extended systems. Ad-

ditional requirements could be necessary to stabi-

lize solutions at moderate GGC.

5. Conclusions

The main result of our work is the effectiveextension of a detailed kinetic model to the de-

scription of pattern formation. This model con-

tains all the available experimental information

about the kinetic reaction mechanism, including

experimentally determined sticking probabilities

and rate constants. Diffusion was also introduced

through a general mass-balance equation which

couples migration processes with the phase tran-sition. We note that this equation can only be

formulated if:

(i) local instead of global coverages are used in

the formulation of the model,

(ii) both trapping and untrapping processes are

included in the mechanism,

(iii) both trapping and untrapping processes de-

pend only on the local coverage.

General aspects of pattern formation in dissi-pative systems can be reproduced. Though the

numerical analysis presented in this Letter is still

preliminary and additional work is necessary to

explore spatial behaviour, we believe that the in-

clusion of experimentally determined kinetic pa-

rameters is an advance with respect to previous

works, since a quantitative comparison with ex-

periments is now possible. In this regard, evenwhen experimental measurements on this system

exist [8,30,31] which are qualitatively in agreement

with our results, they were made at higher total

pressures. Because this strongly affects the values

of the propagation front velocity we avoid further

comparisons here, but note that new measure-

ments are then desirable. We also note that im-

provements to the original kinetic model should bemade in order to describe pattern formation more

accurately. The adequate representation of both

temporal and spatial behaviours is thus a challenge

to the development of kinetic mechanisms of sur-

face chemical reactions.

Additional improvements to our own model

should include transformation of CO from hex to1� 1 phase (according to step (1.8)), in Eqs. (2.1)and (2.2), and consumption of freed sites during

the phase transition in Eq. (2.4). This inclusion

implies however, additional mathematical diffi-

culties which are still under study.

Acknowledgements

This work was supported by the Leverhulme

Trust. We acknowledge Prof. Eduardo E. Mola

and Prof. Charles T. Campbell for very helpful

discussions.

References

[1] Reviewed by P. Gray, S.K. Scott, Chemical Oscillations

and Instabilities: Non-Linear Chemical Kinetics, Oxford

University Press, Oxford, 1990.

278 I.M. Irurzun et al. / Chemical Physics Letters 377 (2003) 269–278

[2] S.K. Scott, Chemical Chaos, Oxford University Press,

Oxford, 1991.

[3] D.P. Tse, A.M. Rucklidge, R.B. Hoyle, M. Silber, Phys. D

146 (2000) 367.

[4] R.B. Hoyle, Phys. Rev. E 58 (1998) 7315.

[5] S.B. Schwartz, L.D. Schmidt, Surf. Sci. 206 (1988) 169.

[6] R. Imbihl, G. Ertl, Chem. Rev. 95 (1995) 697.

[7] H.H. Rotermund, in: Pattern Formation in Continuous

and Coupled System, IMA-Series, vol. 115, Springer, New

York, 1999, p. 231.

[8] Chaos, Focus Issue: Non linear Pattern Formation in

Surface Science 12 (2002) 107.

[9] K. Krischer, M. Eiswirth, G. Erlt, J. Chem. Phys. 96 (1992)

9161.

[10] A. von Oertzen, H.H. Rotermund, A.S. Mikhailov, G.

Erlt, J. Phys. Chem. B 104 (2000) 3155.

[11] M. Falcke, H. Engel, J. Chem. Phys. 101 (1994) 6295.

[12] M. Falcke, H. Engel, Phys. Rev. E 50 (1994) 1353.

[13] M. Falcke, H. Engel, Phys. Rev. E 56 (2) (1997) 635.

[14] M. Bertram, Ph.D. Thesis, Technische Universit€aat Berlin

(2002).

[15] M. Kim, M. Bertram, M. Pollmann, A. von Oertzen, A.S.

Mikhailov, H.H. Rotermund, G. Ertl, Science 292 (2001)

1357.

[16] A.L. Lin, M. Bertram, K. Martinez, H.L. Swinney, A.

Ardelea, G.F. Carey, Phys. Rev. Lett. 84 (2000) 4240.

[17] V.P Zhdanov, Phys. Rev. E. 60 (1999) 7554;

V.P Zhdanov, Surf. Sci. 426 (1999) 345.

[18] R. Imbihl, M.P. Cox, G. Ertl, J. Chem. Phys. 83 (1985)

1578.

[19] A. Hopkinson, X.C. Guo, J.M. Bradley, D.A. King, J.

Chem. Phys. 99 (1993) 8262.

[20] A. Hopkinson, J.M. Bradley, X.C. Guo, D.A. King, Phys.

Rev. Lett. 71 (1993) 1597.

[21] D.A. King, Surf. Rev. Lett. 4 (1994) 435.

[22] M. Gruyters, T. Ali, D.A. King, Chem. Phys. Lett. 232

(1995) 1.

[23] M. Gruyters, T. Ali, D.A. King, J. Phys. Chem. 100 (1996)

14417.

[24] A.T. Pasteur, X.C. Guo, T. Ali, M. Gruyters, D.A. King,

Surf. Sci. 366 (1996) 564.

[25] M. Gruyters, D.A. King, J. Chem. Soc., Faraday Trans. 93

(1997) 2947.

[26] E.E. Mola, D.A. King, I.M. Irurzun, J.L. Vicente, Surf.

Rev. Lett. 10 (2003) 23.

[27] H.H. Rotermund, M. Pollmann, I.G. Kevrekidis, Chaos 12

(1) (2002) 157.

[28] A. von Oertzen, A.S. Mikhailov, H.H. Rotermund, G.

Ertl, J. Phys. Chem. B. 102 (1998) 4966.

[29] T. Lele, J. Lauterbach, Chaos 12 (1) (2002) 164.

[30] J. Lauterbach, PhD. Thesis, Free University, Berlin (1994).

[31] J. Lauterbach, H.H. Rotermund, Surf. Sci. 311 (1994) 231.