(362f) Accelerated Kinetic Monte Carlo (KMC) Algorithm for off-Lattice Particle- Based Reaction-Diffusion Systems

Accelerated
Kinetic Monte Carlo (KMC) Algorithm for Off-Lattice Particle- Based Reaction-Diffusion
Systems

Vikram
Thapar* and Paulette Clancy

School
of Chemical and Biomolecular Engineering

Cornell
University

Ithaca
NY 14853

*vt87@cornell.edu

We present an accelerated stochastic algorithm for simulating
off-lattice/particle-based reaction-diffusion network. Off-lattice Kinetic
Monte Carlo models are widely acknowledged to be valuable for a wide range of
applications (especially in materials growth applications and cell biology) for
which a fixed geometry lattice-based approach is restrictive. The most popular current
particle-based models include ChemCell, MCell, Smoldyn and DADOS. The
first three models (used largely for biological applications) make use of a fixed,
or adaptive-based, time step, whereas the commercial package DADOS (used almost
exclusively for Si-based semiconductor applications) uses a variable time step.
DADOS is based on a Gillespie SSA-like that "fires" just
one reaction or diffusion event in a time step. This, inevitably, makes it
computationally inefficient and limits its scope to complex problems.

Our newly developed algorithm offers a computationally
efficient, open source, code that makes use of a variable time step. The "tau-leaping"
approach developed by Gillespie for "well-mixed" systems lays the foundation
for our algorithm. The idea of "tau leaping" is based on firing multiple reaction-diffusion events in a
time interval chosen such that the rates of those events will not change by a significant
amount. Hence, the most critical and challenging aspect is to select an
appropriate value of tau. The expression for the time step we have implemented
is derived from concepts of how one event affects another. A traditional (single
event, rejection-free) KMC method is implemented for use as a "gold standard"
comparison to our accelerated simulations. Neighbor lists and data structures
are used to increase computational efficiency. The use of constrained binomial
random numbers to calculate the number of executions for each event takes care
of the issue of overlapping of particles and negative populations. Two
interesting reaction-diffusion examples – Fisher's Equation and the Lotka-Volterra (Predator-Prey) Equations -- are used to
check the performance and accuracy of the algorithm. We highlight improvements
in computational efficiency obtained by implementing this approach. Major
challenges applying this technique to a general reaction-diffusion system are
discussed, and possible solutions are presented.