Smoldyn
A spatial stochastic simulator for
chemical reaction networks
Introduction
Smoldyn is a computer program for simulating reactions and diffusion at
a mesoscopic size scale. This scale lies between the more
traditional microscopic size scale, in which all or most atoms are
considered, and the macroscopic size scale, in which molecules are
treated with concentrations rather than as individuals. Smoldyn
is based on and named for Smoluchowski dynamics, which is a very simple
description of chemical reaction dynamics. In Smoldyn, every
molecule of interest is represented as
an individual point, while those that are not of interest (water,
non-reactive molecules, etc.) are not represented.
Simulated space is
continuous, meaning that space is neither a lattice, nor
binned into
voxels. The use of continuous space provides good spatial
resolution and avoids artifacts that can arise from
lattice geometries. Space can be 1, 2, or 3-dimensional.
Some of the algorithms that Smoldyn executes are described below.
Diffusion
Smoldyn does not consider inter-molecular forces, so each molecule
diffuses independently of each other. If diffusion is isotropic,
this is simulated by displacing each molecule, on each axis, by
an amount that is chosen from a Gaussian probability density that has a
standard deviation of s = (2D∆t)1/2;
s is called the rms step
length, D is the diffusion
coefficient, and ∆t is the
simulation time step. Anisotropic diffusion is simulated by
choosing a Gaussian distributed random displacement vector, exactly as
in the isotropic case, and then multiplying it by the appropriate
anisotropy matrix. Anistropic diffusion occurs in liquid
crystals, aligned polymer media, and other ordered systems.
The pictures below show isotropic and anistropic diffusion,
respectively, from Smoldyn simulations. In the former, red,
green, and blue molecules were all started at the system center and
allowed to diffuse outwards. Red molecules were assigned a high
diffusion coefficient, green were assigned an intermediate value, and
blue a very small value. For the latter picture, red and green
molecules were again started at the system center and allowed to
diffuse outwards, but this time both were assigned anisotropic
diffusion matrices.

Surfaces
Surfaces in Smoldyn are built of panels, where each panel is a
rectangle, triangle, sphere, cylinder, hemisphere, or a 1- or
2-dimensional analog of these shapes. Triangular panels are
particularly useful for constructing complex mesh surfaces, while the
others can be more useful for designing surfaces manually.
Included in the Smoldyn download package is a small utility program
called wrl2smol which converts wrl type triangle mesh data to Smoldyn
format input.
When a
molecule collides with a surface, it can be reflected off of the
surface, absorbed by the surface, or transmitted through the
surface. Reflection is performed with ballistic-type reflection,
despite the fact that molecules are assumed to move exclusively through
Brownian motion, because it yields the correct molecule position
probability distribution. For absorption, Smoldyn used to absorb
a molecule into an absorbing surface if the molecule either crossed the
surface during a time step, or it was determined to have both crossed
and recrossed the surface during the time step. It was found that
the latter determination was difficult to work with, computationally
intensive, and made little difference to results, so it is no longer
implemented. Following are some images of molecules that are
diffusing among surfaces in 1, 2, and 3 dimensions, and in a complex
triangulated mesh (which was generously provided by Dylan Morris and
made Smoldyn-readable with wrl2smol).

Reactions
Smoldyn can simulate three types of fundamental reactions:
zeroth-order, unimolecular, and bimolecular. In zeroth-order
reactions, product molecules are produced spontaneously, which is
unphysical, but can be useful for simulating chemical inputs to the
system. They are simulated by adding a Poisson-distributed random
number of product molecules to the system at each time step. For
unimolecular reactions, the probability that a
molecule reacts is constant over time, which means that the probability
that a reactant is still in existance after some time decreases
exponentially. To simulate this, a reaction is considered to have
happened during a time step if a uniformly distributed random number is
greater than exp(-k1∆t), where k1 is the first-order
reaction rate constant.
A particular strong point of Smoldyn is its ability to accurately and
efficiently simulate bimolecular reactions. Following is a
typical sequence of events.

A) Red and green molecules are diffusing.
B) Through diffusion, the reactants approach each other to within a binding radius, σb, which causes a reaction to
occur. The size of the binding radius depends on the diffusion
coefficients, the reaction rate constant, and the simulation time
step. This radius is almost always smaller than the sum of the
physical radii of the two reactants.
C) The reaction product is a blue molecule, which is placed at the
reaction location.
D) The blue molecule diffuses.
E) The blue molecule dissociates into a red and a green molecule which
are initially separated by the unbinding
radius, σu. This radius is
typically slightly larger than the sum of the physical radii of the two
reactants.
F) The red and green product molecules of the dissociation reaction
diffuse away from each other. Alternatively, they could diffuse
back together and react again, which is called a geminate recombination.
The binding and unbinding radii are chosen by Smoldyn so that the
simulated steady-state reaction rates and geminate recombination
probabilities will equal values that the user requests.
The
following
images show simulations of the Lotka-Volterra reaction
system in a thin 3-D system and in a crowded 2-D system,
respectively. Both systems have periodic boundaries. Using
R for red and G for green, the reactions are simply:
R → 2 R, R + G → 2 G, and
G → nothing.

