A hybrid smoothed dissipative particle dynamics (SDPD) spatial stochastic simulation algorithm (sSSA) for advection-diffusion-reaction problems

Abstract

We have developed a new algorithm which merges discrete stochastic simulation, using the spatial stochastic simulation algorithm (sSSA), with the particle based fluid dynamics simulation framework of smoothed dissipative particle dynamics (SDPD). This hybrid algorithm enables discrete stochastic simulation of spatially resolved chemically reacting systems on a mesh-free dynamic domain with a Lagrangian frame of reference. SDPD combines two popular mesoscopic techniques: smoothed particle hydrodynamics and dissipative particle dynamics (DPD), linking the macroscopic and mesoscopic hydrodynamics effects of these two methods. We have implemented discrete stochastic simulation using the reaction-diffusion master equations (RDME) formalism, and deterministic reaction-diffusion equations based on the SDPD method. We validate the new method by comparing our results to four canonical models, and demonstrate the versatility of our method by simulating a flow containing a chemical gradient past a yeast cell in a microfluidics chamber.

Keywords: Discrete Stochastic Simulation; Particle Based Fluid Dynamics; Reaction-Diffusion Master Equation.

Fig. 1.

Example of biological system possessing reaction, diffusion, advection, and stochastic dynamics. In an artery, a neutrophil (N) chemotaxes against the blood flow toward a bacterium (B) by sensing the gradient of individual bacterial peptides (blue triangles) that bind receptors (green Y-shapes) on the neutrophil cell surface.

Fig. 2.

Comparison between mathematical and physical models of (a) a typical deterministic approach, (b) SDPD, (c) sSSA, for the idealized case of a channel flow with diluted species. In the hybrid method, advection is solved using SDPD, while the mass transport can be solved using either a spatial stochastic algorithm or SDPD, depending if stochastic effects are relevant or not. In Figs. (a) and (b), the dotted parts denote magnified regions, illustrating the differences between the methods. In Fig (b), v_ij, Q_ij and R_ij denote the pairwise velocity, flux of mass and reaction term, respectively. In Fig. (c), p in the Kolmogorov equation denotes the probability that the system can be found in state X at time t.

Fig. 3.

(a) Schematic of degradation reaction of two arbitrary species, with initial and boundary conditions. (b) Plot of the concentration of the chemical species along the length of the cylindrical domain, at t = 0.2, 0.4, 0.6, 0.8 and 1.0[s] for N = 33 particles.

Fig. 4.

(a) Schematic of natural convection in a cylinder inside a square enclosure. (b), (c) and (d) show isocontour plots of dimensionless concentration C*at steady-state, for (b) Ra = 10⁴, (c) Ra = 10⁵ and (d) Ra = 10⁶. Results shown in (b), (c) and (d) were interpolated using a Gaussian kernel, for clarity purposes.

Fig. 5.

(a) Comparison of the mean dimensionless concentration (C*) and (b) dimensionless (re-scaled by $\sqrt{Ra/Sc}$) y-velocity $(v_y^{*})$ centerline profiles, for Ra = 10⁴, 10⁵ and 10⁶, computed at t = 50, over Nr = 100 realizations.

Fig. 6.

(a) Upper: schematic of our model of a yeast cell in a microfluidics channel flow containing a pheromone gradient. The working fluid is initially at rest, and at time t > 0, two inlets start to inject fluid into the chamber; the upper inlet consists of a diluted mixture of working fluid and alpha-factor L, and the lower inlet injects working fluid in the system. (a) Lower: magnified region around the cell, illustrating the particles and their initial velocities. (b) Contour plots of the dimensionless concentration of alpha-factor, $C_L^{*}$, at t* = 2, obtained using the proposed hybrid method (upper) and FEM (lower). (c) Comparison of the dimensionless concentration profiles of alpha-factor, $C_L^{*}$, between hybrid sSSA-SDPD, SDPD and finite element method (FEM), measured across the y* centerline, at time t* = 2. Results shown in (b) were interpolated using a Gaussian kernel, for clarity purposes.

Fig. 7.

Top: dimensionless concentration of species RL at the cell wall, $C_{RL}^{*}$, at t* = 1 (equivalent to physical time t = 0.4[s]), obtained by: (a) SDPD, (b) hybrid sSSA-SDPD method and (c) finite element method. Bottom: dimensionless concentration of species RL at the cell wall, $C_{RL}^{*}$, at t* = 2 (equivalent to physical time t = 0.8[s]), obtained by: (d) SDPD, (e) sSSA-SDPD method and (f) finite element method. Note that the deterministic SDPD simulations in (a) and (d) match closely the finite element results in (c) and (f). However, the sSSA-SDPD simulations in (b) and (e) are able to capture the discrete stochastic dynamics of this subcellular biological process.