| 281 |
|
(q,p,t)dq_1 } ...dq_f dp_1 } ...dp_f }}{{\int { \ldots \int {\rho |
| 282 |
|
(q,p,t)dq_1 } ...dq_f dp_1 } ...dp_f }}. |
| 283 |
|
\label{introEquation:ensembelAverage} |
| 284 |
– |
\end{equation} |
| 285 |
– |
|
| 286 |
– |
There are several different types of ensembles with different |
| 287 |
– |
statistical characteristics. As a function of macroscopic |
| 288 |
– |
parameters, such as temperature \textit{etc}, the partition function |
| 289 |
– |
can be used to describe the statistical properties of a system in |
| 290 |
– |
thermodynamic equilibrium. As an ensemble of systems, each of which |
| 291 |
– |
is known to be thermally isolated and conserve energy, the |
| 292 |
– |
Microcanonical ensemble (NVE) has a partition function like, |
| 293 |
– |
\begin{equation} |
| 294 |
– |
\Omega (N,V,E) = e^{\beta TS}. \label{introEquation:NVEPartition} |
| 295 |
– |
\end{equation} |
| 296 |
– |
A canonical ensemble (NVT) is an ensemble of systems, each of which |
| 297 |
– |
can share its energy with a large heat reservoir. The distribution |
| 298 |
– |
of the total energy amongst the possible dynamical states is given |
| 299 |
– |
by the partition function, |
| 300 |
– |
\begin{equation} |
| 301 |
– |
\Omega (N,V,T) = e^{ - \beta A}. |
| 302 |
– |
\label{introEquation:NVTPartition} |
| 303 |
– |
\end{equation} |
| 304 |
– |
Here, $A$ is the Helmholtz free energy which is defined as $ A = U - |
| 305 |
– |
TS$. Since most experiments are carried out under constant pressure |
| 306 |
– |
condition, the isothermal-isobaric ensemble (NPT) plays a very |
| 307 |
– |
important role in molecular simulations. The isothermal-isobaric |
| 308 |
– |
ensemble allow the system to exchange energy with a heat bath of |
| 309 |
– |
temperature $T$ and to change the volume as well. Its partition |
| 310 |
– |
function is given as |
| 311 |
– |
\begin{equation} |
| 312 |
– |
\Delta (N,P,T) = - e^{\beta G}. |
| 313 |
– |
\label{introEquation:NPTPartition} |
| 284 |
|
\end{equation} |
| 315 |
– |
Here, $G = U - TS + PV$, and $G$ is called Gibbs free energy. |
| 285 |
|
|
| 286 |
|
\subsection{\label{introSection:liouville}Liouville's theorem} |
| 287 |
|
|
| 505 |
|
\subsection{\label{introSection:exactFlow}Exact Propagator} |
| 506 |
|
|
| 507 |
|
Let $x(t)$ be the exact solution of the ODE |
| 508 |
< |
system,$\frac{{dx}}{{dt}} = f(x) \label{introEquation:ODE}$, we can |
| 509 |
< |
define its exact propagator(solution) $\varphi_\tau$ |
| 508 |
> |
system, |
| 509 |
> |
\begin{equation} |
| 510 |
> |
\frac{{dx}}{{dt}} = f(x), \label{introEquation:ODE} |
| 511 |
> |
\end{equation} we can |
| 512 |
> |
define its exact propagator $\varphi_\tau$: |
| 513 |
|
\[ x(t+\tau) |
| 514 |
|
=\varphi_\tau(x(t)) |
| 515 |
|
\] |
| 581 |
|
\] |
| 582 |
|
Using the chain rule, one may obtain, |
| 583 |
|
\[ |
| 584 |
< |
\sum\limits_i {\frac{{dG}}{{dx_i }}} f_i (x) = f \dot \nabla G, |
| 584 |
> |
\sum\limits_i {\frac{{dG}}{{dx_i }}} f_i (x) = f \cdot \nabla G, |
| 585 |
|
\] |
| 586 |
|
which is the condition for conserved quantities. For a canonical |
| 587 |
|
Hamiltonian system, the time evolution of an arbitrary smooth |
| 706 |
|
\end{align} |
| 707 |
|
where $F(t)$ is the force at time $t$. This integration scheme is |
| 708 |
|
known as \emph{velocity verlet} which is |
| 709 |
< |
symplectic(\ref{introEquation:SymplecticFlowComposition}), |
| 710 |
< |
time-reversible(\ref{introEquation:timeReversible}) and |
| 711 |
< |
volume-preserving (\ref{introEquation:volumePreserving}). These |
| 709 |
> |
symplectic(Eq.~\ref{introEquation:SymplecticFlowComposition}), |
| 710 |
> |
time-reversible(Eq.~\ref{introEquation:timeReversible}) and |
| 711 |
> |
volume-preserving (Eq.~\ref{introEquation:volumePreserving}). These |
| 712 |
|
geometric properties attribute to its long-time stability and its |
| 713 |
|
popularity in the community. However, the most commonly used |
| 714 |
|
velocity verlet integration scheme is written as below, |
| 750 |
|
|
| 751 |
|
The Baker-Campbell-Hausdorff formula can be used to determine the |
| 752 |
|
local error of a splitting method in terms of the commutator of the |
| 753 |
< |
operators(\ref{introEquation:exponentialOperator}) associated with |
| 753 |
> |
operators(Eq.~\ref{introEquation:exponentialOperator}) associated with |
| 754 |
|
the sub-propagator. For operators $hX$ and $hY$ which are associated |
| 755 |
|
with $\varphi_1(t)$ and $\varphi_2(t)$ respectively , we have |
| 756 |
|
\begin{equation} |
| 834 |
|
These three individual steps will be covered in the following |
| 835 |
|
sections. Sec.~\ref{introSec:initialSystemSettings} deals with the |
| 836 |
|
initialization of a simulation. Sec.~\ref{introSection:production} |
| 837 |
< |
will discuss issues of production runs. |
| 837 |
> |
discusses issues of production runs. |
| 838 |
|
Sec.~\ref{introSection:Analysis} provides the theoretical tools for |
| 839 |
|
analysis of trajectories. |
| 840 |
|
|
| 902 |
|
properties \textit{etc}, become independent of time. Strictly |
| 903 |
|
speaking, minimization and heating are not necessary, provided the |
| 904 |
|
equilibration process is long enough. However, these steps can serve |
| 905 |
< |
as a means to arrive at an equilibrated structure in an effective |
| 905 |
> |
as a mean to arrive at an equilibrated structure in an effective |
| 906 |
|
way. |
| 907 |
|
|
| 908 |
|
\subsection{\label{introSection:production}Production} |
| 972 |
|
V(r_{ij})= \frac{q_i q_j \textrm{erfc}(\alpha |
| 973 |
|
r_{ij})}{r_{ij}}-\lim_{r_{ij}\rightarrow |
| 974 |
|
R_\textrm{c}}\left\{\frac{q_iq_j \textrm{erfc}(\alpha |
| 975 |
< |
r_{ij})}{r_{ij}}\right\}. \label{introEquation:shiftedCoulomb} |
| 975 |
> |
r_{ij})}{r_{ij}}\right\}, \label{introEquation:shiftedCoulomb} |
| 976 |
|
\end{equation} |
| 977 |
|
where $\alpha$ is the convergence parameter. Due to the lack of |
| 978 |
|
inherent periodicity and rapid convergence,this method is extremely |
| 989 |
|
|
| 990 |
|
\subsection{\label{introSection:Analysis} Analysis} |
| 991 |
|
|
| 992 |
< |
Recently, advanced visualization technique have become applied to |
| 992 |
> |
Recently, advanced visualization techniques have been applied to |
| 993 |
|
monitor the motions of molecules. Although the dynamics of the |
| 994 |
|
system can be described qualitatively from animation, quantitative |
| 995 |
|
trajectory analysis is more useful. According to the principles of |
| 1059 |
|
\label{introEquation:timeCorrelationFunction} |
| 1060 |
|
\end{equation} |
| 1061 |
|
If $A$ and $B$ refer to same variable, this kind of correlation |
| 1062 |
< |
function is called an \emph{autocorrelation function}. One example |
| 1063 |
< |
of an auto correlation function is the velocity auto-correlation |
| 1062 |
> |
functions are called \emph{autocorrelation functions}. One example |
| 1063 |
> |
of auto correlation function is the velocity auto-correlation |
| 1064 |
|
function which is directly related to transport properties of |
| 1065 |
|
molecular liquids: |
| 1066 |
|
\[ |
| 1111 |
|
computational penalty and the loss of angular momentum conservation |
| 1112 |
|
still remain. A singularity-free representation utilizing |
| 1113 |
|
quaternions was developed by Evans in 1977\cite{Evans1977}. |
| 1114 |
< |
Unfortunately, this approach uses a nonseparable Hamiltonian |
| 1115 |
< |
resulting from the quaternion representation, which prevents the |
| 1114 |
> |
Unfortunately, this approach used a nonseparable Hamiltonian |
| 1115 |
> |
resulting from the quaternion representation, which prevented the |
| 1116 |
|
symplectic algorithm from being utilized. Another different approach |
| 1117 |
|
is to apply holonomic constraints to the atoms belonging to the |
| 1118 |
|
rigid body. Each atom moves independently under the normal forces |
| 1155 |
|
Q^T Q = 1, \label{introEquation:orthogonalConstraint} |
| 1156 |
|
\end{equation} |
| 1157 |
|
which is used to ensure the rotation matrix's unitarity. Using |
| 1158 |
< |
Equation (\ref{introEquation:motionHamiltonianCoordinate}, |
| 1159 |
< |
\ref{introEquation:motionHamiltonianMomentum}), one can write down |
| 1158 |
> |
Eq.~\ref{introEquation:motionHamiltonianCoordinate} and Eq.~ |
| 1159 |
> |
\ref{introEquation:motionHamiltonianMomentum}, one can write down |
| 1160 |
|
the equations of motion, |
| 1161 |
|
\begin{eqnarray} |
| 1162 |
|
\frac{{dq}}{{dt}} & = & \frac{p}{m}, \label{introEquation:RBMotionPosition}\\ |
| 1188 |
|
\[ |
| 1189 |
|
\tilde Q = Q,\tilde P = \frac{1}{2}\left( {P - QP^T Q} \right), |
| 1190 |
|
\] |
| 1191 |
< |
the mechanical system subject to a holonomic constraint manifold $M$ |
| 1191 |
> |
the mechanical system subjected to a holonomic constraint manifold $M$ |
| 1192 |
|
can be re-formulated as a Hamiltonian system on the cotangent space |
| 1193 |
|
\[ |
| 1194 |
|
T^* SO(3) = \left\{ {(\tilde Q,\tilde P):\tilde Q^T \tilde Q = |
| 1377 |
|
Splitting for Rigid Body} |
| 1378 |
|
|
| 1379 |
|
The Hamiltonian of rigid body can be separated in terms of kinetic |
| 1380 |
< |
energy and potential energy,$H = T(p,\pi ) + V(q,Q)$. The equations |
| 1380 |
> |
energy and potential energy, $H = T(p,\pi ) + V(q,Q)$. The equations |
| 1381 |
|
of motion corresponding to potential energy and kinetic energy are |
| 1382 |
< |
listed in the below table, |
| 1382 |
> |
listed in Table~\ref{introTable:rbEquations} |
| 1383 |
|
\begin{table} |
| 1384 |
|
\caption{EQUATIONS OF MOTION DUE TO POTENTIAL AND KINETIC ENERGIES} |
| 1385 |
+ |
\label{introTable:rbEquations} |
| 1386 |
|
\begin{center} |
| 1387 |
|
\begin{tabular}{|l|l|} |
| 1388 |
|
\hline |