Papers that used or discuss Smoldyn:
Key: •• = research paper that used Smoldyn, • = Smoldyn referenced,
but not used.
2004
•• Andrews, Steven S. and Dennis Bray. "Stochastic simulation of
chemical reactions with spatial resolution and single molecule detail."
Phys. Biol. 1:137-151, 2004.
(Description and analysis of algorithms).
• Ander, M., P. Beltrao, B. Di Ventura, J. Ferkinghoff-Borg, M.
Foglierini, A. Kaplan, C. Lemerle, I. Tomás-Oliveira, and L.
Serrano. "SmartCell, a framework to simulate cellular processes
that combines stochastic approximation with diffusion and localisation:
analysis of simple networks" Syst.
Biol., 1:129-138, 2004. (Alternate software for spatial
stochastic simulation.)
• Rivas, Germán, Frank Ferrone, and Judith Herzfeld. "Life in a
crowded world." EMBO Reports
5:23-27, 2004. (Meeting report.)
2005
•• Andrews, Steven S. “Serial rebinding
of ligands to clustered receptors as exemplified by bacterial
chemotaxis.” Phys. Biol.
2:111-122, 2005. (Different code, but identical algorithm.)
•• Lipkow, Karen, Steven S. Andrews, and Dennis Bray. “Simulated
diffusion of phosphorylated CheY through the cytoplasm of Escherichia coli.” J. Bact., 187:45-53, 2005.
(Used Smoldyn.)
• Hazelbauer, Gerald L. “Myriad molecules in motion: Simulated
diffusion as a new tool to study molecular movement and interaction in
a living cell.” J. Bact.,
187:23-25, 2005. (Commentary.)
• Takahashi, Kouichi, Satya Nanda Vel Arjunan, Masaru Tomita.
“Space in systems biology of signaling pathways – towards intracellular
molecular crowding in silico” FEBS
Lett., 579:1783-1788, 2005. (Review.)
• Plimpton, Steven J. and Alex Slepoy. “Microbial cell modeling
via reacting diffusive particles” J. Phys.: Conf. Ser., 16:305-309,
2005. (About ChemCell.)
•• DePristo, Mark A., Lynne Chang, Karen Lipkow, Ronald D. Vale, and
Shahid Khan. “FRAP kinetics of the Escherichia coli chemoreceptor
cluster analyzed by simulated diffusion” Submitted, 2005.
(Used Smoldyn.)
• Le Novère, Nicolas and Dominic Tolle. “Particle-based
stochastic simulations” Proceedings of the 4th
Workshop on Computation of Biochemical Pathways and Genetic
Networks. Logos, Berlin, pp. 41-45. (Review.)
• Lemerle, Caroline, Barbara Di Ventura, and Luis Serrano. "Space as
the final frontier in stochastic simulations of biological systems." FEBS Lett. 579:1789-1794,
2005. (Review.)
2006
•• Lipkow, Karen. "Changing cellular location of CheZ predicted
by
molecular simulations." PLoS
Comp. Biol. 2:0301-0310, 2006. (Used Smoldyn).
•• Andrews, Steven S. and Adam P. Arkin. "Simulating cell biology" Curr. Biol. 16:R523-R527, 2006
(Used Smoldyn).
• Tolle, Dominic P. and Nicolas Le Novère. "Particle-based
stochastic simulation in systems biology" Curr. Bioinf. 1:315-320, 2006.
(Review.)
• Ridgway, Douglas, Gordon Broderick, and Michael J. Ellison.
"Accomodating space, time and randomness in network simulation." Curr. Opinion Biotechnol.
17:493-498, 2006. (Review.)
• Tournier, Alexander L., Paul W. Fitzjohn, and Paul A. Bates.
"Probability-based model of protein-protein interactions on biological
timescales." Algorithms Molec. Biol.
1:25, 2006. (Investigated an alternative reaction algorithm to the one
that Smoldyn uses.)
• Fange, David and Johan Elf. "Noise-induced Min phenotypes in E. coli." PLoS Comput. Biol. 2:637-648, 2006.
(Used Meso-RD simulation program.)
2007
• Grima, Ramon and Santiago Schnell. "A mesoscopic simulation approach
for modeling intracellular reactions." J. Stat. Phys. 128:139-164, 2007.
(Developed an alternative algorithm.)
• Bray, Dennis, Matthew D. Levin, and Karen Lipkow. "The chemotactic
behavior of computer-based surrogate bacteria." Curr. Biol. 17:12-19, 2007. (Used a
non-spatial simulation method.)
Smoldyn download
Platform: Compiles fairly easily for Macintosh or Linux
(instructions are in the user manual and Makefiles are provided).
An executable application is included in the download package for
Windows. Because Smoldyn uses the OpenGL library for graphics and
the libtiff library for saving tiff files, these need to be supported
by your computer. Mac users should already have OpenGL installed
but will need libtiff, Linux users may need both, and Windows users
don't need either but may need several .dll files instead. These
upgrade instructions are in the manual. Here's the libtiff
link: libtiff.
Current version: 1.85.
Documentation: The manual is in two parts. Part I
(pdf)
is for Smoldyn users and Part II (in the download package) is for
Smoldyn programmers. Some algorithms are
presented in a research paper that is available on my publications
page.
Download: A downloadable Smoldyn package with the
source
code, documentation, makefiles, test files, binaries, and the wrl2smol
utility is available as a .tgz format for Mac and Linux users, Smoldyn.tgz, or as a .zip format for
Windows users, Smoldyn.zip. There is a
sourceforge
website for Smoldyn here
but it's not active.
Old versions:
1.72
1.73
1.74
1.75
1.76
1.77
1.78
1.79
1.80
1.81
1.82
1.83
1.84
License: The
rxnparam.c file is in the public domain. The rest
of the code is Copyright 2003-2008 by Steven Andrews and distributed
under the Gnu Lesser General Public License.
Collaborators: The Computational Research Laboratories,
located in Pune, India, are adding parallel processing capability to
Smoldyn. Also, Upi Bhalla is
adding the Smoldyn simulation core into the MOOSE program.
History:
1.5 (7/2003) - First public release.
1.51 (9/5/03) - Unimolecular reactions can now have much lower reaction
rates, a memory
error upon termination was fixed, parameter checking was improved,
several minor bugs were fixed, and two new interpreter commands were
added.
1.52 (10/24/03) - A significant bug was fixed in the equation for
absorbing walls (thanks to Dan Gillespie for catching it), it is now
possible to pause simulations while using graphics, and a few minor
bugs were fixed.
1.53 (2/9/04) - Cleaned up code
for the runtime command interpreter, it is now possible to save file
stacks, minor display and user interface improvements, and several new
commands for text output and excluded volumes.
1.55 (8/20/04) - Minor graphics changes, some
improvements in the design of the command queue for the runtime command
interpreter, and it is now possible for commands to return error
messages.
1.56 (1/14/05) - Improved
graphics and a few other minor changes.
1.57 (2/17/05) - A new runtime command.
1.58 (7/22/05) - A new runtime command, a randomized Brownian dynamics
lookup table to remove artifacts that can arise from random number
generator problems, and a more versitile mol statement.
1.59 (?) - Displays the random number seed.
1.70 (5/17/06) - Reasonably good surface support, where surfaces can be
reflective, absorbing or transparent, and they are comprised of panels
that can be rectangular, triangular, or spherical.
1.71 (12/8/06) - Better surface support with cylinder and
hemisphere panels, and improved reflective surfaces.
1.72 (2/26/07) - Surfaces are finally totally functional for
both diffusion and reactions (or so I thought for a while).
1.73 (9/25/07) - Minor surface improvements, added exponential
time-stepping option for run-time commands, and improved multi-platform
compiling.
1.74 (10/22/07) - Fixed yet more bugs in surfaces, added disk-shaped
surface panels, and made some changes in surface input statement
formats. Only the code got updated for this release, so the
documentation is slightly out of date.
1.75 (11/6/07) - Added support of surface-bound molecules, which still
need some debugging. Also, changed several surface input
statement formats. The documentation is current, but incomplete.
1.76 (11/7/07) - Main source file was split to two files and assorted
function rearranging. More command-line options are now available.
1.77 (11/18/07) - Fixed a few minor bugs and modified several functions
to make Smoldyn easier to parallelize.
1.78 (11/30/07) - Fixed lots of bugs with surface-bound molecules and
commands, and improved input formats of several statements.
1.79 (12/6/07) - Fixed more bugs with commands. Also added basic
compartments.
1.80 (12/22/07) - 'Q' now quits the program, time steps can be changed
mid-simulation, and basic porting between Smoldyn and MOOSE was added.
1.81 (1/23/08) - Surface-bound molecules now work, lots of minor bugs
were fixed, multiple molecule lists were implemented, porting works
well, and documentation was improved.
1.82 (2/28/08) - Fixed a bug in which reaction products ended up across
surfaces. Also, reactions were completely overhauled and can now
be entered with a vastly improved format. Allosteric reactions
were added.
1.83 (3/14/08) - Improved unimolecular reaction rates for multiple
reaction channels, improved allosteric reactions, and added molecular
drift.
1.84 (4/11/08) - First Windows release. Also, added reactions
that are only active in specified compartments, added several commands,
and fixed the change time step command.
1.85 (6/3/08) - The Mersenne Twister is now used for random numbers,
added define statements to parser, renamed allosteric reactions to
conformational spread reactions, and fixed lots of minor bugs and one
significant one (reaction products ended up across surfaces).
Known bugs: No significant ones; a few very minor ones
are listed at end of part II of the documentation. If
you have a configuration file that causes crashes or incorrect
behavior, please send it to me so I can find the bugs. My e-mail
address is steven.s.andrews@gmail.com.
© Steven
Andrews, 2008. All rights reserved.
Last modified 6/3/08.