Stochastic simulation of reaction-diffusion systems: A fluctuating-hydrodynamics approach

Journal Article · Fri Mar 24 00:00:00 EDT 2017 · Journal of Chemical Physics

DOI:https://doi.org/10.1063/1.4978775· OSTI ID:1379783

Kim, Changho ^[1]; Nonaka, Andy ^[1]; Bell, John B. ^[1]; Garcia, Alejandro L. ^[2]; Donev, Aleksandar ^[3]

---

1. Lawrence Berkeley National Lab. (LBNL), Berkeley, CA (United States). Computational Research Division
2. San Jose State Univ., CA (United States). Dept. of Physics and Astronomy
3. New York Univ. (NYU), NY (United States). Courant Inst. of Mathematical Sciences

Here, we develop numerical methods for stochastic reaction-diffusion systems based on approaches used for fluctuating hydrodynamics (FHD). For hydrodynamic systems, the FHD formulation is formally described by stochastic partial differential equations (SPDEs). In the reaction-diffusion systems we consider, our model becomes similar to the reaction-diffusion master equation (RDME) description when our SPDEs are spatially discretized and reactions are modeled as a source term having Poisson fluctuations. However, unlike the RDME, which becomes prohibitively expensive for an increasing number of molecules, our FHD-based description naturally extends from the regime where fluctuations are strong, i.e., each mesoscopic cell has few (reactive) molecules, to regimes with moderate or weak fluctuations, and ultimately to the deterministic limit. By treating diffusion implicitly, we avoid the severe restriction on time step size that limits all methods based on explicit treatments of diffusion and construct numerical methods that are more efficient than RDME methods, without compromising accuracy. Guided by an analysis of the accuracy of the distribution of steady-state fluctuations for the linearized reaction-diffusion model, we construct several two-stage (predictor-corrector) schemes, where diffusion is treated using a stochastic Crank-Nicolson method, and reactions are handled by the stochastic simulation algorithm of Gillespie or a weakly second-order tau leaping method. We find that an implicit midpoint tau leaping scheme attains second-order weak accuracy in the linearized setting and gives an accurate and stable structure factor for a time step size of an order of magnitude larger than the hopping time scale of diffusing molecules. We study the numerical accuracy of our methods for the Schlögl reaction-diffusion model both in and out of thermodynamic equilibrium. We demonstrate and quantify the importance of thermodynamic fluctuations to the formation of a two-dimensional Turing-like pattern and examine the effect of fluctuations on three-dimensional chemical front propagation. Furthermore, by comparing stochastic simulations to deterministic reaction-diffusion simulations, we show that fluctuations accelerate pattern formation in spatially homogeneous systems and lead to a qualitatively different disordered pattern behind a traveling wave.

• Research Organization: Lawrence Berkeley National Laboratory (LBNL), Berkeley, CA (United States)
• Sponsoring Organization: USDOE Office of Science (SC), Advanced Scientific Computing Research (ASCR)
• Grant/Contract Number: AC02-05CH11231; SC0008271
• OSTI ID: 1379783
• Journal Information: Journal of Chemical Physics, Vol. 146, Issue 12; ISSN 0021-9606
• Publisher: American Institute of Physics (AIP)
• Country of Publication: United States
• Language: English

Related Subjects

37 INORGANIC, ORGANIC, PHYSICAL, AND ANALYTICAL CHEMISTRY
71 CLASSICAL AND QUANTUM MECHANICS
GENERAL PHYSICS