Technische Universität München

Microkinetic Modeling

Karsten Reuter
Chemistry Department and Catalysis Research Center
Technische Universität München
 Catalysis Research Triangle

Testing
(Kinetics)

Preparation Characterization
 Macrokinetic studies

Fin [mol/min] Flow Fout [mol/min]

reactor

Batch reactor

creactants(t)
cproducts(t)

r = rate of reaction = f(initial gas composition, T, pi)

f(catalyst, active sites)
f(catalyst meso/macro structure)
f(reactor, flow (contact time), heat, dilution,…)

Power law kinetics:

r = k ∏ piα i
i
 From macro- to microkinetics

Reaction engineering: „Remove“ all reactor/catalyst macrostructure effects

to get intrinsic rate
Normalize to rate per active site = turnover frequency
TOF = r/N [molecules site-1 s-1]

→ intrinsic steady-state TOF = TOF(T, pi )

Microkinetic modeling aims to relate this

to detailed mechanistic understanding:

Top-down: deduce reaction mechanism from measured TOF data

Bottom-up: generate TOF data starting from first-principles calculations
 Disclaimer: Mind the gap!

Model catalysts

Supported „real“ catalyst

Courtesy: G. Rupprechter and Ch. Weiland, NanoToday 2, 20 (2007).

 I. The top-down approach:
Langmuir-Hinshelwood-Hougen-Watson
(LHHW) kinetics

Chemical kinetics and catalysis,

R.A. van Santen and J.W. Niemantsverdriet,
Plenum Press (New York, 1995)

Microkinetic simulation of catalytic reactions,

P. Stoltze, Prog. Surf. Sci. 65, 65 (2000)
 Step 1: Measured steady-state TOF data

TOF = r/N = f( pi , T )

CO oxidation over Rh(111)

TOF ~
10-2 – 102 molecules site-1 s-1 S.B. Schwartz, L.D. Schmidt, and G.B.
Fisher, J. Phys. Chem. 90, 6194 (1986)
 Step 2: Formulate a reaction mechanism and rate equation

CO + ½ O2 → CO2

O2(gas) + 2 (*) ↔ 2O(ad) (1) equilibrium

CO(gas) + (*) ↔ CO(ad) (2) equilibrium
CO(ad) + O(ad) → CO2(ad) + (*) (3) rate determining step
CO2(ad) → CO2(gas) + (*) (4) fast

r (CO2) ~ d/dt [CO2(gas)] = k4 θCO2(ad) first-order dependence (4)

d/dt θCO2(ad) = k3 θO(ad) θCO(ad) - k4 θCO2(ad) ≡ 0 steady-state approximation

r (CO2) ~ k3 θO(ad) θCO(ad)

⇒ determine adsorbate concentrations at the surface

in equilibrium with the gas phase
 Step 3: Derive adsorption isotherms

θeq = θeq( pi , T )

e.g. Langmuir: i) finite number of equivalent sites to hold reactants

ii) each site can hold at most one adsorbate
iii) no interaction between adsorbed particles

CO(gas) + (*) ↔ CO(ad) r adsorption = kads pCO(gas) θ(*)

r desorption = kdes θCO(ad)

r adsorption = r desorption (equilibrium)

θCO(ad) θCO(ad)
⇒ KeqCO = kads / kdes = =
pCO(gas) θ(*) pCO(gas) (1 - θ CO(ad))

Langmuir isotherm for

KeqCO pCO(gas)
θCO(ad) = non-dissociative,
1 + KeqCO pCO(gas) non-competitive adsorption
 Step 3: Adsorption isotherms cont’d

KeqCO pCO(gas)
θCO(ad) =
1 + KeqCO pCO(gas)

i) saturation due to finite number saturation

of sites linear
ii) KeqCO is an equilibrium property
→ thermodynamics

Similarly: molecular adsorption on several sites

dissociative adsorption (e.g. 2 sites)
competitive adsorption

Refinement: multiple sites → Freundlich/Toth isotherms

adsorbate interactions → Tempkin/Fowler isotherms
(only averaged interactions, still analytical)
→ Lattice gas Hamiltonians
(explicit interactions, no longer analytical)
 Step 4: Plug isotherms into rate equation

CO oxidation over Rh(111)

KeqCO pCO(gas)
r (CO2) ~ k3 θO(ad) θCO(ad) θCO(ad) =
1 + KeqCO pCO(gas) + KeqO2 pO2(gas)

KeqO2 pO2(gas)
θO(ad) =
1 + KeqCO pCO(gas) + KeqO2 pO2(gas)

Langmuir isotherms
competitive adsorption of O and CO
rapid O2 dissociation ( first order in θO(ad) )

Fits perfectly…

S.B. Schwartz, L.D. Schmidt, and G.B.

Fisher, J. Phys. Chem. 90, 6194 (1986)
 II: The bottom-up approach:
First-principles kinetic Monte Carlo simulations

