Efficient simulation of semiflexible polymers

Using a recently developed bead-spring model for semiflexible polymers that takes into account their natural extensibility, we report an efficient algorithm to simulate the dynamics for polymers like double-stranded DNA (dsDNA) in the absence of hydrodynamic interactions. The dsDNA is modeled with one bead-spring element per base pair, and the polymer dynamics is described by the Langevin equation. The key to efficiency is that we describe the equations of motion for the polymer in terms of the amplitudes of the polymer's fluctuation modes, as opposed to the use of the physical positions of the beads. We show that, within an accuracy tolerance level of 5% of several key observables, the model allows for single Langevin time steps of ≈1.6, 8, 16, and 16 ps for a dsDNA model chain consisting of 64, 128, 256, and 512 base pairs (i.e., chains of 0.55, 1.11, 2.24, and 4.48 persistence lengths), respectively. Correspondingly, in 1 h, a standard desktop computer can simulate 0.23, 0.56, 0.56, and 0.26 ms of these dsDNA chains, respectively. We compare our results to those obtained from other methods, in particular, the (inextensible discretized) wormlike chain (WLC) model. Importantly, we demonstrate that at the same level of discretization, i.e., when each discretization element is one base pair long, our algorithm gains about five to six orders of magnitude in the size of time steps over the inextensible WLC model. Further, we show that our model can be mapped one on one to a discretized version of the extensible WLC model, implying that the speed-up we achieve in our model must hold equally well for the latter. We also demonstrate the use of the method by simulating efficiently the tumbling behavior of a dsDNA segment in a shear flow.