JOURNAL OF COMPUTATIONAL PHYSICS 127, 196–207 (1996) ARTICLE NO. 0168 Stochastic Simulation of Coupled Reaction–Diffusion Processes

actants is used to compute the time evolution of reactant concentrations. The stochastic algorithm is rigorous in theThe stochastic time evolution method has been used previously to study non-linear chemical reaction processes in well-stirred ho- sense that it provides an exact solution to the correspond-mogeneous systems. We present the first treatment of diffusion, in ing master equation for chemical reaction in a homoge-the stochastic method, for non-linear reaction–diffusion processes. neous, well-stirred reaction volume [6]. Because the Gilles-The derivation introduces mesoscopic rates of diffusion that are pie method follows unit-by-unit changes in the totalformally analogous to reaction rates. We map, using Green’s func-numbers of each reactant species, it is especially well suitedtion, the bulk diffusion coefficient D in Fick’s differential law to the to the study of systems in which reactant densities arecorresponding transition rate probability for diffusion of a particle between finite volume elements. This generalized stochastic algo- low and the application of methods based on continuum rithm enables us to numerically calculate the time evolution of a approximations, such as the traditional ordinary differen-spatially inhomogeneous mixture of reaction–diffusion species in tial equations of chemical kinetics, is questionable. This a finite volume. The algorithm is equivalent to solving the time capability is especially relevant to biophysics and cell biol-evolution of the spatially inhomogeneous master equation. A ogy. Within the intact living cell, number densities of keyunique feature of our method is that the time step is stochastic an