First-principles kinetic Monte Carlo simulations for heterogeneous catalysis: Concepts, status and frontiers
K. Reuter, in “Modeling Heterogeneous Catalytic Reactions: From the Molecular Process to the Technical System”,
(Ed.) O. Deutschmann, Wiley-VCH, Weinheim (2009). http://www.fhi-berlin.mpg.de/th/paper.html
Elementary processes and catalytic function

- Continuous bond making and breaking

Challenge I: predictive-quality QM energetics

- Rare event time scale: > ~msec

Challenge II: long time-scale simulations

First-principles modeling of surface reactions

0.9 eV

COcus + Ocus → CO2

- Active site model

- „Level of theory“

M. Sabbe, M.F. Reyniers, and K. Reuter,

Catal. Sci. Technol. 2, 2010 (2012); invited perspective
 Markovian state dynamics: Kinetic Monte Carlo simulations

kB→A

kA→B
ΔEA→B ΔEB→A
kinetic
Transition State
Monte Carlo
Theory TS
N
 k T  Z TS ( i → j ) B
ki → j = B 
 h  Zi B A
 − ∆Ei → j 
= Γ exp   Molecular
 k BT  A Dynamics t

dPi (t )
= −∑ ki → j Pi (t ) + ∑ k j →i Pj (t )
dt j j
Kinetic Monte Carlo: essentially „coarse-grained MD“

Molecular Dynamics: Kinetic Monte Carlo:

the whole trajectory coarse-grained hops

ab initio MD: ab initio kMC:

up to 50 ps up to minutes
Building a first-principles kinetic Monte Carlo model

CO oxidation @ RuO2(110)

CO

O
O
Ru

26 elementary processes (site-specific):

- O2 adsorption/desorption (dissociative/associative)
- CO adsorption/desorption (unimolecular)
- O and CO diffusion
- CO + O reaction
K. Reuter, Oil&Gas Sci. Technol. 61, 471 (2006)
K. Reuter and M. Scheffler, Phys. Rev. B 73, 045433 (2006)
Surface structure and composition in the reactive environment

CO oxidation at RuO2(110)

600 K pO2 (atm)

10-15 10-10 10-5 1 10+5
105
CO br/CO
br/CO cus cus
CO

pCO (atm)
T = 600 K, pO 2 = 1 atm, pCO = 7 atm
1

CObr/-
10-5
Obr/ - Obr
O/O cus
br/O cus

K. Reuter, D. Frenkel and M. Scheffler,

Phys. Rev. Lett. 93, 116105 (2004)
 Steady-state and transient parameter-free turnover frequencies

350 K TPR
pCO (10-9 atm)
0.0 1.0 2.0 3.0
TOFCO2 (1012 mol/cm2 s)

6 pO2 = 10-10 atm

xx Exp.
4 x
xx

x
x
2 x
x
x Theory
x
x x
xx xxx x
0

K. Reuter and M. Scheffler, M. Rieger, J. Rogal, and K. Reuter,

Phys. Rev. B 73, 045433 (2006) Phys. Rev. Lett. 100, 016105 (2008)
 The dawn of a new era

J.K. Nørskov et al.,

Nature Chem. 1, 37 (2009)
 III. Towards error-controlled
first-principles microkinetic models
 Key ingredients to „predictive-quality“ microkinetic modeling

Accurate rate constants:  − ∆E i → j 

k i→ j = Γ exp  
 k BT 
Transition state theory and beyond
DFT functionals: „self-interaction“
van der Waals interactions

Reaction mechanism:
CO + * ↔ CO*
Process identification O2 + 2* ↔ O* + O*
Lattice mapping / spatial distributions CO* + O* ↔ CO2
…
„Hot chemistry“ beyond Markov
Mean-field approximation: Phenomenological rate equations

dPi (t )
= −∑ k i → j Pi (t ) + ∑ k j →i Pj (t )
dt j j

θ(Ocus) P (O br , COcus , t )
θ(Obr)
↓
θ(COcus)
θ(CObr) θ (O br , t ) ⋅θ (COcus , t )

dθ (O cus , t )
dt
{
= f 1 k i → j , θ (O cus , t ), θ (O br , t ), θ (CO cus , t ), θ (CO br , t ) }
dθ (O br , t )
dt
{
= f 2 k i → j , θ (O cus , t ), θ (O br , t ), θ (CO cus , t ), θ (CO br , t ) }

 Spatial arrangement and catalytic activity

T = 600 K → Inhomogeneous spatial

pO2 = 1 atm distribution even in
absence of
lateral interactions
in kMC model
This image cannot currently be displayed.

106
TOF (site-1s-1)

104
O-poisoned CO-poisoned

102

100

1 10 100
pCO (atm)
 The „power“ of fitting

CO oxidation
at RuO2(110)

TOF (site-1s-1) 104

102 MFA(kDFT) MFA(kFit)

kMC

pO2 = 1 atm, T = 600 K

100
0.5 5 100
pCO (atm)
Fitted rate constants deviate from “real” rate constants
by up to two orders in magnitude for dominant processes
Effective parameters without microscopic meaning B. Temel et al.,
J. Chem. Phys. 126, 204711 (2007)
 „E pluribus unum“: Water-gas-shift at Rh(111)

