Bacterial Chemotaxis in Silico
Home Biology Software Data Publications Links

BCT | StochSim | Smoldyn



The computer program Smoldyn (for Smoluchowski dynamics) was written by Steve Andrews when he was a postdoc in our group 2001-2003. It was developed as part of a study of bacterial chemotaxis as a more realistic way to represent the diffusion of signalling molecules through the cytoplasm - especially CheY. The program now provides a general purpose biochemical simulator in which the diffusion of molecules in 1-, 2-, or 3- dimensions can be computed and displayed graphically. Molecules move randomly at rates specified by their diffusion coefficients and bind or react with other molecules according to probabilities derived from known rate constants. The ability to represent reactions in solution is so far as we know a unique feature of Smoldyn.

Smoldyn consists of a core simulation engine encapsulating the algorithm described below. It is platform-independent and written in C and OpenGL. The compiled program, documentation and source code are available for download from

Description of the algorithm

Smoldyn employs the Smoluchowski level of detail. That is, molecules have an identity and an exact position in continuous space, but no volume, shape or inertia. They diffuse in random directions by distances calculated from Fick's second law rewritten as a stochastic master equation.

Stochastic master equation

The product pB(r,t)dr is the probability that a specific B molecule is within a volume dr about position r at time t. Solving the above equation shows that the probability density for the displacement of a molecule after a time step has a Gaussian profile on each Cartesian coordinate (Andrews & Bray, 2004). These results form the basis of the simulation method called Brownian dynamics in which diffusion is simulated by picking a normally distributed random displacement for each molecule at each step. Since space is continuous, not compartmentalised, the level of detail can be adjusted by a suitable choice of step time dt.

Unimolecular reactions occur every time step at a probability that is derived from the user-defined rate-constant. Bimolecular reactions occur if the two molecules come within a binding radius, which is dependent on time-step length and rate constant. The algorithms employed are based on precise statistical-mechanical calculations of the rates of uni- and bimolecular reactions (Andrews & Bray, 2004).

Iteration step of the Smoldyn algorithm

Typical sequence of steps in Smoldyn. Two molecules, red and green, shown in the first box diffuse to new positions in the second box. If they come within a binding radius, shown in the third box, then they may react, depending on a reaction probability. If the reaction does occur, then the two molecules are replaced by a product molecule (yellow, as seen in the last box).

Running Smoldyn

The user writes a configuration file that defines all molecular species, their diffusion coefficients and initial positions in the cell, followed by runtime commands and specification of all reactions. Each molecule and molecular complex is a separate software species, and this includes variants created by conformational states and posttranslational modifications. In each time step, all mobile molecules (those with a non-zero diffusion coefficient) are moved to a new position. If any molecule enters an 'excluded volume' it is returned to its previous position. Reactions are then carried out where appropriate and the list of molecules updated. Products are positioned according to an unbinding radius, which prevents automatic re-binding in the next time step for reversible reactions. Before the next diffusion step, runtime commands are performed, such as change of molecule identities and output of molecule numbers and positions to a text file. If required, Smoldyn simulations can be displayed simultaneously as graphical animations (in OpenGL) and saved as a series of TIFF files.

Membrane molecules such as the receptor cluster are placed at fixed positions just inside the cell walls. (Their distance from the cell walls must be at least as large as their largest binding radius for correct rates.) With a diffusion coefficient of zero, the membrane proteins remain at their position, but can change identity spontaneously or in reactions. Currently, Smoldyn does not allow interactions with anything outside the cell volume, such as ligands. Instead, the user can generate input to the pathway by explicitly changing the identity of molecules from one activity state to another, either directly or probabilistically.

Sample output

Sample output from the Smoldyn program

Sample output from the Smoldyn program. This is a snapshot taken from an OpenGL movie showing the receptor cluster at one end of the cell and signalling molecules diffusing in the cytoplasm. Colour code: brown CheY, red CheYp, green CheZ, blue motors. The same data can be analysed in many ways, for example, to show the loci of individual molecules and the lifetimes of phosphorylated species (Lipkow et al., 2005).

Timecourses of numbers of molecules in response to a pulse of aspartate

Timecourses showing changes in numbers of molecules in response to a pulse of the attractant aspartate. Data produced by Smoldyn were in this case analysed by MATLAB.

Spatial distribution of CheYp molecules

Distribution of CheYp molecules in a cell containing freely diffusing CheZ molecules. Since phosphorylation occurs mainly at one end of the cell (the receptor cluster), a gradient of phosphorylated CheY has been set up.


Andrews, S. S. & Bray, D. (2004) Stochastic simulation of chemical reactions with spatial resolution and single molecular detail. Phys. Biol. 1, 137-151.

Lipkow, K., Andrews, S. S. & Bray, D. (2005) Simulated diffusion of phosphorylated CheY through the cytoplasm of Escherichia coli. J. Bacteriol. 187, 45-53.

Smoluchowsky, M. v. (1916) Versuch einer mathematischen Theorie der Koagulationskinetik kolloider Lösungen. Z. Phys. Chem. 92, 129-168.

Site designed by Matthew Levin. This page updated 11 January 2011.
Valid HTML 4.01 Transitional