ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/tengDissertation/Introduction.tex
(Generate patch)

Comparing trunk/tengDissertation/Introduction.tex (file contents):
Revision 2720 by tim, Wed Apr 19 19:46:53 2006 UTC vs.
Revision 2725 by tim, Fri Apr 21 05:45:14 2006 UTC

# Line 892 | Line 892 | T(t) = \sum\limits_{i = 1}^N {\frac{{m_i v_i^2 }}{{fk_
892   simulations. For instance, instantaneous temperature of an
893   Hamiltonian system of $N$ particle can be measured by
894   \[
895 < T(t) = \sum\limits_{i = 1}^N {\frac{{m_i v_i^2 }}{{fk_B }}}
895 > T = \sum\limits_{i = 1}^N {\frac{{m_i v_i^2 }}{{fk_B }}}
896   \]
897   where $m_i$ and $v_i$ are the mass and velocity of $i$th particle
898   respectively, $f$ is the number of degrees of freedom, and $k_B$ is
# Line 913 | Line 913 | discusses issues in production run, including the forc
913   These three individual steps will be covered in the following
914   sections. Sec.~\ref{introSec:initialSystemSettings} deals with the
915   initialization of a simulation. Sec.~\ref{introSec:production} will
916 < discusses issues in production run, including the force evaluation
917 < and the numerical integration schemes of the equations of motion .
918 < Sec.~\ref{introSection:Analysis} provides the theoretical tools for
919 < trajectory analysis.
916 > discusses issues in production run. Sec.~\ref{introSection:Analysis}
917 > provides the theoretical tools for trajectory analysis.
918  
919   \subsection{\label{introSec:initialSystemSettings}Initialization}
920  
# Line 986 | Line 984 | way.
984  
985   \subsection{\label{introSection:production}Production}
986  
987 < \subsubsection{\label{introSec:forceCalculation}The Force Calculation}
987 > Production run is the most important steps of the simulation, in
988 > which the equilibrated structure is used as a starting point and the
989 > motions of the molecules are collected for later analysis. In order
990 > to capture the macroscopic properties of the system, the molecular
991 > dynamics simulation must be performed in correct and efficient way.
992  
993 < \subsubsection{\label{introSection:integrationSchemes} Integration
994 < Schemes}
993 > The most expensive part of a molecular dynamics simulation is the
994 > calculation of non-bonded forces, such as van der Waals force and
995 > Coulombic forces \textit{etc}. For a system of $N$ particles, the
996 > complexity of the algorithm for pair-wise interactions is $O(N^2 )$,
997 > which making large simulations prohibitive in the absence of any
998 > computation saving techniques.
999 >
1000 > A natural approach to avoid system size issue is to represent the
1001 > bulk behavior by a finite number of the particles. However, this
1002 > approach will suffer from the surface effect. To offset this,
1003 > \textit{Periodic boundary condition} is developed to simulate bulk
1004 > properties with a relatively small number of particles. In this
1005 > method, the simulation box is replicated throughout space to form an
1006 > infinite lattice. During the simulation, when a particle moves in
1007 > the primary cell, its image in other cells move in exactly the same
1008 > direction with exactly the same orientation. Thus, as a particle
1009 > leaves the primary cell, one of its images will enter through the
1010 > opposite face.
1011 > %\begin{figure}
1012 > %\centering
1013 > %\includegraphics[width=\linewidth]{pbcFig.eps}
1014 > %\caption[An illustration of periodic boundary conditions]{A 2-D
1015 > %illustration of periodic boundary conditions. As one particle leaves
1016 > %the right of the simulation box, an image of it enters the left.}
1017 > %\label{introFig:pbc}
1018 > %\end{figure}
1019  
1020 + %cutoff and minimum image convention
1021 + Another important technique to improve the efficiency of force
1022 + evaluation is to apply cutoff where particles farther than a
1023 + predetermined distance, are not included in the calculation
1024 + \cite{Frenkel1996}. The use of a cutoff radius will cause a
1025 + discontinuity in the potential energy curve
1026 + (Fig.~\ref{introFig:shiftPot}). Fortunately, one can shift the
1027 + potential to ensure the potential curve go smoothly to zero at the
1028 + cutoff radius. Cutoff strategy works pretty well for Lennard-Jones
1029 + interaction because of its short range nature. However, simply
1030 + truncating the electrostatic interaction with the use of cutoff has
1031 + been shown to lead to severe artifacts in simulations. Ewald
1032 + summation, in which the slowly conditionally convergent Coulomb
1033 + potential is transformed into direct and reciprocal sums with rapid
1034 + and absolute convergence, has proved to minimize the periodicity
1035 + artifacts in liquid simulations. Taking the advantages of the fast
1036 + Fourier transform (FFT) for calculating discrete Fourier transforms,
1037 + the particle mesh-based methods are accelerated from $O(N^{3/2})$ to
1038 + $O(N logN)$. An alternative approach is \emph{fast multipole
1039 + method}, which treats Coulombic interaction exactly at short range,
1040 + and approximate the potential at long range through multipolar
1041 + expansion. In spite of their wide acceptances at the molecular
1042 + simulation community, these two methods are hard to be implemented
1043 + correctly and efficiently. Instead, we use a damped and
1044 + charge-neutralized Coulomb potential method developed by Wolf and
1045 + his coworkers. The shifted Coulomb potential for particle $i$ and
1046 + particle $j$ at distance $r_{rj}$ is given by:
1047 + \begin{equation}
1048 + V(r_{ij})= \frac{q_i q_j \textrm{erfc}(\alpha
1049 + r_{ij})}{r_{ij}}-\lim_{r_{ij}\rightarrow
1050 + R_\textrm{c}}\left\{\frac{q_iq_j \textrm{erfc}(\alpha
1051 + r_{ij})}{r_{ij}}\right\}. \label{introEquation:shiftedCoulomb}
1052 + \end{equation}
1053 + where $\alpha$ is the convergence parameter. Due to the lack of
1054 + inherent periodicity and rapid convergence,this method is extremely
1055 + efficient and easy to implement.
1056 + %\begin{figure}
1057 + %\centering
1058 + %\includegraphics[width=\linewidth]{pbcFig.eps}
1059 + %\caption[An illustration of shifted Coulomb potential]{An illustration of shifted Coulomb potential.}
1060 + %\label{introFigure:shiftedCoulomb}
1061 + %\end{figure}
1062 +
1063 + %multiple time step
1064 +
1065   \subsection{\label{introSection:Analysis} Analysis}
1066  
1067 + Recently, advanced visualization technique are widely applied to
1068 + monitor the motions of molecules. Although the dynamics of the
1069 + system can be described qualitatively from animation, quantitative
1070 + trajectory analysis are more appreciable. According to the
1071 + principles of Statistical Mechanics,
1072 + Sec.~\ref{introSection:statisticalMechanics}, one can compute
1073 + thermodynamics properties, analyze fluctuations of structural
1074 + parameters, and investigate time-dependent processes of the molecule
1075 + from the trajectories.
1076 +
1077 + \subsubsection{\label{introSection:thermodynamicsProperties}Thermodynamics Properties}
1078 +
1079 + Thermodynamics properties, which can be expressed in terms of some
1080 + function of the coordinates and momenta of all particles in the
1081 + system, can be directly computed from molecular dynamics. The usual
1082 + way to measure the pressure is based on virial theorem of Clausius
1083 + which states that the virial is equal to $-3Nk_BT$. For a system
1084 + with forces between particles, the total virial, $W$, contains the
1085 + contribution from external pressure and interaction between the
1086 + particles:
1087 + \[
1088 + W =  - 3PV + \left\langle {\sum\limits_{i < j} {r{}_{ij} \cdot
1089 + f_{ij} } } \right\rangle
1090 + \]
1091 + where $f_{ij}$ is the force between particle $i$ and $j$ at a
1092 + distance $r_{ij}$. Thus, the expression for the pressure is given
1093 + by:
1094 + \begin{equation}
1095 + P = \frac{{Nk_B T}}{V} - \frac{1}{{3V}}\left\langle {\sum\limits_{i
1096 + < j} {r{}_{ij} \cdot f_{ij} } } \right\rangle
1097 + \end{equation}
1098 +
1099 + \subsubsection{\label{introSection:structuralProperties}Structural Properties}
1100 +
1101 + Structural Properties of a simple fluid can be described by a set of
1102 + distribution functions. Among these functions,\emph{pair
1103 + distribution function}, also known as \emph{radial distribution
1104 + function}, is of most fundamental importance to liquid-state theory.
1105 + Pair distribution function can be gathered by Fourier transforming
1106 + raw data from a series of neutron diffraction experiments and
1107 + integrating over the surface factor \cite{Powles73}. The experiment
1108 + result can serve as a criterion to justify the correctness of the
1109 + theory. Moreover, various equilibrium thermodynamic and structural
1110 + properties can also be expressed in terms of radial distribution
1111 + function \cite{allen87:csl}.
1112 +
1113 + A pair distribution functions $g(r)$ gives the probability that a
1114 + particle $i$ will be located at a distance $r$ from a another
1115 + particle $j$ in the system
1116 + \[
1117 + g(r) = \frac{V}{{N^2 }}\left\langle {\sum\limits_i {\sum\limits_{j
1118 + \ne i} {\delta (r - r_{ij} )} } } \right\rangle.
1119 + \]
1120 + Note that the delta function can be replaced by a histogram in
1121 + computer simulation. Figure
1122 + \ref{introFigure:pairDistributionFunction} shows a typical pair
1123 + distribution function for the liquid argon system. The occurrence of
1124 + several peaks in the plot of $g(r)$ suggests that it is more likely
1125 + to find particles at certain radial values than at others. This is a
1126 + result of the attractive interaction at such distances. Because of
1127 + the strong repulsive forces at short distance, the probability of
1128 + locating particles at distances less than about 2.5{\AA} from each
1129 + other is essentially zero.
1130 +
1131 + %\begin{figure}
1132 + %\centering
1133 + %\includegraphics[width=\linewidth]{pdf.eps}
1134 + %\caption[Pair distribution function for the liquid argon
1135 + %]{Pair distribution function for the liquid argon}
1136 + %\label{introFigure:pairDistributionFunction}
1137 + %\end{figure}
1138 +
1139 + \subsubsection{\label{introSection:timeDependentProperties}Time-dependent
1140 + Properties}
1141 +
1142 + Time-dependent properties are usually calculated using \emph{time
1143 + correlation function}, which correlates random variables $A$ and $B$
1144 + at two different time
1145 + \begin{equation}
1146 + C_{AB} (t) = \left\langle {A(t)B(0)} \right\rangle.
1147 + \label{introEquation:timeCorrelationFunction}
1148 + \end{equation}
1149 + If $A$ and $B$ refer to same variable, this kind of correlation
1150 + function is called \emph{auto correlation function}. One example of
1151 + auto correlation function is velocity auto-correlation function
1152 + which is directly related to transport properties of molecular
1153 + liquids:
1154 + \[
1155 + D = \frac{1}{3}\int\limits_0^\infty  {\left\langle {v(t) \cdot v(0)}
1156 + \right\rangle } dt
1157 + \]
1158 + where $D$ is diffusion constant. Unlike velocity autocorrelation
1159 + function which is averaging over time origins and over all the
1160 + atoms, dipole autocorrelation are calculated for the entire system.
1161 + The dipole autocorrelation function is given by:
1162 + \[
1163 + c_{dipole}  = \left\langle {u_{tot} (t) \cdot u_{tot} (t)}
1164 + \right\rangle
1165 + \]
1166 + Here $u_{tot}$ is the net dipole of the entire system and is given
1167 + by
1168 + \[
1169 + u_{tot} (t) = \sum\limits_i {u_i (t)}
1170 + \]
1171 + In principle, many time correlation functions can be related with
1172 + Fourier transforms of the infrared, Raman, and inelastic neutron
1173 + scattering spectra of molecular liquids. In practice, one can
1174 + extract the IR spectrum from the intensity of dipole fluctuation at
1175 + each frequency using the following relationship:
1176 + \[
1177 + \hat c_{dipole} (v) = \int_{ - \infty }^\infty  {c_{dipole} (t)e^{ -
1178 + i2\pi vt} dt}
1179 + \]
1180 +
1181   \section{\label{introSection:rigidBody}Dynamics of Rigid Bodies}
1182  
1183   Rigid bodies are frequently involved in the modeling of different
# Line 1854 | Line 2039 | joining center of resistance $R$ and origin $O$.
2039   \]
2040   where $x_OR$, $y_OR$, $z_OR$ are the components of the vector
2041   joining center of resistance $R$ and origin $O$.
1857
1858 %\section{\label{introSection:correlationFunctions}Correlation Functions}

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines