Efficient Brownian dynamics simulation of DNA molecules with hydrodynamic interactions in linear flows

The coarse-grained molecular dynamics (MD) or Brownian dynamics (BD) simulation is a particle-based approach that has been applied to a wide range of biological problems that involve interactions with surrounding fluid molecules or the so-called hydrodynamic interactions (HIs). In this paper, an efficient algorithm is proposed to simulate the motion of a single DNA molecule in linear flows. The algorithm utilizes the integrating factor to cope with the effect of the linear flow of the surrounding fluid and applies the Metropolis method (MM) by Bou-Rabee, Donev, and Vanden-Eijnden [Multiscale Model. Simul. 12, 781 (2014)] to achieve more efficient BD simulation. Thus our method permits much larger time step size than previous methods while still maintaining the stability of the BD simulation, which is advantageous for long-time BD simulation. Our numerical results on $\lambda$-DNA agree very well with both experimental data and previous simulation results. Finally, when combined with fast algorithms such as the fast multipole method which has nearly optimal complexity in the total number of beads, the resulting method is parallelizable, scalable to large systems, and stable for large time step size, thus making the long-time large-scale BD simulation within practical reach. This will be useful for the study of membranes, long-chain molecules, and a large collection of molecules in the fluids.

Figure 1

Transient fractional extension for a seven-lambda ($L=150$ $\mu \mathrm{m}$) DNA in a planar extensional flow. 60 trajectories from simulations are used for ensemble average.

Figure 2

Comparison between experiments [11] (thin curves for individual trajectories and filled circles for the average) and our numerical simulations (empty circles). The vertical dashed line in (a) shows the point below which continuous data could not be collected in some experiments. The horizontal dashed line in (b) shows the steady-state of the stretched $\sim 22\mu \mathrm{m}$ $\lambda$-DNA.

Figure 3

Comparison between numerical results [6] (filled and empty circles for the assemble averages of FD and HI) and our numerical simulations (empty triangles). Solid curves are single trajectories of HI simulations from [6].

Figure 4

Fractional extensions for different De and shear rate $\dot{\gamma }$ from our simulations with $(\mathrm{De},\dot{\gamma })=(3.2,0.5)$, $(6.3,1.0)$, $(76.0,4.0)$, respectively. The relaxation time $\tau$ is $6.3\text{s}$ for the first two cases and $19.0\text{s}$ for the third case. The time steps are: $\mathrm{\Delta }t={10}^{-3}\text{s}$ for $\mathrm{De}=3.2$, $\mathrm{\Delta }t=5\times {10}^{-4}\text{s}$ for $\mathrm{De}=6.3$, and $\mathrm{\Delta }t=2.5\times {10}^{-4}\text{s}$ for $\mathrm{De}=76.0$.

Figure 5

Mean fractional extensions for shear flow and extensional flow. Experimental data [2] are symbols with error bars, bead model with and without HI (see legend) are from [6] and our results as red symbols (empty circles for the extensional flow; triangles and crosses for the shear flow).

Figure 6

Molecular extensions for 1.3 mm DNA in an extensional flow for $\mathrm{De}=0.30$ (top) and $\mathrm{De}=0.57$ (bottom). Filled symbols are experimental data from [10]. Blue circles are our simulation results. Time steps $\mathrm{\Delta }t={10}^{-2}\text{s}$ for $\mathrm{De}=0.30$ (top), and $\mathrm{\Delta }t=5\times {10}^{-3}\text{s}$ for $\mathrm{De}=0.57$ (bottom).

Figure 7

Numerical experiments of many-DNA in an oscillatory shear flow at $t=0$ and $t=25.6$, when shear flow velocity is zero.

Figure 8

Numerical experiments of many-DNA in an oscillatory shear flow at $t=38.4$ and $t=49.024$.

Figure 9

Numerical simulations for 25 DNA molecules in an oscillatory shear flow. Blue trajectories are fractional extensions of each molecule, dashed green curve is the magnitude of periodic shear flow and solid red curve is the assemble average.

Figure 10

Energies versus time $t$ for 25 DNA molecules in an oscillatory shear flow. Total energy (top) is the sum of WLC spring energy (middle) and EV potential energy (bottom). Red dotted vertical lines represent half period of flow oscillation.