3.3 Monte Carlo (MC), quantum Monte Carlo (QMC), and
kinetic Monte Carlo (KMC) methods
In the limit of steady-state, the master equation in (Eq. 1) can be
written in matrix form as WP = P or in the familiar Einstein’s
convention of wijPej =
Pei , where the summation over the
repeated index is implicitly assumed. It is important to recognize is
that W is the entire matrix and wij is the
ijth element of the matrix. Note that
wij is the transition probability of migrating from
microstate j to I, consistent with the definition of w in (Eq. 1).
Similarly, P is the entire vector of probabilities of each microstate,
and the ith element of P is Pi, the
probability to access microstate i. Note that here,
Pei is the equilibrium distribution.
More generally, if we start with a non-equilibrium state P(1), here, (1)
is the initial time and the system transitions to later times and is
tracked by (2); (3), etc., then: WP(1) = P(2);WP(2) = P(3);WP(3) = P(4);
…; WP(n) = P(n+1), and as n becomes large, P(n) = P(n+1) =
Pe. In a Monte Carlo simulation, we simulate the
system by sampling accessible microstates according to their equilibrium
distribution Pe, e.g., as given by the Boltzmann
distribution or the appropriate equivalent distribution in different
ensembles (for thermodynamic systems at equilibrium). To achieve this
task, we need to choose a W or all wij that make up the
W, such that Pei=
exp(-Ei/kBT)/[Σjexp(-Ej/kBT)] is satisfied. Metropolis
et al. [30] recognized that this could be achieved by choosing
wij that satisfy equation
wijPej =
Pei by imposing a stronger criterion,
namely,
Pemwnm=Penwmn,
which leads to the Metropolis Monte Carlo method for sampling
microstates of a classical system.
Quantum Monte Carlo (QMC) techniques provide a direct and potentially
efficient means for solving the many-body Schrödinger equation of
quantum mechanics [31]. The simplest quantum Monte Carlo technique,
variational Monte Carlo (VMC), is based on a direct application of Monte
Carlo integration to calculate multi-dimensional integrals of
expectation values such as the total energy. Monte Carlo methods are
statistical, and a key result is that the value of integrals computed
using Monte Carlo converges faster than by using conventional methods of
numerical quadrature, once the problem involves more than a few
dimensions. Therefore, statistical methods provide a practical means of
solving the full many-body Schrödinger equation by direct integration,
making only limited and well-controlled approximations.
The kinetic Monte Carlo (KMC) method is a Monte Carlo method computer
simulation intended to simulate the time evolution of processes that
occur with known transition rates among states (such as chemical
reactions or diffusion transport). It is essential 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 [32]. From a
mathematical standpoint, solving the master equation for such systems is
impossible owing to the combinatorially large number of accessible
microstates, even considering that the limited number of accessible
states renders the transition probability matrix sparse. The Gillespie
algorithm provides an ingenious way out of this issue. The practical
idea behind KMC is not to attempt to deal with the entire matrix, but
instead to generate stochastic trajectories that propagate the system
from state to state (i.e., a Markovian sequence of discrete hops to
random states happening at random times). From this, the correct time
evolution of the probabilities Pi(t) is then obtained by
ensemble averaging over these trajectories. The KMC algorithm does so by
selecting elementary processes according to their probabilities to fire,
followed by an updating of the time.