Computational Research Division, Lawrence Berkeley National Laboratory, 1 Cyclotron Road, Berkeley, California 94720, USA.
Department of Physics and Astronomy, San Jose State University, 1 Washington Square, San Jose, California 95192, USA.
J Chem Phys. 2017 Mar 28;146(12):124110. doi: 10.1063/1.4978775.
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.
我们基于用于随机流体动力学（FHD）的方法开发了用于随机反应扩散系统的数值方法。对于流体动力学系统，FHD 公式通过随机偏微分方程（SPDE）正式描述。在我们考虑的反应扩散系统中，当我们的 SPDE 被空间离散化并且反应被建模为具有泊松波动的源项时，我们的模型变得类似于反应扩散主方程（RDME）描述。然而，与 RDME 不同，随着分子数量的增加，RDME 变得非常昂贵，我们基于 FHD 的描述自然从波动强烈的区域扩展，即每个介观细胞具有少量（反应性）分子，到具有中等或弱波动的区域，最终扩展到确定性极限。通过隐式处理扩散，我们避免了限制所有基于扩散显式处理的方法的时间步长的严重限制，并构建了比 RDME 方法更有效的数值方法，而不会牺牲准确性。通过对线性化反应扩散模型的稳态波动分布的准确性进行分析，我们构建了几种两阶段（预测校正）方案，其中扩散使用随机 Crank-Nicolson 方法处理，反应由 Gillespie 的随机模拟算法或弱二阶 tau 跳跃方法处理。我们发现，在线性化设置中，隐式中点 tau 跳跃方案达到二阶弱精度，并为扩散分子的跳跃时间尺度量级大一个数量级的时间步长提供准确稳定的结构因子。我们研究了我们的方法在热力学平衡内外对 Schlögl 反应扩散模型的数值准确性。我们展示并量化了热力学波动对二维类似 Turing 模式形成的重要性，并研究了波动对三维化学前沿传播的影响。通过将随机模拟与确定性反应扩散模拟进行比较，我们表明波动会加速均匀系统中的模式形成，并导致在传播波后面形成定性不同的无序模式。
