Subject: Summary - Diffusion Equation from Kinetic Lattice Monte Carlo

From: Andrew Horsfield <[email-masked]>
Date: Fri, 17 Sep 1999 17:28:34 +0100 (BST)

First I would like to thank all those who responded to my question.

I have had several requests to let people know what I have learned in
response to the question about generating (approximate) diffusion /
reaction equations from Kinetic Lattice Monte Carlo models.

The original posting was:

> I have a question about kinetic Monte Carlo simulations on a lattice, 
> and how they relate to traditional diffusion and reaction equations.
>
> I understand Monte Carlo in the following way:
>
> 1) You define a set of allowed configurations which you can index
> 
> 2) You evaluate by some means a set of rates for transitions between
>     these configurations. Let the rate for a transition from
>     configuration i to configuration j be W_{ji} (note the order of the
>     indices).
> 
> 3) You start from some initial configuration, and make a transition to a
>    new configuration with a probability proportional to the rate. This
>    process is repeated for a number of transitions.
> 
> I want to derive continuum equations from this algorithm. To do this I
> begin by defining the probability that my system is in configuration i
> to be P_i. I then ask how P_i varies with time. The answer (I think) is:
> 
> d P_i
> ----- = \sum_j L_{ij} P_j
>  dt
> 
> where
> 
> L_{ij} = W_{ij}             i \= j
>        = -\sum_k W_{ki}     i  = j
> 
> This all looks fine until one transforms this equation into a
> corresponding equation for the rate of change of concentration (C):
> 
> d C(r)_n
> -------- = \sum_{sm} R_{rn,sm} C(s)_m
>    dt
> 
> where r and s are indices for the type of atom (C, O, H or whatever), n
> and m are indices for spatial position (lattice site), and R_{rn,sm} is
> a concentration independent rate matrix that can be derived from L_{ij}.
> The feature that worries me about this is that it is *linear* in C (as
> it has to be, given the starting point).
> 
> I believe that this cannot be correct if there are chemical reactions
> taking place, since then we would expect the rate to depend on the
> *product* of two concentrations (at least).
> 
> Can anyone tell me where I have gone wrong? Is there a standard
> reference for this kind of approximation?

I have found the mistake: I assumed that given a concentration profile it
was possible to determine a *unique* probability distribution of
configurations. This is true if you have only one particle (in which case
the diffusion equation is indeed linear), but not true if you have more
than one.

On the question of how to do things properly, I was pointed by Wolfgang
Huber of the IBM Almaden Research Center to a very comprehensive text. The
book is

> "Stochastic processes in physics and chemistry", by N. G. van Kampen,
> Elsevier Science Publishers B.V., 2n ed. 1992. He extensively treats
> the problem of how to derive macroscopic rate equations from master
> equations, the method is a systematic expansion being called
> "Omega-expansion". Higher-order terms in the expansion can even be
> related to Langevin-like fluctuations.

This is not a book for those who are shy of mathematics! But it is
comprehensive. Rather than displaying my lack of understanding again in
public, I refer you to this book.

Cheers,

Andrew

  +----------------------------------------------------+
   Andrew Horsfield       e-mail: [email-masked] 
     FECIT, 2 Longwalk Road, Stockley Park, Uxbridge,   
          Middlesex UB11 1AB, United Kingdom.           
   phone: +44-(0)181-606-4653  FAX: +44-(0)181-606-4422 
  +----------------------------------------------------+