Microscale [CO]Rh + [OH]Rh CO Macroscale

COOH*+* → CO*+OH*
+ [H]Rh 1100
T catalyst

Temperature [°C]
1000
CO*+OH* → COOH*+*
COOH*+* → CO2*+H* 900

CO2*+H* → COOH*+* 800

CO2* + H2O* → COOH* + OH*
700
COOH* + OH* → CO2* + H2O*
CO2*+H* → HCOO** 600 T gas
HCOO** → CO2*+H*
[CO2]Rh +2[H]Rh [H]Rh +[OH]Rh + Rh 500
CO2* + OH* + *→ HCOO** + O*
400
HCOO** + OH* → CO2* + H2O* 2
CH* + H* → CH2* + *
CH* + * → C* + H*
C* + H* → CH* + *
RDS 1.5
CH
H
2
4

mmol/s
CH3* + O* → CH2* + OH* O CO
CH2* + OH* → CH3* + O*
CO2 1 2
CH* + OH* → CH2* + O*
CH2* + O* → CH* + OH* [H2O]Rh + 2Rh 0.5
HO
2
CO
2
…

~ 102 potential steps 0

0 2 4 6 8 10
@ different coverages Axial length [mm]
H2
H2O

M. Maestri and K. Reuter, Chem. Eng. Sci. 74, 296 (2012)

Diffusion at metal surfaces: surprises…

Hopping mechanism

Ag(100) ∆E = 0.45 eV
Au(100) ∆E = 0.83 eV

Exchange mechanism

Ag(100) ∆E = 0.73 eV
Au(100) ∆E = 0.65 eV

B.D. Yu and M. Scheffler, Phys. Rev. B 56, R15569 (1997)

 Automatized process identification

Accelerated molecular dynamics:

Hyperdynamics TAD

Other approaches: - metadynamics

- dimer method
…
Extending the Time Scale in Atomistic Simulation of Materials,
A.F. Voter, F. Montalenti and T.C. Germann,
Annu. Rev. Mater. Res. 32, 321 (2002)
 Error propagation through rate-determining steps

CO oxidation
at RuO2(110)

O-poisoned CO-poisoned
Sensitivity analysis:
pO2 = 1 atm
 k  ∂TOF  COcus + Ocus
X rc,i =  i   T = 600 K
 TOF  ∂ki  k j , K i Ocus/Ocus ad

C.T. Campbell,
Χrc,i

J. Catal. 204, 520 (2001); CObr ad

Nature 432, 282 (2004) COcus +Obr

Obr/Ocus ad

H. Meskine et al.,
Ertl Special Issue
pCO (atm) Surf. Sci. 603, 1724 (2009)
 Source for „rough“ rate constants: Hybrid UBI-QEP ?!

Pt(111),
Rh(111)

Unity Bond-Index Quadratic Exponential Potential

[ [ 
] − ] 
2
∆E∆E = φ =Eφ ( P − −Emin
b, A − − +
UBI-QEP ( P D )
UBI-QEP
A +→
AB→AB B A+B
TS
AB→ A+ B
E ( E ) ) (
b , BABD
=
MEP
φ
AB
E b ,
AB
AB ) 
+
x AB
 ( P D AB )
E. Shustorovich and H. Sellers, Surf. Sci. Rep. 31, 5 (1998)

M. Maestri and K. Reuter, Angew. Chemie Int. Ed. 123, 1226 (2011)
 More than Markov?

Heat dissipation during dissociative adsorption: O2/Pd(100)

k = S  (T )
~ pAuc
2πmk B T

Z=1.5 Å

J. Meyer and K. Reuter (submitted)

 First-Principles Microkinetic Modeling: Where do we stand?

State-of-the-art in catalysis modeling:

- Prevalence of highly coarse-grained models
based on effective parameters without true
microscopic meaning
rate equation theory
based on empirical rate constants
- Emergence of ad-hoc 1p-microkinetic models
kMC and mean-field
model catalysts, show case reactions

Steps towards a predictive character multiscale catalysis modeling:

- Replace effective parameters by clean first-principles data
fitted vs. DFT-based rate constants
battle the curse of complexity (off-lattice, complex networks)
electronic non-adiabaticity, heat dissipation
- Refined modeling at each individual level
reliable and efficient 1p-rate constants (where needed)
necessity to resolve spatial arrangement at surface
integrate 1p-surface chemistry into reactor models
- Robust links between theories that enable reverse-mapping
sensitivity analysis to control flow of error across scales
 Technische Universität München

www.th4.ch.tum.de

www.th4.ch.tum.de

Present members:
Max Hoffmann, Sebastian Matera, Jörg Meyer

Past members: Collaborations:

Matteo Maestri (→ U Milan, I) Daan Frenkel (Cambridge, UK)
Hakim Meskine (→ Wiley, D) Edvin Lundgren (U Lund, SE)
Michael Rieger (→ BASF, D) Horia Metiu (UCSB, USA)
Jutta Rogal (→ RU Bochum, D) Matthias Scheffler (FHI Berlin, D)