48 |
|
\section{Introduction} |
49 |
|
|
50 |
|
%applications of langevin dynamics |
51 |
< |
As alternative to Newtonian dynamics, Langevin dynamics, which |
52 |
< |
mimics a simple heat bath with stochastic and dissipative forces, |
53 |
< |
has been applied in a variety of studies. The stochastic treatment |
54 |
< |
of the solvent enables us to carry out substantially longer time |
55 |
< |
simulations. Implicit solvent Langevin dynamics simulations of |
56 |
< |
met-enkephalin not only outperform explicit solvent simulations for |
57 |
< |
computational efficiency, but also agrees very well with explicit |
58 |
< |
solvent simulations for dynamical properties.\cite{Shen2002} |
59 |
< |
Recently, applying Langevin dynamics with the UNRES model, Liow and |
60 |
< |
his coworkers suggest that protein folding pathways can be possibly |
61 |
< |
explored within a reasonable amount of time.\cite{Liwo2005} The |
62 |
< |
stochastic nature of the Langevin dynamics also enhances the |
63 |
< |
sampling of the system and increases the probability of crossing |
64 |
< |
energy barriers.\cite{Banerjee2004, Cui2003} Combining Langevin |
65 |
< |
dynamics with Kramers's theory, Klimov and Thirumalai identified |
66 |
< |
free-energy barriers by studying the viscosity dependence of the |
67 |
< |
protein folding rates.\cite{Klimov1997} In order to account for |
68 |
< |
solvent induced interactions missing from implicit solvent model, |
69 |
< |
Kaya incorporated desolvation free energy barrier into implicit |
70 |
< |
coarse-grained solvent model in protein folding/unfolding studies |
71 |
< |
and discovered a higher free energy barrier between the native and |
72 |
< |
denatured states. Because of its stability against noise, Langevin |
73 |
< |
dynamics is very suitable for studying remagnetization processes in |
74 |
< |
various systems.\cite{Palacios1998,Berkov2002,Denisov2003} For |
51 |
> |
Langevin dynamics, which mimics a simple heat bath with stochastic and |
52 |
> |
dissipative forces, has been applied in a variety of situations as an |
53 |
> |
alternative to molecular dynamics with explicit solvent molecules. |
54 |
> |
The stochastic treatment of the solvent allows the use of simulations |
55 |
> |
with substantially longer time and length scales. In general, the |
56 |
> |
dynamic and structural properties obtained from Langevin simulations |
57 |
> |
agree quite well with similar properties obtained from explicit |
58 |
> |
solvent simulations. |
59 |
> |
|
60 |
> |
Recent examples of the usefulness of Langevin simulations include a |
61 |
> |
study of met-enkephalin in which Langevin simulations predicted |
62 |
> |
dynamical properties that were largely in agreement with explicit |
63 |
> |
solvent simulations.\cite{Shen2002} By applying Langevin dynamics with |
64 |
> |
the UNRES model, Liow and his coworkers suggest that protein folding |
65 |
> |
pathways can be explored within a reasonable amount of |
66 |
> |
time.\cite{Liwo2005} |
67 |
> |
|
68 |
> |
The stochastic nature of Langevin dynamics also enhances the sampling |
69 |
> |
of the system and increases the probability of crossing energy |
70 |
> |
barriers.\cite{Cui2003,Banerjee2004} Combining Langevin dynamics with |
71 |
> |
Kramers's theory, Klimov and Thirumalai identified free-energy |
72 |
> |
barriers by studying the viscosity dependence of the protein folding |
73 |
> |
rates.\cite{Klimov1997} In order to account for solvent induced |
74 |
> |
interactions missing from the implicit solvent model, Kaya |
75 |
> |
incorporated a desolvation free energy barrier into protein |
76 |
> |
folding/unfolding studies and discovered a higher free energy barrier |
77 |
> |
between the native and denatured states.\cite{XXX} |
78 |
> |
|
79 |
> |
Because of its stability against noise, Langevin dynamics has also |
80 |
> |
proven useful for studying remagnetization processes in various |
81 |
> |
systems.\cite{Palacios1998,Berkov2002,Denisov2003} [Check: For |
82 |
|
instance, the oscillation power spectrum of nanoparticles from |
83 |
< |
Langevin dynamics simulation has the same peak frequencies for |
84 |
< |
different wave vectors, which recovers the property of magnetic |
85 |
< |
excitations in small finite structures.\cite{Berkov2005a} |
83 |
> |
Langevin dynamics has the same peak frequencies for different wave |
84 |
> |
vectors, which recovers the property of magnetic excitations in small |
85 |
> |
finite structures.\cite{Berkov2005a}] |
86 |
|
|
87 |
< |
%review rigid body dynamics |
88 |
< |
Rigid bodies are frequently involved in the modeling of different |
89 |
< |
areas, from engineering, physics, to chemistry. For example, |
90 |
< |
missiles and vehicle are usually modeled by rigid bodies. The |
91 |
< |
movement of the objects in 3D gaming engine or other physics |
92 |
< |
simulator is governed by the rigid body dynamics. In molecular |
93 |
< |
simulation, rigid body is used to simplify the model in |
94 |
< |
protein-protein docking study{\cite{Gray2003}}. |
87 |
> |
In typical LD simulations, the friction and random forces on |
88 |
> |
individual atoms are taken from the Stokes-Einstein hydrodynamic |
89 |
> |
approximation, |
90 |
> |
\begin{eqnarray} |
91 |
> |
m \dot{v}(t) & = & -\nabla U(x) - \xi m v(t) + R(t) \\ |
92 |
> |
\langle R(t) \rangle & = & 0 \\ |
93 |
> |
\langle R(t) R(t') \rangle & = & 2 k_B T \xi m \delta(t - t') |
94 |
> |
\end{eqnarray} |
95 |
> |
where $\xi \approx 6 \pi \eta a$. Here $\eta$ is the viscosity of the |
96 |
> |
implicit solvent, and $a$ is the hydrodynamic radius of the atom. |
97 |
|
|
98 |
< |
It is very important to develop stable and efficient methods to |
99 |
< |
integrate the equations of motion for orientational degrees of |
100 |
< |
freedom. Euler angles are the natural choice to describe the |
101 |
< |
rotational degrees of freedom. However, due to $\frac {1}{sin |
102 |
< |
\theta}$ singularities, the numerical integration of corresponding |
103 |
< |
equations of these motion is very inefficient and inaccurate. |
104 |
< |
Although an alternative integrator using multiple sets of Euler |
105 |
< |
angles can overcome this difficulty\cite{Barojas1973}, the |
106 |
< |
computational penalty and the loss of angular momentum conservation |
107 |
< |
still remain. A singularity-free representation utilizing |
108 |
< |
quaternions was developed by Evans in 1977.\cite{Evans1977} |
109 |
< |
Unfortunately, this approach used a nonseparable Hamiltonian |
110 |
< |
resulting from the quaternion representation, which prevented the |
111 |
< |
symplectic algorithm from being utilized. Another different approach |
112 |
< |
is to apply holonomic constraints to the atoms belonging to the |
104 |
< |
rigid body. Each atom moves independently under the normal forces |
105 |
< |
deriving from potential energy and constraint forces which are used |
106 |
< |
to guarantee the rigidness. However, due to their iterative nature, |
107 |
< |
the SHAKE and Rattle algorithms also converge very slowly when the |
108 |
< |
number of constraints increases.\cite{Ryckaert1977, Andersen1983} |
98 |
> |
The use of rigid substructures,\cite{???} |
99 |
> |
coarse-graining,\cite{Ayton,Sun,Zannoni} and ellipsoidal |
100 |
> |
representations of protein side chains~\cite{Schulten} has made the |
101 |
> |
use of the Stokes-Einstein approximation problematic. A rigid |
102 |
> |
substructure moves as a single unit with orientational as well as |
103 |
> |
translational degrees of freedom. This requires a more general |
104 |
> |
treatment of the hydrodynamics than the spherical approximation |
105 |
> |
provides. The atoms involved in a rigid or coarse-grained structure |
106 |
> |
should properly have solvent-mediated interactions with each |
107 |
> |
other. The theory of interactions {\it between} bodies moving through |
108 |
> |
a fluid has been developed over the past century and has been applied |
109 |
> |
to simulations of Brownian |
110 |
> |
motion.\cite{MarshallNewton,GarciaDeLaTorre} There a need to have a |
111 |
> |
more thorough treatment of hydrodynamics included in the library of |
112 |
> |
methods available for performing Langevin simulations. |
113 |
|
|
114 |
< |
A break-through in geometric literature suggests that, in order to |
114 |
> |
\subsection{Rigid Body Dynamics} |
115 |
> |
Rigid bodies are frequently involved in the modeling of large |
116 |
> |
collections of particles that move as a single unit. In molecular |
117 |
> |
simulations, rigid bodies have been used to simplify protein-protein |
118 |
> |
docking,\cite{Gray2003} and lipid bilayer simulations.\cite{Sun2008} |
119 |
> |
Many of the water models in common use are also rigid-body |
120 |
> |
models,\cite{TIPs,SPC/E} although they are typically evolved using |
121 |
> |
constraints rather than rigid body equations of motion. |
122 |
> |
|
123 |
> |
Euler angles are a natural choice to describe the rotational |
124 |
> |
degrees of freedom. However, due to $1 \over \sin \theta$ |
125 |
> |
singularities, the numerical integration of corresponding equations of |
126 |
> |
these motion can become inaccurate (and inefficient). Although an |
127 |
> |
alternative integrator using multiple sets of Euler angles can |
128 |
> |
overcome this problem,\cite{Barojas1973} the computational penalty and |
129 |
> |
the loss of angular momentum conservation remain. A singularity-free |
130 |
> |
representation utilizing quaternions was developed by Evans in |
131 |
> |
1977.\cite{Evans1977} Unfortunately, this approach uses a nonseparable |
132 |
> |
Hamiltonian resulting from the quaternion representation, which |
133 |
> |
prevented symplectic algorithms from being utilized until very |
134 |
> |
recently.\cite{Miller2002} Another approach is the application of |
135 |
> |
holonomic constraints to the atoms belonging to the rigid body. Each |
136 |
> |
atom moves independently under the normal forces deriving from |
137 |
> |
potential energy and constraint forces which are used to guarantee the |
138 |
> |
rigidness. However, due to their iterative nature, the SHAKE and |
139 |
> |
Rattle algorithms also converge very slowly when the number of |
140 |
> |
constraints increases.\cite{Ryckaert1977,Andersen1983} |
141 |
> |
|
142 |
> |
A breakthrough in geometric literature suggests that, in order to |
143 |
|
develop a long-term integration scheme, one should preserve the |
144 |
|
symplectic structure of the propagator. By introducing a conjugate |
145 |
|
momentum to the rotation matrix $Q$ and re-formulating Hamiltonian's |
146 |
< |
equation, a symplectic integrator, RSHAKE\cite{Kol1997}, was |
147 |
< |
proposed to evolve the Hamiltonian system in a constraint manifold |
148 |
< |
by iteratively satisfying the orthogonality constraint $Q^T Q = 1$. |
149 |
< |
An alternative method using the quaternion representation was |
150 |
< |
developed by Omelyan.\cite{Omelyan1998} However, both of these |
151 |
< |
methods are iterative and inefficient. In this section, we descibe a |
152 |
< |
symplectic Lie-Poisson integrator for rigid bodies developed by |
153 |
< |
Dullweber and his coworkers\cite{Dullweber1997} in depth. |
146 |
> |
equation, a symplectic integrator, RSHAKE,\cite{Kol1997} was proposed |
147 |
> |
to evolve the Hamiltonian system in a constraint manifold by |
148 |
> |
iteratively satisfying the orthogonality constraint $Q^T Q = 1$. An |
149 |
> |
alternative method using the quaternion representation was developed |
150 |
> |
by Omelyan.\cite{Omelyan1998} However, both of these methods are |
151 |
> |
iterative and suffer from some related inefficiencies. A symplectic |
152 |
> |
Lie-Poisson integrator for rigid bodies developed by Dullweber {\it et |
153 |
> |
al.}\cite{Dullweber1997} gets around most of the limitations mentioned |
154 |
> |
above and has become the basis for our Langevin integrator. |
155 |
|
|
156 |
< |
%review langevin/browninan dynamics for arbitrarily shaped rigid body |
157 |
< |
Combining Langevin or Brownian dynamics with rigid body dynamics, |
158 |
< |
one can study slow processes in biomolecular systems. Modeling DNA |
159 |
< |
as a chain of rigid beads, which are subject to harmonic potentials |
160 |
< |
as well as excluded volume potentials, Mielke and his coworkers |
161 |
< |
discovered rapid superhelical stress generations from the stochastic |
162 |
< |
simulation of twin supercoiling DNA with response to induced |
163 |
< |
torques.\cite{Mielke2004} Membrane fusion is another key biological |
164 |
< |
process which controls a variety of physiological functions, such as |
165 |
< |
release of neurotransmitters \textit{etc}. A typical fusion event |
166 |
< |
happens on the time scale of a millisecond, which is impractical to |
167 |
< |
study using atomistic models with newtonian mechanics. With the help |
168 |
< |
of coarse-grained rigid body model and stochastic dynamics, the |
169 |
< |
fusion pathways were explored by many |
170 |
< |
researchers.\cite{Noguchi2001,Noguchi2002,Shillcock2005} Due to the |
171 |
< |
difficulty of numerical integration of anisotropic rotation, most of |
172 |
< |
the rigid body models are simply modeled using spheres, cylinders, |
173 |
< |
ellipsoids or other regular shapes in stochastic simulations. In an |
174 |
< |
effort to account for the diffusion anisotropy of arbitrary |
175 |
< |
particles, Fernandes and de la Torre improved the original Brownian |
176 |
< |
dynamics simulation algorithm\cite{Ermak1978,Allison1991} by |
177 |
< |
incorporating a generalized $6\times6$ diffusion tensor and |
178 |
< |
introducing a simple rotation evolution scheme consisting of three |
179 |
< |
consecutive rotations.\cite{Fernandes2002} Unfortunately, unexpected |
180 |
< |
errors and biases are introduced into the system due to the |
156 |
> |
|
157 |
> |
\subsection{The Hydrodynamic tensor and Brownian dynamics} |
158 |
> |
Combining Brownian dynamics with rigid substructures, one can study |
159 |
> |
slow processes in biomolecular systems. Modeling DNA as a chain of |
160 |
> |
beads which are subject to harmonic potentials as well as excluded |
161 |
> |
volume potentials, Mielke and his coworkers discovered rapid |
162 |
> |
superhelical stress generations from the stochastic simulation of twin |
163 |
> |
supercoiling DNA with response to induced torques.\cite{Mielke2004} |
164 |
> |
Membrane fusion is another key biological process which controls a |
165 |
> |
variety of physiological functions, such as release of |
166 |
> |
neurotransmitters \textit{etc}. A typical fusion event happens on the |
167 |
> |
time scale of a millisecond, which is impractical to study using |
168 |
> |
atomistic models with newtonian mechanics. With the help of |
169 |
> |
coarse-grained rigid body model and stochastic dynamics, the fusion |
170 |
> |
pathways were explored by Noguchi and others.\cite{Noguchi2001,Noguchi2002,Shillcock2005} |
171 |
> |
|
172 |
> |
Due to the difficulty of numerically integrating anisotropic |
173 |
> |
rotational motion, most of the coarse-grained rigid body models are |
174 |
> |
treated as spheres, cylinders, ellipsoids or other regular shapes in |
175 |
> |
stochastic simulations. In an effort to account for the diffusion |
176 |
> |
anisotropy of arbitrarily-shaped particles, Fernandes and Garc\'{i}a |
177 |
> |
de la Torre improved the original Brownian dynamics simulation |
178 |
> |
algorithm~\cite{Ermak1978,Allison1991} by incorporating a generalized |
179 |
> |
$6\times6$ diffusion tensor and introducing a rotational evolution |
180 |
> |
scheme consisting of three consecutive rotations.\cite{Fernandes2002} |
181 |
> |
Unfortunately, biases are introduced into the system due to the |
182 |
|
arbitrary order of applying the noncommuting rotation |
183 |
|
operators.\cite{Beard2003} Based on the observation the momentum |
184 |
|
relaxation time is much less than the time step, one may ignore the |
185 |
< |
inertia in Brownian dynamics. However, the assumption of zero |
186 |
< |
average acceleration is not always true for cooperative motion which |
187 |
< |
is common in protein motion. An inertial Brownian dynamics (IBD) was |
188 |
< |
proposed to address this issue by adding an inertial correction |
185 |
> |
inertia in Brownian dynamics. However, the assumption of zero average |
186 |
> |
acceleration is not always true for cooperative motion which is common |
187 |
> |
in proteins. An inertial Brownian dynamics (IBD) was proposed to |
188 |
> |
address this issue by adding an inertial correction |
189 |
|
term.\cite{Beard2000} As a complement to IBD which has a lower bound |
190 |
|
in time step because of the inertial relaxation time, long-time-step |
191 |
|
inertial dynamics (LTID) can be used to investigate the inertial |
192 |
|
behavior of the polymer segments in low friction |
193 |
|
regime.\cite{Beard2000} LTID can also deal with the rotational |
194 |
|
dynamics for nonskew bodies without translation-rotation coupling by |
195 |
< |
separating the translation and rotation motion and taking advantage |
196 |
< |
of the analytical solution of hydrodynamics properties. However, |
197 |
< |
typical nonskew bodies like cylinders and ellipsoids are inadequate |
198 |
< |
to represent most complex macromolecule assemblies. These intricate |
165 |
< |
molecules have been represented by a set of beads and their |
166 |
< |
hydrodynamic properties can be calculated using variants on the |
167 |
< |
standard hydrodynamic interaction tensors. |
195 |
> |
separating the translation and rotation motion and taking advantage of |
196 |
> |
the analytical solution of hydrodynamics properties. However, typical |
197 |
> |
nonskew bodies like cylinders and ellipsoids are inadequate to |
198 |
> |
represent most complex macromolecular assemblies. |
199 |
|
|
200 |
|
The goal of the present work is to develop a Langevin dynamics |
201 |
|
algorithm for arbitrary-shaped rigid particles by integrating the |
202 |
< |
accurate estimation of friction tensor from hydrodynamics theory |
203 |
< |
into the sophisticated rigid body dynamics algorithms. |
202 |
> |
accurate estimation of friction tensor from hydrodynamics theory into |
203 |
> |
a symplectic rigid body dynamics propagator. In the sections below, |
204 |
> |
we review some of the theory of hydrodynamic tensors developed for |
205 |
> |
Brownian simulations of rigid multi-particle systems, we then present |
206 |
> |
our integration method for a set of generalized Langevin equations of |
207 |
> |
motion, and we compare the behavior of the new Langevin integrator to |
208 |
> |
dynamical quantities obtained via explicit solvent molecular dynamics. |
209 |
|
|
210 |
< |
\subsection{\label{introSection:frictionTensor}Friction Tensor} |
211 |
< |
Theoretically, the friction kernel can be determined using the |
210 |
> |
\subsection{\label{introSection:frictionTensor}The Friction Tensor} |
211 |
> |
Theoretically, a complete friction kernel can be determined using the |
212 |
|
velocity autocorrelation function. However, this approach becomes |
213 |
< |
impractical when the system becomes more and more complicated. |
214 |
< |
Instead, various approaches based on hydrodynamics have been |
215 |
< |
developed to calculate the friction coefficients. In general, the |
216 |
< |
friction tensor $\Xi$ is a $6\times 6$ matrix given by |
217 |
< |
\[ |
213 |
> |
impractical when the solute becomes complex. Instead, various |
214 |
> |
approaches based on hydrodynamics have been developed to calculate the |
215 |
> |
friction coefficients. In general, the friction tensor $\Xi$ is a |
216 |
> |
$6\times 6$ matrix given by |
217 |
> |
\begin{equation} |
218 |
|
\Xi = \left( {\begin{array}{*{20}c} |
219 |
|
{\Xi _{}^{tt} } & {\Xi _{}^{rt} } \\ |
220 |
|
{\Xi _{}^{tr} } & {\Xi _{}^{rr} } \\ |
221 |
|
\end{array}} \right). |
222 |
< |
\] |
223 |
< |
Here, $ {\Xi^{tt} }$ and $ {\Xi^{rr} }$ are $3 \times 3$ |
224 |
< |
translational friction tensor and rotational resistance (friction) |
225 |
< |
tensor respectively, while ${\Xi^{tr} }$ is translation-rotation |
226 |
< |
coupling tensor and $ {\Xi^{rt} }$ is rotation-translation coupling |
227 |
< |
tensor. When a particle moves in a fluid, it may experience friction |
228 |
< |
force or torque along the opposite direction of the velocity or |
229 |
< |
angular velocity, |
230 |
< |
\[ |
222 |
> |
\end{equation} |
223 |
> |
Here, $\Xi^{tt}$ and $\Xi^{rr}$ are $3 \times 3$ translational and |
224 |
> |
rotational resistance (friction) tensors respectively, while |
225 |
> |
$\Xi^{tr}$ is translation-rotation coupling tensor and $\Xi^{rt}$ is |
226 |
> |
rotation-translation coupling tensor. When a particle moves in a |
227 |
> |
fluid, it may experience friction force ($\mathbf{F}_f$) and torque |
228 |
> |
($\mathbf{\tau}_f$) in opposition to the directions of the velocity |
229 |
> |
($\mathbf{v}$) and body-fixed angular velocity ($\mathbf{\omega}$), |
230 |
> |
\begin{equation} |
231 |
|
\left( \begin{array}{l} |
232 |
< |
F_R \\ |
233 |
< |
\tau _R \\ |
232 |
> |
\mathbf{F}_f \\ |
233 |
> |
\mathbf{\tau}_f \\ |
234 |
|
\end{array} \right) = - \left( {\begin{array}{*{20}c} |
235 |
< |
{\Xi ^{tt} } & {\Xi ^{rt} } \\ |
236 |
< |
{\Xi ^{tr} } & {\Xi ^{rr} } \\ |
235 |
> |
\Xi ^{tt} & \Xi ^{rt} \\ |
236 |
> |
\Xi ^{tr} & \Xi ^{rr} \\ |
237 |
|
\end{array}} \right)\left( \begin{array}{l} |
238 |
< |
v \\ |
239 |
< |
w \\ |
240 |
< |
\end{array} \right) |
241 |
< |
\] |
206 |
< |
where $F_r$ is the friction force and $\tau _R$ is the friction |
207 |
< |
torque. |
238 |
> |
\mathbf{v} \\ |
239 |
> |
\mathbf{\omega} \\ |
240 |
> |
\end{array} \right). |
241 |
> |
\end{equation} |
242 |
|
|
243 |
|
\subsubsection{\label{introSection:resistanceTensorRegular}\textbf{The Resistance Tensor for Regular Shapes}} |
244 |
< |
|
245 |
< |
For a spherical particle with slip boundary conditions, the |
246 |
< |
translational and rotational friction constant can be calculated |
247 |
< |
from Stoke's law, |
248 |
< |
\[ |
215 |
< |
\Xi ^{tt} = \left( {\begin{array}{*{20}c} |
244 |
> |
For a spherical particle under ``stick'' boundary conditions, the |
245 |
> |
translational and rotational friction tensors can be calculated from |
246 |
> |
Stoke's law, |
247 |
> |
\begin{equation} |
248 |
> |
\Xi^{tt} = \left( {\begin{array}{*{20}c} |
249 |
|
{6\pi \eta R} & 0 & 0 \\ |
250 |
|
0 & {6\pi \eta R} & 0 \\ |
251 |
|
0 & 0 & {6\pi \eta R} \\ |
252 |
|
\end{array}} \right) |
253 |
< |
\] |
253 |
> |
\end{equation} |
254 |
|
and |
255 |
< |
\[ |
255 |
> |
\begin{equation} |
256 |
|
\Xi ^{rr} = \left( {\begin{array}{*{20}c} |
257 |
|
{8\pi \eta R^3 } & 0 & 0 \\ |
258 |
|
0 & {8\pi \eta R^3 } & 0 \\ |
259 |
|
0 & 0 & {8\pi \eta R^3 } \\ |
260 |
|
\end{array}} \right) |
261 |
< |
\] |
261 |
> |
\end{equation} |
262 |
|
where $\eta$ is the viscosity of the solvent and $R$ is the |
263 |
|
hydrodynamic radius. |
264 |
|
|
265 |
|
Other non-spherical shapes, such as cylinders and ellipsoids, are |
266 |
< |
widely used as references for developing new hydrodynamics theory, |
266 |
> |
widely used as references for developing new hydrodynamics theories, |
267 |
|
because their properties can be calculated exactly. In 1936, Perrin |
268 |
|
extended Stokes's law to general ellipsoids, also called a triaxial |
269 |
|
ellipsoid, which is given in Cartesian coordinates |
270 |
< |
by\cite{Perrin1934, Perrin1936} |
271 |
< |
\[ |
272 |
< |
\frac{{x^2 }}{{a^2 }} + \frac{{y^2 }}{{b^2 }} + \frac{{z^2 }}{{c^2 |
273 |
< |
}} = 1 |
274 |
< |
\] |
275 |
< |
where the semi-axes are of lengths $a$, $b$, and $c$. Unfortunately, |
276 |
< |
due to the complexity of the elliptic integral, only the ellipsoid |
277 |
< |
with the restriction of two axes being equal, \textit{i.e.} |
278 |
< |
prolate($ a \ge b = c$) and oblate ($ a < b = c $), can be solved |
279 |
< |
exactly. Introducing an elliptic integral parameter $S$ for prolate |
280 |
< |
ellipsoids : |
281 |
< |
\[ |
282 |
< |
S = \frac{2}{{\sqrt {a^2 - b^2 } }}\ln \frac{{a + \sqrt {a^2 - b^2 |
283 |
< |
} }}{b}, |
284 |
< |
\] |
285 |
< |
and oblate ellipsoids: |
286 |
< |
\[ |
287 |
< |
S = \frac{2}{{\sqrt {b^2 - a^2 } }}arctg\frac{{\sqrt {b^2 - a^2 } |
288 |
< |
}}{a}, |
289 |
< |
\] |
290 |
< |
one can write down the translational and rotational resistance |
258 |
< |
tensors |
259 |
< |
\begin{eqnarray*} |
260 |
< |
\Xi _a^{tt} & = & 16\pi \eta \frac{{a^2 - b^2 }}{{(2a^2 - b^2 )S - 2a}}. \\ |
261 |
< |
\Xi _b^{tt} & = & \Xi _c^{tt} = 32\pi \eta \frac{{a^2 - b^2 }}{{(2a^2 - 3b^2 )S + |
262 |
< |
2a}}, |
263 |
< |
\end{eqnarray*} |
264 |
< |
and |
265 |
< |
\begin{eqnarray*} |
266 |
< |
\Xi _a^{rr} & = & \frac{{32\pi }}{3}\eta \frac{{(a^2 - b^2 )b^2 }}{{2a - b^2 S}}, \\ |
267 |
< |
\Xi _b^{rr} & = & \Xi _c^{rr} = \frac{{32\pi }}{3}\eta \frac{{(a^4 - b^4 )}}{{(2a^2 - b^2 )S - 2a}}. |
270 |
> |
by\cite{Perrin1934,Perrin1936} |
271 |
> |
\begin{equation} |
272 |
> |
\frac{x^2 }{a^2} + \frac{y^2}{b^2} + \frac{z^2 }{c^2} = 1 |
273 |
> |
\end{equation} |
274 |
> |
where the semi-axes are of lengths $a$, $b$, and $c$. Due to the |
275 |
> |
complexity of the elliptic integral, only uniaxial ellipsoids, |
276 |
> |
{\it i.e.} prolate ($ a \ge b = c$) and oblate ($ a < b = c $), can |
277 |
> |
be solved exactly. Introducing an elliptic integral parameter $S$ for |
278 |
> |
prolate ellipsoids : |
279 |
> |
\begin{equation} |
280 |
> |
S = \frac{2}{\sqrt{a^2 - b^2}} \ln \frac{a + \sqrt{a^2 - b^2}}{b}, |
281 |
> |
\end{equation} |
282 |
> |
and oblate ellipsoids: |
283 |
> |
\begin{equation} |
284 |
> |
S = \frac{2}{\sqrt {b^2 - a^2 }} \arctan \frac{\sqrt {b^2 - a^2}}{a}, |
285 |
> |
\end{equation} |
286 |
> |
one can write down the translational and rotational resistance |
287 |
> |
tensors for oblate, |
288 |
> |
\begin{eqnarray*} |
289 |
> |
\Xi_a^{tt} & = & 16\pi \eta \frac{a^2 - b^2}{(2a^2 - b^2 )S - 2a}. \\ |
290 |
> |
\Xi_b^{tt} = \Xi_c^{tt} & = & 32\pi \eta \frac{a^2 - b^2 }{(2a^2 - 3b^2 )S + 2a}, |
291 |
|
\end{eqnarray*} |
292 |
+ |
and prolate, |
293 |
+ |
\begin{eqnarray*} |
294 |
+ |
\Xi_a^{rr} & = & \frac{32\pi}{3} \eta \frac{(a^2 - b^2 )b^2}{2a - b^2 S}, \\ |
295 |
+ |
\Xi_b^{rr} = \Xi_c^{rr} & = & \frac{32\pi}{3} \eta \frac{(a^4 - b^4)}{(2a^2 - b^2 )S - 2a} |
296 |
+ |
\end{eqnarray*} |
297 |
+ |
ellipsoids. For both spherical and ellipsoidal particles, the |
298 |
+ |
translation-rotation and rotation-translation coupling tensors are |
299 |
+ |
zero. |
300 |
|
|
301 |
|
\subsubsection{\label{introSection:resistanceTensorRegularArbitrary}\textbf{The Resistance Tensor for Arbitrary Shapes}} |
302 |
|
|
304 |
|
analytical solution for the friction tensor for arbitrarily shaped |
305 |
|
rigid molecules. The ellipsoid of revolution model and general |
306 |
|
triaxial ellipsoid model have been used to approximate the |
307 |
< |
hydrodynamic properties of rigid bodies. However, since the mapping |
308 |
< |
from all possible ellipsoidal spaces, $r$-space, to all possible |
309 |
< |
combination of rotational diffusion coefficients, $D$-space, is not |
310 |
< |
unique\cite{Wegener1979} as well as the intrinsic coupling between |
311 |
< |
translational and rotational motion of rigid bodies, general |
312 |
< |
ellipsoids are not always suitable for modeling arbitrarily shaped |
313 |
< |
rigid molecules. A number of studies have been devoted to |
307 |
> |
hydrodynamic properties of rigid bodies. However, the mapping from all |
308 |
> |
possible ellipsoidal spaces, $r$-space, to all possible combination of |
309 |
> |
rotational diffusion coefficients, $D$-space, is not |
310 |
> |
unique.\cite{Wegener1979} Additionally, because there is intrinsic |
311 |
> |
coupling between translational and rotational motion of rigid bodies, |
312 |
> |
general ellipsoids are not always suitable for modeling arbitrarily |
313 |
> |
shaped rigid molecules. A number of studies have been devoted to |
314 |
|
determining the friction tensor for irregularly shaped rigid bodies |
315 |
< |
using more advanced methods where the molecule of interest was |
316 |
< |
modeled by a combinations of spheres\cite{Carrasco1999} and the |
317 |
< |
hydrodynamics properties of the molecule can be calculated using the |
318 |
< |
hydrodynamic interaction tensor. Let us consider a rigid assembly of |
319 |
< |
$N$ beads immersed in a continuous medium. Due to hydrodynamic |
320 |
< |
interaction, the ``net'' velocity of $i$th bead, $v'_i$ is different |
321 |
< |
than its unperturbed velocity $v_i$, |
315 |
> |
using more advanced methods where the molecule of interest was modeled |
316 |
> |
by a combinations of spheres\cite{Carrasco1999} and the hydrodynamics |
317 |
> |
properties of the molecule can be calculated using the hydrodynamic |
318 |
> |
interaction tensor. Let us consider a rigid assembly of $N$ beads |
319 |
> |
immersed in a continuous medium. Due to hydrodynamic interaction, the |
320 |
> |
``net'' velocity of $i$th bead, $v'_i$ is different than its |
321 |
> |
unperturbed velocity $v_i$, |
322 |
|
\[ |
323 |
|
v'_i = v_i - \sum\limits_{j \ne i} {T_{ij} F_j } |
324 |
|
\] |
1092 |
|
|
1093 |
|
\subsection{Composite sphero-ellipsoids} |
1094 |
|
Spherical heads perched on the ends of Gay-Berne ellipsoids have been |
1095 |
< |
used recently as models for lipid molecules.\cite{SunGezelter08,Ayton01} |
1096 |
< |
|
1095 |
> |
used recently as models for lipid |
1096 |
> |
molecules.\cite{SunGezelter08,Ayton01} |
1097 |
|
MORE DETAILS |
1098 |
|
|
1099 |
+ |
A reference system composed of a single lipid rigid body embedded in a |
1100 |
+ |
sea of 1929 solvent particles was created and run under standard |
1101 |
+ |
(microcanonical) molecular dynamics. The resulting viscosity of this |
1102 |
+ |
mixture was 0.349 centipoise (as estimated using |
1103 |
+ |
Eq. (\ref{eq:shear})). To calculate the hydrodynamic properties of |
1104 |
+ |
the lipid rigid body model, we created a rough shell (see |
1105 |
+ |
Fig.~\ref{fig:roughShell}), in which the lipid is represented as a |
1106 |
+ |
``shell'' made of 3550 identical beads (0.25 \AA\ in diameter) |
1107 |
+ |
distributed on the surface. Applying the procedure described in |
1108 |
+ |
Sec.~\ref{introEquation:ResistanceTensorArbitraryOrigin}, we |
1109 |
+ |
identified the center of resistance, ${\bf r} = $(0 \AA, 0 \AA, 1.46 |
1110 |
+ |
\AA). |
1111 |
|
|
1112 |
+ |
|
1113 |
|
\subsection{Summary} |
1114 |
|
According to our simulations, the langevin dynamics is a reliable |
1115 |
|
theory to apply to replace the explicit solvents, especially for the |
1116 |
|
translation properties. For large molecules, the rotation properties |
1117 |
|
are also mimiced reasonablly well. |
1118 |
+ |
|
1119 |
+ |
\begin{figure} |
1120 |
+ |
\centering |
1121 |
+ |
\includegraphics[width=\linewidth]{graph} |
1122 |
+ |
\caption[Mean squared displacements and orientational |
1123 |
+ |
correlation functions for each of the model rigid bodies.]{The |
1124 |
+ |
mean-squared displacements ($\langle r^2(t) \rangle$) and |
1125 |
+ |
orientational correlation functions ($C_2(t)$) for each of the model |
1126 |
+ |
rigid bodies studied. The circles are the results for microcanonical |
1127 |
+ |
simulations with explicit solvent molecules, while the other data sets |
1128 |
+ |
are results for Langevin dynamics using the different hydrodynamic |
1129 |
+ |
tensor approximations. The Perrin model for the ellipsoids is |
1130 |
+ |
considered the ``exact'' hydrodynamic behavior (this can also be said |
1131 |
+ |
for the translational motion of the dumbbell operating under the bead |
1132 |
+ |
model). In most cases, the various hydrodynamics models reproduce |
1133 |
+ |
each other quantitatively.} |
1134 |
+ |
\label{fig:results} |
1135 |
+ |
\end{figure} |
1136 |
|
|
1137 |
|
\begin{table*} |
1138 |
|
\begin{minipage}{\linewidth} |
1152 |
|
\cline{2-3} \cline{5-7} |
1153 |
|
model & $\eta$ (centipoise) & D & & Analytical & method & Hydrodynamics & simulation \\ |
1154 |
|
\hline |
1155 |
< |
sphere & 0.261 & ? & & 2.59 & exact & 2.59 & 2.56 \\ |
1155 |
> |
sphere & 0.279 & 3.06 & & 2.42 & exact & 2.42 & 2.33 \\ |
1156 |
|
ellipsoid & 0.255 & 2.44 & & 2.34 & exact & 2.34 & 2.37 \\ |
1157 |
|
& 0.255 & 2.44 & & 2.34 & rough shell & 2.36 & 2.28 \\ |
1158 |
< |
dumbbell & 0.322 & ? & & 1.57 & bead model & 1.57 & 1.57 \\ |
1159 |
< |
& 0.322 & ? & & 1.57 & rough shell & ? & ? \\ |
1158 |
> |
dumbbell & 0.308 & 2.06 & & 1.64 & bead model & 1.65 & 1.62 \\ |
1159 |
> |
& 0.308 & 2.06 & & 1.64 & rough shell & 1.59 & 1.62 \\ |
1160 |
|
banana & 0.298 & 1.53 & & & rough shell & 1.56 & 1.55 \\ |
1161 |
|
lipid & 0.349 & 0.96 & & & rough shell & 1.33 & 1.32 \\ |
1162 |
|
\end{tabular} |
1181 |
|
\cline{2-3} \cline{5-7} |
1182 |
|
model & $\eta$ (centipoise) & $\tau$ & & Perrin & method & Hydrodynamic & simulation \\ |
1183 |
|
\hline |
1184 |
< |
sphere & 0.261 & & & 9.06 & exact & 9.06 & 9.11 \\ |
1184 |
> |
sphere & 0.279 & & & 9.69 & exact & 9.69 & 9.64 \\ |
1185 |
|
ellipsoid & 0.255 & 46.7 & & 22.0 & exact & 22.0 & 22.2 \\ |
1186 |
|
& 0.255 & 46.7 & & 22.0 & rough shell & 22.6 & 22.2 \\ |
1187 |
< |
dumbbell & 0.322 & 14.0 & & & bead model & 52.3 & 52.8 \\ |
1188 |
< |
& 0.322 & 14.0 & & & rough shell & ? & ? \\ |
1187 |
> |
dumbbell & 0.308 & 14.1 & & & bead model & 50.0 & 50.1 \\ |
1188 |
> |
& 0.308 & 14.1 & & & rough shell & 41.5 & 41.3 \\ |
1189 |
|
banana & 0.298 & 63.8 & & & rough shell & 70.9 & 70.9 \\ |
1190 |
|
lipid & 0.349 & 78.0 & & & rough shell & 76.9 & 77.9 \\ |
1191 |
|
\hline |