11 |
|
simulations. Implicit solvent Langevin dynamics simulations of |
12 |
|
met-enkephalin not only outperform explicit solvent simulations for |
13 |
|
computational efficiency, but also agrees very well with explicit |
14 |
< |
solvent simulations for dynamical properties\cite{Shen2002}. |
14 |
> |
solvent simulations for dynamical properties.\cite{Shen2002} |
15 |
|
Recently, applying Langevin dynamics with the UNRES model, Liow and |
16 |
|
his coworkers suggest that protein folding pathways can be possibly |
17 |
< |
explored within a reasonable amount of time\cite{Liwo2005}. The |
17 |
> |
explored within a reasonable amount of time.\cite{Liwo2005} The |
18 |
|
stochastic nature of the Langevin dynamics also enhances the |
19 |
|
sampling of the system and increases the probability of crossing |
20 |
< |
energy barriers\cite{Banerjee2004, Cui2003}. Combining Langevin |
20 |
> |
energy barriers.\cite{Banerjee2004, Cui2003} Combining Langevin |
21 |
|
dynamics with Kramers's theory, Klimov and Thirumalai identified |
22 |
|
free-energy barriers by studying the viscosity dependence of the |
23 |
< |
protein folding rates\cite{Klimov1997}. In order to account for |
23 |
> |
protein folding rates.\cite{Klimov1997} In order to account for |
24 |
|
solvent induced interactions missing from implicit solvent model, |
25 |
|
Kaya incorporated desolvation free energy barrier into implicit |
26 |
|
coarse-grained solvent model in protein folding/unfolding studies |
27 |
|
and discovered a higher free energy barrier between the native and |
28 |
|
denatured states. Because of its stability against noise, Langevin |
29 |
|
dynamics is very suitable for studying remagnetization processes in |
30 |
< |
various systems\cite{Palacios1998,Berkov2002,Denisov2003}. For |
30 |
> |
various systems.\cite{Palacios1998,Berkov2002,Denisov2003} For |
31 |
|
instance, the oscillation power spectrum of nanoparticles from |
32 |
|
Langevin dynamics simulation has the same peak frequencies for |
33 |
|
different wave vectors, which recovers the property of magnetic |
34 |
< |
excitations in small finite structures\cite{Berkov2005a}. |
34 |
> |
excitations in small finite structures.\cite{Berkov2005a} |
35 |
|
|
36 |
|
%review langevin/browninan dynamics for arbitrarily shaped rigid body |
37 |
|
Combining Langevin or Brownian dynamics with rigid body dynamics, |
40 |
|
as well as excluded volume potentials, Mielke and his coworkers |
41 |
|
discovered rapid superhelical stress generations from the stochastic |
42 |
|
simulation of twin supercoiling DNA with response to induced |
43 |
< |
torques\cite{Mielke2004}. Membrane fusion is another key biological |
43 |
> |
torques.\cite{Mielke2004} Membrane fusion is another key biological |
44 |
|
process which controls a variety of physiological functions, such as |
45 |
|
release of neurotransmitters \textit{etc}. A typical fusion event |
46 |
< |
happens on the time scale of millisecond, which is impractical to |
46 |
> |
happens on the time scale of a millisecond, which is impractical to |
47 |
|
study using atomistic models with newtonian mechanics. With the help |
48 |
|
of coarse-grained rigid body model and stochastic dynamics, the |
49 |
|
fusion pathways were explored by many |
50 |
< |
researchers\cite{Noguchi2001,Noguchi2002,Shillcock2005}. Due to the |
50 |
> |
researchers.\cite{Noguchi2001,Noguchi2002,Shillcock2005} Due to the |
51 |
|
difficulty of numerical integration of anisotropic rotation, most of |
52 |
|
the rigid body models are simply modeled using spheres, cylinders, |
53 |
|
ellipsoids or other regular shapes in stochastic simulations. In an |
56 |
|
dynamics simulation algorithm\cite{Ermak1978,Allison1991} by |
57 |
|
incorporating a generalized $6\times6$ diffusion tensor and |
58 |
|
introducing a simple rotation evolution scheme consisting of three |
59 |
< |
consecutive rotations\cite{Fernandes2002}. Unfortunately, unexpected |
59 |
> |
consecutive rotations.\cite{Fernandes2002} Unfortunately, unexpected |
60 |
|
errors and biases are introduced into the system due to the |
61 |
|
arbitrary order of applying the noncommuting rotation |
62 |
< |
operators\cite{Beard2003}. Based on the observation the momentum |
62 |
> |
operators.\cite{Beard2003} Based on the observation the momentum |
63 |
|
relaxation time is much less than the time step, one may ignore the |
64 |
|
inertia in Brownian dynamics. However, the assumption of zero |
65 |
|
average acceleration is not always true for cooperative motion which |
66 |
|
is common in protein motion. An inertial Brownian dynamics (IBD) was |
67 |
|
proposed to address this issue by adding an inertial correction |
68 |
< |
term\cite{Beard2000}. As a complement to IBD which has a lower bound |
68 |
> |
term.\cite{Beard2000} As a complement to IBD which has a lower bound |
69 |
|
in time step because of the inertial relaxation time, long-time-step |
70 |
|
inertial dynamics (LTID) can be used to investigate the inertial |
71 |
|
behavior of the polymer segments in low friction |
72 |
< |
regime\cite{Beard2000}. LTID can also deal with the rotational |
72 |
> |
regime.\cite{Beard2000} LTID can also deal with the rotational |
73 |
|
dynamics for nonskew bodies without translation-rotation coupling by |
74 |
|
separating the translation and rotation motion and taking advantage |
75 |
|
of the analytical solution of hydrodynamics properties. However, |
172 |
|
one can write down the translational and rotational resistance |
173 |
|
tensors |
174 |
|
\begin{eqnarray*} |
175 |
< |
\Xi _a^{tt} = 16\pi \eta \frac{{a^2 - b^2 }}{{(2a^2 - b^2 )S - 2a}}. \\ |
176 |
< |
\Xi _b^{tt} = \Xi _c^{tt} = 32\pi \eta \frac{{a^2 - b^2 }}{{(2a^2 - 3b^2 )S + |
175 |
> |
\Xi _a^{tt} & = & 16\pi \eta \frac{{a^2 - b^2 }}{{(2a^2 - b^2 )S - 2a}}. \\ |
176 |
> |
\Xi _b^{tt} & = & \Xi _c^{tt} = 32\pi \eta \frac{{a^2 - b^2 }}{{(2a^2 - 3b^2 )S + |
177 |
|
2a}}, |
178 |
|
\end{eqnarray*} |
179 |
|
and |
180 |
|
\begin{eqnarray*} |
181 |
< |
\Xi _a^{rr} = \frac{{32\pi }}{3}\eta \frac{{(a^2 - b^2 )b^2 }}{{2a - b^2 S}}, \\ |
182 |
< |
\Xi _b^{rr} = \Xi _c^{rr} = \frac{{32\pi }}{3}\eta \frac{{(a^4 - b^4 )}}{{(2a^2 - b^2 )S - 2a}}. |
183 |
< |
\end{eqnarray} |
181 |
> |
\Xi _a^{rr} & = & \frac{{32\pi }}{3}\eta \frac{{(a^2 - b^2 )b^2 }}{{2a - b^2 S}}, \\ |
182 |
> |
\Xi _b^{rr} & = & \Xi _c^{rr} = \frac{{32\pi }}{3}\eta \frac{{(a^4 - b^4 )}}{{(2a^2 - b^2 )S - 2a}}. |
183 |
> |
\end{eqnarray*} |
184 |
|
|
185 |
|
\subsubsection{\label{introSection:resistanceTensorRegularArbitrary}\textbf{The Resistance Tensor for Arbitrary Shapes}} |
186 |
|
|
223 |
|
Here $R_{ij}$ is the distance vector between bead $i$ and bead $j$. |
224 |
|
A second order expression for element of different size was |
225 |
|
introduced by Rotne and Prager\cite{Rotne1969} and improved by |
226 |
< |
Garc\'{i}a de la Torre and Bloomfield\cite{Torre1977}, |
226 |
> |
Garc\'{i}a de la Torre and Bloomfield,\cite{Torre1977} |
227 |
|
\begin{equation} |
228 |
|
T_{ij} = \frac{1}{{8\pi \eta R_{ij} }}\left[ {\left( {I + |
229 |
|
\frac{{R_{ij} R_{ij}^T }}{{R_{ij}^2 }}} \right) + R\frac{{\sigma |
278 |
|
bead $i$ and origin $O$, the elements of resistance tensor at |
279 |
|
arbitrary origin $O$ can be written as |
280 |
|
\begin{eqnarray} |
281 |
< |
\Xi _{}^{tt} = \sum\limits_i {\sum\limits_j {C_{ij} } } \notag , \\ |
282 |
< |
\Xi _{}^{tr} = \Xi _{}^{rt} = \sum\limits_i {\sum\limits_j {U_i C_{ij} } } , \\ |
283 |
< |
\Xi _{}^{rr} = - \sum\limits_i {\sum\limits_j {U_i C_{ij} } } U_j. \notag \\ |
281 |
> |
\Xi _{}^{tt} & = & \sum\limits_i {\sum\limits_j {C_{ij} } } \notag , \\ |
282 |
> |
\Xi _{}^{tr} & = & \Xi _{}^{rt} = \sum\limits_i {\sum\limits_j {U_i C_{ij} } } , \\ |
283 |
> |
\Xi _{}^{rr} & = & - \sum\limits_i {\sum\limits_j {U_i C_{ij} } } U_j. \notag \\ |
284 |
|
\label{introEquation:ResistanceTensorArbitraryOrigin} |
285 |
|
\end{eqnarray} |
286 |
|
The resistance tensor depends on the origin to which they refer. The |
528 |
|
\end{center} |
529 |
|
\end{table} |
530 |
|
|
531 |
< |
Another set of calculation were performed to study the efficiency of |
531 |
> |
Another set of calculations were performed to study the efficiency of |
532 |
|
temperature control using different temperature coupling schemes. |
533 |
|
The starting configuration is cooled to 173~K and evolved using NVE, |
534 |
|
NVT, and Langevin dynamic with time step of 2 fs. |
538 |
|
Nos\'e-Hoover temperature scaling scheme with thermostat of 5 ps |
539 |
|
which gives reasonable tight coupling, while the blue one from |
540 |
|
Langevin dynamics with viscosity of 0.1 poise demonstrates a faster |
541 |
< |
scaling to the desire temperature. In extremely lower friction |
542 |
< |
regime (when $ \eta \approx 0$), Langevin dynamics becomes normal |
541 |
> |
scaling to the desire temperature. When $ \eta = 0$, Langevin dynamics becomes normal |
542 |
|
NVE (see orange curve in Fig.~\ref{langevin:temperature}) which |
543 |
|
loses the temperature control ability. |
544 |
|
|
564 |
|
made of 2266 small identical beads with size of 0.3 \AA on the |
565 |
|
surface. Applying the procedure described in |
566 |
|
Sec.~\ref{introEquation:ResistanceTensorArbitraryOrigin}, we |
567 |
< |
identified the center of resistance at $(0\AA, 0.7482\AA, |
568 |
< |
-0.1988\AA)$, as well as the resistance tensor, |
567 |
> |
identified the center of resistance at (0 $\rm{\AA}$, 0.7482 $\rm{\AA}$, |
568 |
> |
-0.1988 $\rm{\AA}$), as well as the resistance tensor, |
569 |
|
\[ |
570 |
|
\left( {\begin{array}{*{20}c} |
571 |
|
0.9261 & 0 & 0&0&0.08585&0.2057\\ |
572 |
|
0& 0.9270&-0.007063& 0.08585&0&0\\ |
573 |
< |
0&0.007063&0.7494&0.2057&0&0\\ |
574 |
< |
0&0.0858&0.2057& 58.64& 0&-8.5736\\ |
573 |
> |
0&-0.007063&0.7494&0.2057&0&0\\ |
574 |
> |
0&0.0858&0.2057& 58.64& 0&0\\ |
575 |
|
0.08585&0&0&0&48.30&3.219&\\ |
576 |
|
0.2057&0&0&0&3.219&10.7373\\ |
577 |
|
\end{array}} \right). |
578 |
|
\] |
579 |
< |
%\[ |
581 |
< |
%\left( {\begin{array}{*{20}c} |
582 |
< |
%0.9261 & 1.310e-14 & -7.292e-15&5.067e-14&0.08585&0.2057\\ |
583 |
< |
%3.968e-14& 0.9270&-0.007063& 0.08585&6.764e-14&4.846e-14\\ |
584 |
< |
%-6.561e-16&-0.007063&0.7494&0.2057&4.846e-14&1.5036e-14\\ |
585 |
< |
%5.067e-14&0.0858&0.2057& 58.64& 8.563e-13&-8.5736\\ |
586 |
< |
%0.08585&6.764e-14&4.846e-14&1.555e-12&48.30&3.219&\\ |
587 |
< |
%0.2057&4.846e-14&1.5036e-14&-3.904e-13&3.219&10.7373\\ |
588 |
< |
%\end{array}} \right). |
589 |
< |
%\] |
579 |
> |
where the units for translational, translation-rotation coupling and rotational tensors are $\frac{kcal \cdot fs}{mol \cdot \rm{\AA}^2}$, $\frac{kcal \cdot fs}{mol \cdot \rm{\AA} \cdot rad}$ and $\frac{kcal \cdot fs}{mol \cdot rad^2}$ respectively. |
580 |
|
|
581 |
|
Curves of the velocity auto-correlation functions in |
582 |
|
Fig.~\ref{langevin:vacf} were shown to match each other very well. |
583 |
|
However, because of the stochastic nature, simulation using Langevin |
584 |
|
dynamics was shown to decay slightly faster than MD. In order to |
585 |
|
study the rotational motion of the molecules, we also calculated the |
586 |
< |
auto- correlation function of the principle axis of the second GB |
586 |
> |
auto-correlation function of the principle axis of the second GB |
587 |
|
particle, $u$. The discrepancy shown in Fig.~\ref{langevin:uacf} was |
588 |
< |
probably due to the reason that the viscosity using in the |
599 |
< |
simulations only partially preserved the dynamics of the system. |
588 |
> |
probably due to the reason that we used the experimental viscosity directly instead of calculating bulk viscosity from simulation. |
589 |
|
|
590 |
|
\begin{figure} |
591 |
|
\centering |
606 |
|
\centering |
607 |
|
\includegraphics[width=\linewidth]{vacf.eps} |
608 |
|
\caption[Plots of Velocity Auto-correlation Functions]{Velocity |
609 |
< |
auto-correlation functions of NVE (explicit solvent) in blue) and |
609 |
> |
auto-correlation functions of NVE (explicit solvent) in blue and |
610 |
|
Langevin dynamics (implicit solvent) in red.} \label{langevin:vacf} |
611 |
|
\end{figure} |
612 |
|
|