Kinetic Monte Carlo

The kinetic Monte Carlo (KMC) method is a Monte Carlo method computer simulation intended to simulate the time evolution of some processes occurring in nature. Typically these are processes that occur with a given known rate. It is important to understand that these rates are inputs to the KMC algorithm, the method itself cannot predict them.

The KMC method is essentially the same as the dynamic Monte Carlo method and the Gillespie algorithm, the main difference seems to be in terminology and usage areas: KMC is used mainly in physics while the "dynamic" method is mostly used in chemistry.

Algorithm
The KMC algorithm for simulating the time evolution of a system where some processes can occur with known rates r can be written for instance as follows:

0. Set the time t = 0.

1. Form a list of all possible rates in the system ri

2. Calculate the cumulative function $$R_i=\sum_{j=1}^i r_j$$ for $$i=1,\ldots,N$$ where N is the total number of transitions. Denote R = RN

3. Get a uniform random number u' &isin; (0, 1]

4. Find the event to carry out i by finding the i for which Ri &minus; 1 < u'R &le; Ri (this can be achieved efficiently using binary search).

5. Carry out event i.

6. Recalculate all rates ri which may have changed due to the transition. If appropriate, add new transitions i. Update N accordingly.

7. Get a new uniform random number u &isin; ( 0, 1 ].

8. Update the time with t = t + &Delta;t where
 * $$\Delta t = - {\log u\over R}. $$

9. Return to step 1.

(note that the same average time scale can be obtained also using $$ \Delta t = {1 \over R} $$ in step 8. However, including the random number describes better the stochastic nature of the process).

This algorithm is known in different sources variously as the residence-time algorithm or the n-fold way or the Bortz-Kalos-Liebowitz (BKL) algorithm or just the kinetic Monte Carlo (KMC) algorithm.

Time dependent Algorithms
If the rate ri(t) are time dependent step 8 has to be modified by $$ \int_{0}^{\Delta t} R_i(t') dt' = -\log u.  $$. The reaction (step 4) has to be chosen after this by $$ \int_{0}^{\Delta t} R_{i-1}(t') dt' <    u' \int_{0}^{\Delta t} R_N(t') dt' \leq   \int_{0}^{\Delta t} R_i(t') dt'. $$

Another very similar algorithm is called the First Reaction Method (FRM). It consists of choosing the first occurring reaction, meaning to choose the smallest time $$ \Delta t_i $$, and the corresponding reaction number $$ i $$, from the formula $$ \int_{0}^{\Delta t_i} r_i(t') dt' = -\log u_i $$ where the ui &isin; ( 0, 1 ] are N random numbers.

Comments on the algorithm
The beauty of the KMC algorithm, and of the FRM one, is that if the rates are correct, the processes associated with the rates are of the Poisson process type, and different processes are independent, i.e. not correlated, the KMC algorithm gives the correct time scale for the evolution of the simulated system.

If furthermore the transitions follow detailed balance, the KMC algorithm can be used to simulate thermodynamic equilibrium. However, KMC is widely used to simulate non-equilibrium processes, in which case detailed balance need not be obeyed.

The KMC algorithm is efficient in the sense that every iteration is guaranteed to produce a transition. However, in the form presented above it requires $$N$$ operations for each transition, which is not too efficient. In many cases this can be much improved on by binning the same kinds of transitions into bins, and/or forming a tree data structure of the events.

The major disadvantage with KMC is that all possible rates $$r_i$$ and reactions have to be known in advance. The method itself can do nothing about predicting them.

Examples of use
KMC has been used in simulations of the e.g. the following physical systems:

1. Surface diffusion

2. Surface growth

3. Vacancy diffusion in alloys (this was the original use in (Young 1966))

4. Coarsening of domain evolution

5. Defect mobility and clustering in ion or neutron irradiated solids

To give an idea what the "objects" and "events" may be in practice, here is one concrete simple example, corresponding to example 2 above.

Consider a system where individual atoms are deposited on a surface one at a time (typical of physical vapor deposition), but also may migrate on the surface with some known jump rate $$w$$. In this case the "objects" of the KMC algorithm are simply the individual atoms.

If two atoms come right next to each other, they become immobile. Then the flux of incoming atoms determines a rate rdeposit, and the system can be simulated with KMC considering all deposited mobile atoms which have not (yet) met a counterpart and become immobile. This way there are the following events possible at each KMC step:


 * A new atom comes in with rate 'r''deposit


 * An already deposited atom jumps one step with rate w.

After an event has been selected and carried out with the KMC algorithm, one then needs to check whether the new or just jumped atom has become immediately adjacent to some other atom. If this has happened, the atom(s) which are now adjacent needs to be moved away from the list of mobile atoms, and correspondingly their jump events removed from the list of possible events.

Naturally in applying KMC to problems in physics and chemistry, one has to first consider whether the real system follows the assumptions underlying KMC well enough. Real processes do not necessarily have well-defined rates, the transition processes may be correlated, in case of atom or particle jumps the jumps may not occur in random directions, and so on. When simulating widely disparate time scales one also needs to consider whether new processes may be present at longer time scales. If any of these issues are valid, the time scale and system evolution predicted by KMC may be skewed or even completely wrong.

History
The first publication which described the basic features of the KMC method (namely using a cumulative function to select an event and a time scale calculation of the form 1/R) was by Young and Elcock in 1966 (Young 1966). The residence-time algorithm was also published at about the same time in (Cox 1965).

Apparently independent of the work of Young and Elcock, Bortz, Kalos and Lebowitz (Bortz 1975) developed a KMC algorithm for simulating the Ising model, which they called the n-fold way. The basics of their algorithm is the same as that of (Young 1966), but they do provide much greater detail on the method.

The following year Dan Gillespie published what is now known as the Gillespie algorithm to describe chemical reactions (Gillespie 1976). The algorithm is similar and the time advancement scheme essentially the same as in KMC.

There is as of the writing of this (June 2006) no definitive treatise of the theory of KMC, but Fichthorn and Weinberg have discussed the theory for thermodynamic equilibrium KMC simulations in detail in (Fichthorn 1991). A good introduction is given also by Art Voter (Voter 2005), and by A.P.J. Jansen (Jansen 2003),.