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

Abstract

We develop numerical methods for stochastic reaction-diffusion systems based on approaches used for fluctuatinghydrodynamics (FHD). For hydrodynamicsystems, 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 Poissonfluctuations. 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 thermodynamicfluctuations to the formation of a two-dimensional Turing-like pattern and examine the effect of fluctuations on three-dimensional chemical front propagation. 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.

Fluctuating Hydrodynamics and Debye-Hückel-Onsager Theory for Electrolytes

“…It should be noted that our FHD approach extends closely-related density functional Email address: donev@courant.nyu.edu (Aleksandar Donev) theory calculations of relaxation corrections to the conductivity [12] of binary electrolytes to account for advection, which enables us to also compute the electrophoretic corrections [2].First formulated by Landau and Lifshitz to predict light scattering spectra [13,14], more recently FHD has been applied to study various mesoscopic phenomena in fluid dynamics [10,15,16,17]. The current popularity of fluctuating hydrodynamics is due, in part, to the availability of efficient and accurate numerical schemes for solving the FHD equations [18,19,20,21,22,23,24,25,26]. While in this work we give analytical results for dilute electrolytes under simplifying assumptions, numerical techniques can in principle be used to compute the transport coefficients for moderately dilute solutions.From the work of Onsager et al [7,8,9] we can obtain the Fickian diffusion matrix for sufficiently dilute electrolytes.…”

“…First formulated by Landau and Lifshitz to predict light scattering spectra [13,14], more recently FHD has been applied to study various mesoscopic phenomena in fluid dynamics [10,15,16,17]. The current popularity of fluctuating hydrodynamics is due, in part, to the availability of efficient and accurate numerical schemes for solving the FHD equations [18,19,20,21,22,23,24,25,26]. While in this work we give analytical results for dilute electrolytes under simplifying assumptions, numerical techniques can in principle be used to compute the transport coefficients for moderately dilute solutions.…”

Fluctuating hydrodynamics of reactive liquid mixtures

“…for four examples and discuss various aspects of our numerical method and the effects of thermal fluctuations. First, we verify that for dilute solutions our algorithm preserves the robustness and accuracy properties of our previous method for reaction-diffusion systems [4] by modeling the hydrolysis of sucrose at micrometer scales. Second, to assess the fidelity of our approach in a non-dilute setting, we consider a binary mixture undergoing a dimerization reaction 2A ⇋ A 2 at thermodynamic equilibrium with a small number of molecules per cell.…”

“…The use of the CME enables us to correctly capture large fluctuations of composition, going beyond the Gaussian approximation inherent in the chemical Langevin equation (CLE) used in our prior work [12]. While our previous work on reaction-diffusion systems [4] also employed the CME, it was restricted to dilute solutions. Here we generalize the CME to non-dilute ideal mixtures with a complete Maxwell-Stefan formulation of diffusive transport in multispecies mixtures.…”

“…The key difference is that here we describe chemical fluctuations using a nonlinear FHD description based on a master equation, rather than a linearized Langevin description. This is common in stochastic reaction-diffusion models used in biochemical modeling [20,21], as we have discussed in more detail in prior work [4]. However, traditional RDME descriptions have been restricted to dilute solutions and do not account for velocity (momentum) fluctuations.…”