[ViewVC] Diff of: group/trunk/fennellDissertation/dissertation.tex

Comparing trunk/fennellDissertation/dissertation.tex (file contents):
Revision 2920 by chrisfen, Mon Jul 3 14:51:28 2006 UTC vs.
Revision 2971 by chrisfen, Sun Aug 6 22:29:27 2006 UTC

-<
+\documentclass[12pt]{ndthesis}
->
+\documentclass[11pt]{ndthesis}
+% some packages for things like equations and graphics
-+
+\usepackage[tbtags]{amsmath}
+\usepackage{amsmath,bm}
+\usepackage{amssymb}
+\usepackage{mathrsfs}
+\usepackage{tabularx}
+\usepackage{graphicx}
+\usepackage{booktabs}
-+
+\usepackage{cite}
-+
+\usepackage{enumitem}
+\begin{document}
+structural and dynamic properties of water compared the same
+properties obtained using the Ewald sum.
-<
+\section{Simple Forms for Pairwise Electrostatics}
->
+\section{Simple Forms for Pairwise Electrostatics}\label{sec:PairwiseDerivation}
+The potentials proposed by Wolf \textit{et al.} and Zahn \textit{et
+al.} are constructed using two different (and separable) computational
+tricks:
-<
+\begin{enumerate}
->
+\begin{enumerate}[itemsep=0pt]
+\item shifting through the use of image charges, and
+\item damping the electrostatic interaction.
+\end{enumerate}
-<
+Wolf \textit{et al.} treated the
-<
+development of their summation method as a progressive application of
-<
+these techniques,\cite{Wolf99} while Zahn \textit{et al.} founded
-<
+their damped Coulomb modification (Eq. (\ref{eq:ZahnPot})) on the
-<
+post-limit forces (Eq. (\ref{eq:WolfForces})) which were derived using
-<
+both techniques.  It is possible, however, to separate these
-<
+tricks and study their effects independently.
->
+Wolf \textit{et al.} treated the development of their summation method
->
+as a progressive application of these techniques,\cite{Wolf99} while
->
+Zahn \textit{et al.} founded their damped Coulomb modification
->
+(Eq. (\ref{eq:ZahnPot})) on the post-limit forces
->
+(Eq. (\ref{eq:WolfForces})) which were derived using both techniques.
->
+It is possible, however, to separate these tricks and study their
->
+effects independently.
+Starting with the original observation that the effective range of the
+electrostatic interaction in condensed phases is considerably less
+differences.
+Results and discussion for the individual analysis of each of the
-<
+system types appear in sections \ref{sec:SystemResults}, while the
->
+system types appear in sections \ref{sec:IndividualResults}, while the
+cumulative results over all the investigated systems appear below in
+sections \ref{sec:EnergyResults}.
+simulations (i.e. from liquids of neutral molecules to ionic
+crystals), so the systems studied were:
-<
+\begin{enumerate}
->
+\begin{enumerate}[itemsep=0pt]
+\item liquid water (SPC/E),\cite{Berendsen87}
+\item crystalline water (Ice I$_\textrm{c}$ crystals of SPC/E),
+\item NaCl crystals,
+We compared the following alternative summation methods with results
+from the reference method ({\sc spme}):
-<
+\begin{enumerate}
->
+\begin{enumerate}[itemsep=0pt]
+\item {\sc sp} with damping parameters ($\alpha$) of 0.0, 0.1, 0.2,
+and 0.3\AA$^{-1}$,
+\item {\sc sf} with damping parameters ($\alpha$) of 0.0, 0.1, 0.2,
+.3119, and 0.2476\AA$^{-1}$ for cutoff radii of 9, 12, and 15\AA\
+respectively.
-<
+\section{Combined Configuration Energy Difference Results}\label{sec:EnergyResults}
->
+\section{Configuration Energy Difference Results}\label{sec:EnergyResults}
+In order to evaluate the performance of the pairwise electrostatic
+summation methods for Monte Carlo (MC) simulations, the energy
+differences between configurations were compared to the values
+function.\cite{Adams79,Steinbach94,Leach01} However, we do not see
+significant improvement using the group-switched cutoff because the
+salt and salt solution systems contain non-neutral groups.  Section
-<
+\ref{sec:SystemResults} includes results for systems comprised entirely
->
+\ref{sec:IndividualResults} includes results for systems comprised entirely
+of neutral groups.
+For the {\sc sp} method, inclusion of electrostatic damping improves
+the complementary error function is required).
+The reaction field results illustrates some of that method's
-<
+limitations, primarily that it was developed for use in homogenous
->
+limitations, primarily that it was developed for use in homogeneous
+systems; although it does provide results that are an improvement over
+those from an unmodified cutoff.
-<
+\section{Magnitude of the Force and Torque Vector Results}
->
+\section{Magnitude of the Force and Torque Vector Results}\label{sec:FTMagResults}
+Evaluation of pairwise methods for use in Molecular Dynamics
+simulations requires consideration of effects on the forces and
+molecular bodies. Therefore it is not surprising that reaction field
+performs best of all of the methods on molecular torques.
-<
+\section{Directionality of the Force and Torque Vector Results}
->
+\section{Directionality of the Force and Torque Vector Results}\label{sec:FTDirResults}
+It is clearly important that a new electrostatic method can reproduce
+the magnitudes of the force and torque vectors obtained via the Ewald
+all do equivalently well at capturing the direction of both the force
+and torque vectors.  Using the electrostatic damping improves the
+angular behavior significantly for the {\sc sp} and moderately for the
-<
+{\sc sf} methods.  Overdamping is detrimental to both methods.  Again
->
+{\sc sf} methods.  Over-damping is detrimental to both methods.  Again
+it is important to recognize that the force vectors cover all
+particles in all seven systems, while torque vectors are only
+available for neutral molecular groups.  Damping is more beneficial to
+charged bodies, and this observation is investigated further in
-<
+section \ref{SystemResults}.
->
+section \ref{sec:IndividualResults}.
+Although not discussed previously, group based cutoffs can be applied
+to both the {\sc sp} and {\sc sf} methods.  The group-based cutoffs
+\footnotesize
+\begin{center}
-<
+\begin{tabular}{@{} ccrrrrrrrr @{}} \\
->
+\begin{tabular}{@{} ccrrrrrrrr @{}}
+\toprule
+\toprule
-–
+& &\multicolumn{4}{c}{Shifted Potential} & \multicolumn{4}{c}{Shifted
+Force} \\
+\cmidrule(lr){3-6} \cmidrule(l){7-10} $R_\textrm{c}$ & Groups &
+increases, something that is more obvious with group-based cutoffs.
+The complimentary error function inserted into the potential weakens
+the electrostatic interaction as the value of $\alpha$ is increased.
-<
+However, at larger values of $\alpha$, it is possible to overdamp the
->
+However, at larger values of $\alpha$, it is possible to over-damp the
+electrostatic interaction and to remove it completely.  Kast
+\textit{et al.}  developed a method for choosing appropriate $\alpha$
+values for these types of electrostatic summation methods by fitting
+observations, empirical damping up to 0.2\AA$^{-1}$ is beneficial,
+but damping may be unnecessary when using the {\sc sf} method.
-<
+\section{Individual System Analysis Results}
->
+\section{Individual System Analysis Results}\label{sec:IndividualResults}
+The combined results of the previous sections show how the pairwise
+methods compare to the Ewald summation in the general sense over all
+vector, and torque vector analyses are presented on an individual
+system basis.
-<
+\subsection{SPC/E Water Results}
->
+\subsection{SPC/E Water Results}\label{sec:WaterResults}
-<
+\subsection{SPC/E Ice I$_\textrm{c}$ Results}
->
+The first system considered was liquid water at 300K using the SPC/E
->
+model of water.\cite{Berendsen87} The results for the energy gap
->
+comparisons and the force and torque vector magnitude comparisons are
->
+shown in table \ref{tab:spce}.  The force and torque vector
->
+directionality results are displayed separately in table
->
+\ref{tab:spceAng}, where the effect of group-based cutoffs and
->
+switching functions on the {\sc sp} and {\sc sf} potentials are also
->
+investigated.  In all of the individual results table, the method
->
+abbreviations are as follows:
-<
+\subsection{NaCl Melt Results}
->
+\begin{itemize}[itemsep=0pt]
->
+\item PC = Pure Cutoff,
->
+\item SP = Shifted Potential,
->
+\item SF = Shifted Force,
->
+\item GSC = Group Switched Cutoff,
->
+\item RF = Reaction Field (where $\varepsilon \approx\infty$),
->
+\item GSSP = Group Switched Shifted Potential, and
->
+\item GSSF = Group Switched Shifted Force.
->
+\end{itemize}
-<
+\subsection{NaCl Crystal Results}
->
+\begin{table}[htbp]
->
+\centering
->
+\caption{REGRESSION RESULTS OF THE LIQUID WATER SYSTEM FOR THE
->
+$\Delta E$ VALUES ({\it upper}), FORCE VECTOR MAGNITUDES ({\it middle})
->
+AND TORQUE VECTOR MAGNITUDES ({\it lower})}
-<
+\subsection{0.1M NaCl Solution Results}
->
+\footnotesize
->
+\begin{tabular}{@{} ccrrrrrr @{}}
->
+\toprule
->
+\toprule
->
+& & \multicolumn{2}{c}{9\AA} & \multicolumn{2}{c}{12\AA} & \multicolumn{2}{c}{15\AA}\\
->
+\cmidrule(lr){3-4}
->
+\cmidrule(lr){5-6}
->
+\cmidrule(l){7-8}
->
+Method & $\alpha$ & slope & $R^2$ & slope & $R^2$ & slope & $R^2$ \\
->
+\midrule
->
+PC  &     & 3.046 & 0.002 & -3.018 & 0.002 & 4.719 & 0.005 \\
->
+SP  & 0.0 & 1.035 & 0.218 & 0.908 & 0.313 & 1.037 & 0.470 \\
->
+    & 0.1 & 1.021 & 0.387 & 0.965 & 0.752 & 1.006 & 0.947 \\
->
+    & 0.2 & 0.997 & 0.962 & 1.001 & 0.994 & 0.994 & 0.996 \\
->
+    & 0.3 & 0.984 & 0.980 & 0.997 & 0.985 & 0.982 & 0.987 \\
->
+SF  & 0.0 & 0.977 & 0.974 & 0.996 & 0.992 & 0.991 & 0.997 \\
->
+    & 0.1 & 0.983 & 0.974 & 1.001 & 0.994 & 0.996 & 0.998 \\
->
+    & 0.2 & 0.992 & 0.989 & 1.001 & 0.995 & 0.994 & 0.996 \\
->
+    & 0.3 & 0.984 & 0.980 & 0.996 & 0.985 & 0.982 & 0.987 \\
->
+GSC &     & 0.918 & 0.862 & 0.852 & 0.756 & 0.801 & 0.700 \\
->
+RF  &     & 0.971 & 0.958 & 0.975 & 0.987 & 0.959 & 0.983 \\
->
+\midrule
->
+PC  &     & -1.647 & 0.000 & -0.127 & 0.000 & -0.979 & 0.000 \\
->
+SP  & 0.0 & 0.735 & 0.368 & 0.813 & 0.537 & 0.865 & 0.659 \\
->
+    & 0.1 & 0.850 & 0.612 & 0.956 & 0.887 & 0.992 & 0.979 \\
->
+    & 0.2 & 0.996 & 0.989 & 1.000 & 1.000 & 1.000 & 1.000 \\
->
+    & 0.3 & 0.996 & 0.998 & 0.997 & 0.998 & 0.996 & 0.998 \\
->
+SF  & 0.0 & 0.998 & 0.995 & 1.000 & 0.999 & 1.000 & 0.999 \\
->
+    & 0.1 & 0.998 & 0.995 & 1.000 & 0.999 & 1.000 & 1.000 \\
->
+    & 0.2 & 0.999 & 0.998 & 1.000 & 1.000 & 1.000 & 1.000 \\
->
+    & 0.3 & 0.996 & 0.998 & 0.997 & 0.998 & 0.996 & 0.998 \\
->
+GSC &     & 0.998 & 0.995 & 1.000 & 0.999 & 1.000 & 1.000 \\
->
+RF  &     & 0.999 & 0.995 & 1.000 & 0.999 & 1.000 & 1.000 \\
->
+\midrule
->
+PC  &     & 2.387 & 0.000 & 0.183 & 0.000 & 1.282 & 0.000 \\
->
+SP  & 0.0 & 0.847 & 0.543 & 0.904 & 0.694 & 0.935 & 0.786 \\
->
+    & 0.1 & 0.922 & 0.749 & 0.980 & 0.934 & 0.996 & 0.988 \\
->
+    & 0.2 & 0.987 & 0.985 & 0.989 & 0.992 & 0.990 & 0.993 \\
->
+    & 0.3 & 0.965 & 0.973 & 0.967 & 0.975 & 0.967 & 0.976 \\
->
+SF  & 0.0 & 0.978 & 0.990 & 0.988 & 0.997 & 0.993 & 0.999 \\
->
+    & 0.1 & 0.983 & 0.991 & 0.993 & 0.997 & 0.997 & 0.999 \\
->
+    & 0.2 & 0.986 & 0.989 & 0.989 & 0.992 & 0.990 & 0.993 \\
->
+    & 0.3 & 0.965 & 0.973 & 0.967 & 0.975 & 0.967 & 0.976 \\
->
+GSC &     & 0.995 & 0.981 & 0.999 & 0.991 & 1.001 & 0.994 \\
->
+RF  &     & 0.993 & 0.989 & 0.998 & 0.996 & 1.000 & 0.999 \\
->
+\bottomrule
->
+\end{tabular}
->
+\label{tab:spce}
->
+\end{table}
-<
+\subsection{1M NaCl Solution Results}
->
+\begin{table}[htbp]
->
+\centering
->
+\caption{VARIANCE RESULTS FROM GAUSSIAN FITS TO ANGULAR
->
+DISTRIBUTIONS OF THE FORCE AND TORQUE VECTORS IN THE LIQUID WATER
->
+SYSTEM}
->
->
+\footnotesize
->
+\begin{tabular}{@{} ccrrrrrr @{}}
->
+\toprule
->
+\toprule
->
+& & \multicolumn{3}{c}{Force $\sigma^2$} & \multicolumn{3}{c}{Torque $\sigma^2$} \\
->
+\cmidrule(lr){3-5}
->
+\cmidrule(l){6-8}
->
+Method & $\alpha$ & 9\AA & 12\AA & 15\AA & 9\AA & 12\AA & 15\AA \\
->
+\midrule
->
+PC  &     & 783.759 & 481.353 & 332.677 & 248.674 & 144.382 & 98.535 \\
->
+SP  & 0.0 & 659.440 & 380.699 & 250.002 & 235.151 & 134.661 & 88.135 \\
->
+    & 0.1 & 293.849 & 67.772 & 11.609 & 105.090 & 23.813 & 4.369 \\
->
+    & 0.2 & 5.975 & 0.136 & 0.094 & 5.553 & 1.784 & 1.536 \\
->
+    & 0.3 & 0.725 & 0.707 & 0.693 & 7.293 & 6.933 & 6.748 \\
->
+SF  & 0.0 & 2.238 & 0.713 & 0.292 & 3.290 & 1.090 & 0.416 \\
->
+    & 0.1 & 2.238 & 0.524 & 0.115 & 3.184 & 0.945 & 0.326 \\
->
+    & 0.2 & 0.374 & 0.102 & 0.094 & 2.598 & 1.755 & 1.537 \\
->
+    & 0.3 & 0.721 & 0.707 & 0.693 & 7.322 & 6.933 & 6.748 \\
->
+GSC &     & 2.431 & 0.614 & 0.274 & 5.135 & 2.133 & 1.339 \\
->
+RF  &     & 2.091 & 0.403 & 0.113 & 3.583 & 1.071 & 0.399 \\
->
+\midrule
->
+GSSP  & 0.0 & 2.431 & 0.614 & 0.274 & 5.135 & 2.133 & 1.339 \\
->
+      & 0.1 & 1.879 & 0.291 & 0.057 & 3.983 & 1.117 & 0.370 \\
->
+      & 0.2 & 0.443 & 0.103 & 0.093 & 2.821 & 1.794 & 1.532 \\
->
+      & 0.3 & 0.728 & 0.694 & 0.692 & 7.387 & 6.942 & 6.748 \\
->
+GSSF  & 0.0 & 1.298 & 0.270 & 0.083 & 3.098 & 0.992 & 0.375 \\
->
+      & 0.1 & 1.296 & 0.210 & 0.044 & 3.055 & 0.922 & 0.330 \\
->
+      & 0.2 & 0.433 & 0.104 & 0.093 & 2.895 & 1.797 & 1.532 \\
->
+      & 0.3 & 0.728 & 0.694 & 0.692 & 7.410 & 6.942 & 6.748 \\
->
+\bottomrule
->
+\end{tabular}
->
+\label{tab:spceAng}
->
+\end{table}
->
->
+The water results parallel the combined results seen in sections
->
+\ref{sec:EnergyResults} through \ref{sec:FTDirResults}.  There is good
->
+agreement with {\sc spme} in both energetic and dynamic behavior when
->
+using the {\sc sf} method with and without damping. The {\sc sp}
->
+method does well with an $\alpha$ around 0.2\AA$^{-1}$, particularly
->
+with cutoff radii greater than 12\AA. Over-damping the electrostatics
->
+reduces the agreement between both these methods and {\sc spme}.
->
->
+The pure cutoff ({\sc pc}) method performs poorly, again mirroring the
->
+observations from the combined results.  In contrast to these results, however, the use of a switching function and group
->
+based cutoffs greatly improves the results for these neutral water
->
+molecules.  The group switched cutoff ({\sc gsc}) does not mimic the
->
+energetics of {\sc spme} as well as the {\sc sp} (with moderate
->
+damping) and {\sc sf} methods, but the dynamics are quite good.  The
->
+switching functions correct discontinuities in the potential and
->
+forces, leading to these improved results.  Such improvements with the
->
+use of a switching function have been recognized in previous
->
+studies,\cite{Andrea83,Steinbach94} and this proves to be a useful
->
+tactic for stably incorporating local area electrostatic effects.
-+
+The reaction field ({\sc rf}) method simply extends upon the results
-+
+observed in the {\sc gsc} case.  Both methods are similar in form
-+
+(i.e. neutral groups, switching function), but {\sc rf} incorporates
-+
+an added effect from the external dielectric. This similarity
-+
+translates into the same good dynamic results and improved energetic
-+
+agreement with {\sc spme}.  Though this agreement is not to the level
-+
+of the moderately damped {\sc sp} and {\sc sf} methods, these results
-+
+show how incorporating some implicit properties of the surroundings
-+
+(i.e. $\epsilon_\textrm{S}$) can improve the solvent depiction.
-+
-+
+As a final note for the liquid water system, use of group cutoffs and a
-+
+switching function leads to noticeable improvements in the {\sc sp}
-+
+and {\sc sf} methods, primarily in directionality of the force and
-+
+torque vectors (table \ref{tab:spceAng}). The {\sc sp} method shows
-+
+significant narrowing of the angle distribution when using little to
-+
+no damping and only modest improvement for the recommended conditions
-+
+($\alpha = 0.2$\AA$^{-1}$ and $R_\textrm{c}~\geqslant~12$\AA).  The
-+
+{\sc sf} method shows modest narrowing across all damping and cutoff
-+
+ranges of interest.  When over-damping these methods, group cutoffs and
-+
+the switching function do not improve the force and torque
-+
+directionalities.
-+
-+
+\subsection{SPC/E Ice I$_\textrm{c}$ Results}\label{sec:IceResults}
-+
-+
+In addition to the disordered molecular system above, the ordered
-+
+molecular system of ice I$_\textrm{c}$ was also considered.  Ice
-+
+polymorph could have been used to fit this role; however, ice
-+
+I$_\textrm{c}$ was chosen because it can form an ideal periodic
-+
+lattice with the same number of water molecules used in the disordered
-+
+liquid state case.  The results for the energy gap comparisons and the
-+
+force and torque vector magnitude comparisons are shown in table
-+
+\ref{tab:ice}.  The force and torque vector directionality results are
-+
+displayed separately in table \ref{tab:iceAng}, where the effect of
-+
+group-based cutoffs and switching functions on the {\sc sp} and {\sc
-+
+sf} potentials are also displayed.
-+
-+
+\begin{table}[htbp]
-+
+\centering
-+
+\caption{REGRESSION RESULTS OF THE ICE I$_\textrm{c}$ SYSTEM FOR
-+
+$\Delta E$ VALUES ({\it upper}), FORCE VECTOR MAGNITUDES ({\it
-+
+middle}) AND TORQUE VECTOR MAGNITUDES ({\it lower})}
-+
-+
+\footnotesize
-+
+\begin{tabular}{@{} ccrrrrrr @{}}
-+
+\toprule
-+
+\toprule
-+
+& & \multicolumn{2}{c}{9\AA} & \multicolumn{2}{c}{12\AA} & \multicolumn{2}{c}{15\AA}\\
-+
+\cmidrule(lr){3-4}
-+
+\cmidrule(lr){5-6}
-+
+\cmidrule(l){7-8}
-+
+Method & $\alpha$ & slope & $R^2$ & slope & $R^2$ & slope & $R^2$ \\
-+
+\midrule
-+
+PC  &     & 19.897 & 0.047 & -29.214 & 0.048 & -3.771 & 0.001 \\
-+
+SP  & 0.0 & -0.014 & 0.000 & 2.135 & 0.347 & 0.457 & 0.045 \\
-+
+    & 0.1 & 0.321 & 0.017 & 1.490 & 0.584 & 0.886 & 0.796 \\
-+
+    & 0.2 & 0.896 & 0.872 & 1.011 & 0.998 & 0.997 & 0.999 \\
-+
+    & 0.3 & 0.983 & 0.997 & 0.992 & 0.997 & 0.991 & 0.997 \\
-+
+SF  & 0.0 & 0.943 & 0.979 & 1.048 & 0.978 & 0.995 & 0.999 \\
-+
+    & 0.1 & 0.948 & 0.979 & 1.044 & 0.983 & 1.000 & 0.999 \\
-+
+    & 0.2 & 0.982 & 0.997 & 0.969 & 0.960 & 0.997 & 0.999 \\
-+
+    & 0.3 & 0.985 & 0.997 & 0.961 & 0.961 & 0.991 & 0.997 \\
-+
+GSC &     & 0.983 & 0.985 & 0.966 & 0.994 & 1.003 & 0.999 \\
-+
+RF  &     & 0.924 & 0.944 & 0.990 & 0.996 & 0.991 & 0.998 \\
-+
+\midrule
-+
+PC  &     & -4.375 & 0.000 & 6.781 & 0.000 & -3.369 & 0.000 \\
-+
+SP  & 0.0 & 0.515 & 0.164 & 0.856 & 0.426 & 0.743 & 0.478 \\
-+
+    & 0.1 & 0.696 & 0.405 & 0.977 & 0.817 & 0.974 & 0.964 \\
-+
+    & 0.2 & 0.981 & 0.980 & 1.001 & 1.000 & 1.000 & 1.000 \\
-+
+    & 0.3 & 0.996 & 0.998 & 0.997 & 0.999 & 0.997 & 0.999 \\
-+
+SF  & 0.0 & 0.991 & 0.995 & 1.003 & 0.998 & 0.999 & 1.000 \\
-+
+    & 0.1 & 0.992 & 0.995 & 1.003 & 0.998 & 1.000 & 1.000 \\
-+
+    & 0.2 & 0.998 & 0.998 & 0.981 & 0.962 & 1.000 & 1.000 \\
-+
+    & 0.3 & 0.996 & 0.998 & 0.976 & 0.957 & 0.997 & 0.999 \\
-+
+GSC &     & 0.997 & 0.996 & 0.998 & 0.999 & 1.000 & 1.000 \\
-+
+RF  &     & 0.988 & 0.989 & 1.000 & 0.999 & 1.000 & 1.000 \\
-+
+\midrule
-+
+PC  &     & -6.367 & 0.000 & -3.552 & 0.000 & -3.447 & 0.000 \\
-+
+SP  & 0.0 & 0.643 & 0.409 & 0.833 & 0.607 & 0.961 & 0.805 \\
-+
+    & 0.1 & 0.791 & 0.683 & 0.957 & 0.914 & 1.000 & 0.989 \\
-+
+    & 0.2 & 0.974 & 0.991 & 0.993 & 0.998 & 0.993 & 0.998 \\
-+
+    & 0.3 & 0.976 & 0.992 & 0.977 & 0.992 & 0.977 & 0.992 \\
-+
+SF  & 0.0 & 0.979 & 0.997 & 0.992 & 0.999 & 0.994 & 1.000 \\
-+
+    & 0.1 & 0.984 & 0.997 & 0.996 & 0.999 & 0.998 & 1.000 \\
-+
+    & 0.2 & 0.991 & 0.997 & 0.974 & 0.958 & 0.993 & 0.998 \\
-+
+    & 0.3 & 0.977 & 0.992 & 0.956 & 0.948 & 0.977 & 0.992 \\
-+
+GSC &     & 0.999 & 0.997 & 0.996 & 0.999 & 1.002 & 1.000 \\
-+
+RF  &     & 0.994 & 0.997 & 0.997 & 0.999 & 1.000 & 1.000 \\
-+
+\bottomrule
-+
+\end{tabular}
-+
+\label{tab:ice}
-+
+\end{table}
-+
-+
+\begin{table}[htbp]
-+
+\centering
-+
+\caption{VARIANCE RESULTS FROM GAUSSIAN FITS TO ANGULAR DISTRIBUTIONS
-+
+OF THE FORCE AND TORQUE VECTORS IN THE ICE I$_\textrm{c}$ SYSTEM}
-+
-+
+\footnotesize
-+
+\begin{tabular}{@{} ccrrrrrr @{}}
-+
+\toprule
-+
+\toprule
-+
+& & \multicolumn{3}{c}{Force $\sigma^2$} & \multicolumn{3}{c}{Torque
-+
+$\sigma^2$} \\
-+
+\cmidrule(lr){3-5}
-+
+\cmidrule(l){6-8}
-+
+Method & $\alpha$ & 9\AA & 12\AA & 15\AA & 9\AA & 12\AA & 15\AA \\
-+
+\midrule
-+
+PC  &     & 2128.921 & 603.197 & 715.579 & 329.056 & 221.397 & 81.042 \\
-+
+SP  & 0.0 & 1429.341 & 470.320 & 447.557 & 301.678 & 197.437 & 73.840 \\
-+
+    & 0.1 & 590.008 & 107.510 & 18.883 & 118.201 & 32.472 & 3.599 \\
-+
+    & 0.2 & 10.057 & 0.105 & 0.038 & 2.875 & 0.572 & 0.518 \\
-+
+    & 0.3 & 0.245 & 0.260 & 0.262 & 2.365 & 2.396 & 2.327 \\
-+
+SF  & 0.0 & 1.745 & 1.161 & 0.212 & 1.135 & 0.426 & 0.155 \\
-+
+    & 0.1 & 1.721 & 0.868 & 0.082 & 1.118 & 0.358 & 0.118 \\
-+
+    & 0.2 & 0.201 & 0.040 & 0.038 & 0.786 & 0.555 & 0.518 \\
-+
+    & 0.3 & 0.241 & 0.260 & 0.262 & 2.368 & 2.400 & 2.327 \\
-+
+GSC &     & 1.483 & 0.261 & 0.099 & 0.926 & 0.295 & 0.095 \\
-+
+RF  &     & 2.887 & 0.217 & 0.107 & 1.006 & 0.281 & 0.085 \\
-+
+\midrule
-+
+GSSP  & 0.0 & 1.483 & 0.261 & 0.099 & 0.926 & 0.295 & 0.095 \\
-+
+      & 0.1 & 1.341 & 0.123 & 0.037 & 0.835 & 0.234 & 0.085 \\
-+
+      & 0.2 & 0.558 & 0.040 & 0.037 & 0.823 & 0.557 & 0.519 \\
-+
+      & 0.3 & 0.250 & 0.251 & 0.259 & 2.387 & 2.395 & 2.328 \\
-+
+GSSF  & 0.0 & 2.124 & 0.132 & 0.069 & 0.919 & 0.263 & 0.099 \\
-+
+      & 0.1 & 2.165 & 0.101 & 0.035 & 0.895 & 0.244 & 0.096 \\
-+
+      & 0.2 & 0.706 & 0.040 & 0.037 & 0.870 & 0.559 & 0.519 \\
-+
+      & 0.3 & 0.251 & 0.251 & 0.259 & 2.387 & 2.395 & 2.328 \\
-+
+\bottomrule
-+
+\end{tabular}
-+
+\label{tab:iceAng}
-+
+\end{table}
-+
-+
+Highly ordered systems are a difficult test for the pairwise methods
-+
+in that they lack the implicit periodicity of the Ewald summation.  As
-+
+expected, the energy gap agreement with {\sc spme} is reduced for the
-+
+{\sc sp} and {\sc sf} methods with parameters that were ideal for the
-+
+disordered liquid system.  Moving to higher $R_\textrm{c}$ helps
-+
+improve the agreement, though at an increase in computational cost.
-+
+The dynamics of this crystalline system (both in magnitude and
-+
+direction) are little affected. Both methods still reproduce the Ewald
-+
+behavior with the same parameter recommendations from the previous
-+
+section.
-+
-+
+It is also worth noting that {\sc rf} exhibits improved energy gap
-+
+results over the liquid water system.  One possible explanation is
-+
+that the ice I$_\textrm{c}$ crystal is ordered such that the net
-+
+dipole moment of the crystal is zero.  With $\epsilon_\textrm{S} =
-+
+\infty$, the reaction field incorporates this structural organization
-+
+by actively enforcing a zeroed dipole moment within each cutoff
-+
+sphere.
-+
-+
+\subsection{NaCl Melt Results}\label{sec:SaltMeltResults}
-+
-+
+A high temperature NaCl melt was tested to gauge the accuracy of the
-+
+pairwise summation methods in a disordered system of charges. The
-+
+results for the energy gap comparisons and the force vector magnitude
-+
+comparisons are shown in table \ref{tab:melt}.  The force vector
-+
+directionality results are displayed separately in table
-+
+\ref{tab:meltAng}.
-+
-+
+\begin{table}[htbp]
-+
+\centering
-+
+\caption{REGRESSION RESULTS OF THE MOLTEN SODIUM CHLORIDE SYSTEM FOR
-+
+$\Delta E$ VALUES ({\it upper}) AND FORCE VECTOR MAGNITUDES ({\it
-+
+lower})}
-+
-+
+\footnotesize
-+
+\begin{tabular}{@{} ccrrrrrr @{}}
-+
+\toprule
-+
+\toprule
-+
+& & \multicolumn{2}{c}{9\AA} & \multicolumn{2}{c}{12\AA} & \multicolumn{2}{c}{15\AA}\\
-+
+\cmidrule(lr){3-4}
-+
+\cmidrule(lr){5-6}
-+
+\cmidrule(l){7-8}
-+
+Method & $\alpha$ & slope & $R^2$ & slope & $R^2$ & slope & $R^2$ \\
-+
+\midrule
-+
+PC  &     & -0.008 & 0.000 & -0.049 & 0.005 & -0.136 & 0.020 \\
-+
+SP  & 0.0 & 0.928 & 0.996 & 0.931 & 0.998 & 0.950 & 0.999 \\
-+
+    & 0.1 & 0.977 & 0.998 & 0.998 & 1.000 & 0.997 & 1.000 \\
-+
+    & 0.2 & 0.960 & 1.000 & 0.813 & 0.996 & 0.811 & 0.954 \\
-+
+    & 0.3 & 0.671 & 0.994 & 0.439 & 0.929 & 0.535 & 0.831 \\
-+
+SF  & 0.0 & 0.996 & 1.000 & 0.995 & 1.000 & 0.997 & 1.000 \\
-+
+    & 0.1 & 1.021 & 1.000 & 1.024 & 1.000 & 1.007 & 1.000 \\
-+
+    & 0.2 & 0.966 & 1.000 & 0.813 & 0.996 & 0.811 & 0.954 \\
-+
+    & 0.3 & 0.671 & 0.994 & 0.439 & 0.929 & 0.535 & 0.831 \\
-+
+            \midrule
-+
+PC  &     & 1.103 & 0.000 & 0.989 & 0.000 & 0.802 & 0.000 \\
-+
+SP  & 0.0 & 0.973 & 0.981 & 0.975 & 0.988 & 0.979 & 0.992 \\
-+
+    & 0.1 & 0.987 & 0.992 & 0.993 & 0.998 & 0.997 & 0.999 \\
-+
+    & 0.2 & 0.993 & 0.996 & 0.985 & 0.988 & 0.986 & 0.981 \\
-+
+    & 0.3 & 0.956 & 0.956 & 0.940 & 0.912 & 0.948 & 0.929 \\
-+
+SF  & 0.0 & 0.996 & 0.997 & 0.997 & 0.999 & 0.998 & 1.000 \\
-+
+    & 0.1 & 1.000 & 0.997 & 1.001 & 0.999 & 1.000 & 1.000 \\
-+
+    & 0.2 & 0.994 & 0.996 & 0.985 & 0.988 & 0.986 & 0.981 \\
-+
+    & 0.3 & 0.956 & 0.956 & 0.940 & 0.912 & 0.948 & 0.929 \\
-+
+\bottomrule
-+
+\end{tabular}
-+
+\label{tab:melt}
-+
+\end{table}
-+
-+
+\begin{table}[htbp]
-+
+\centering
-+
+\caption{VARIANCE RESULTS FROM GAUSSIAN FITS TO ANGULAR DISTRIBUTIONS
-+
+OF THE FORCE VECTORS IN THE MOLTEN SODIUM CHLORIDE SYSTEM}
-+
-+
+\footnotesize
-+
+\begin{tabular}{@{} ccrrrrrr @{}}
-+
+\toprule
-+
+\toprule
-+
+& & \multicolumn{3}{c}{Force $\sigma^2$} \\
-+
+\cmidrule(lr){3-5}
-+
+\cmidrule(l){6-8}
-+
+Method & $\alpha$ & 9\AA & 12\AA & 15\AA \\
-+
+\midrule
-+
+PC  &     & 13.294 & 8.035 & 5.366 \\
-+
+SP  & 0.0 & 13.316 & 8.037 & 5.385 \\
-+
+    & 0.1 & 5.705 & 1.391 & 0.360 \\
-+
+    & 0.2 & 2.415 & 7.534 & 13.927 \\
-+
+    & 0.3 & 23.769 & 67.306 & 57.252 \\
-+
+SF  & 0.0 & 1.693 & 0.603 & 0.256 \\
-+
+    & 0.1 & 1.687 & 0.653 & 0.272 \\
-+
+    & 0.2 & 2.598 & 7.523 & 13.930 \\
-+
+    & 0.3 & 23.734 & 67.305 & 57.252 \\
-+
+\bottomrule
-+
+\end{tabular}
-+
+\label{tab:meltAng}
-+
+\end{table}
-+
-+
+The molten NaCl system shows more sensitivity to the electrostatic
-+
+damping than the water systems. The most noticeable point is that the
-+
+undamped {\sc sf} method does very well at replicating the {\sc spme}
-+
+configurational energy differences and forces. Light damping appears
-+
+to minimally improve the dynamics, but this comes with a deterioration
-+
+of the energy gap results. In contrast, this light damping improves
-+
+the {\sc sp} energy gaps and forces. Moderate and heavy electrostatic
-+
+damping reduce the agreement with {\sc spme} for both methods. From
-+
+these observations, the undamped {\sc sf} method is the best choice
-+
+for disordered systems of charges.
-+
-+
+\subsection{NaCl Crystal Results}\label{sec:SaltCrystalResults}
-+
-+
+Similar to the use of ice I$_\textrm{c}$ to investigate the role of
-+
+order in molecular systems on the effectiveness of the pairwise
-+
+methods, the 1000K NaCl crystal system was used to investigate the
-+
+accuracy of the pairwise summation methods in an ordered system of
-+
+charged particles. The results for the energy gap comparisons and the
-+
+force vector magnitude comparisons are shown in table \ref{tab:salt}.
-+
+The force vector directionality results are displayed separately in
-+
+table \ref{tab:saltAng}.
-+
-+
+\begin{table}[htbp]
-+
+\centering
-+
+\caption{REGRESSION RESULTS OF THE CRYSTALLINE SODIUM CHLORIDE
-+
+SYSTEM FOR $\Delta E$ VALUES ({\it upper}) AND FORCE VECTOR MAGNITUDES
-+
+({\it lower})}
-+
-+
+\footnotesize
-+
+\begin{tabular}{@{} ccrrrrrr @{}}
-+
+\toprule
-+
+\toprule
-+
+& & \multicolumn{2}{c}{9\AA} & \multicolumn{2}{c}{12\AA} & \multicolumn{2}{c}{15\AA}\\
-+
+\cmidrule(lr){3-4}
-+
+\cmidrule(lr){5-6}
-+
+\cmidrule(l){7-8}
-+
+Method & $\alpha$ & slope & $R^2$ & slope & $R^2$ & slope & $R^2$ \\
-+
+\midrule
-+
+PC  &     & -20.241 & 0.228 & -20.248 & 0.229 & -20.239 & 0.228 \\
-+
+SP  & 0.0 & 1.039 & 0.733 & 2.037 & 0.565 & 1.225 & 0.743 \\
-+
+    & 0.1 & 1.049 & 0.865 & 1.424 & 0.784 & 1.029 & 0.980 \\
-+
+    & 0.2 & 0.982 & 0.976 & 0.969 & 0.980 & 0.960 & 0.980 \\
-+
+    & 0.3 & 0.873 & 0.944 & 0.872 & 0.945 & 0.872 & 0.945 \\
-+
+SF  & 0.0 & 1.041 & 0.967 & 0.994 & 0.989 & 0.957 & 0.993 \\
-+
+    & 0.1 & 1.050 & 0.968 & 0.996 & 0.991 & 0.972 & 0.995 \\
-+
+    & 0.2 & 0.982 & 0.975 & 0.959 & 0.980 & 0.960 & 0.980 \\
-+
+    & 0.3 & 0.873 & 0.944 & 0.872 & 0.945 & 0.872 & 0.944 \\
-+
+\midrule
-+
+PC  &     & 0.795 & 0.000 & 0.792 & 0.000 & 0.793 & 0.000 \\
-+
+SP  & 0.0 & 0.916 & 0.829 & 1.086 & 0.791 & 1.010 & 0.936 \\
-+
+    & 0.1 & 0.958 & 0.917 & 1.049 & 0.943 & 1.001 & 0.995 \\
-+
+    & 0.2 & 0.981 & 0.981 & 0.982 & 0.984 & 0.981 & 0.984 \\
-+
+    & 0.3 & 0.950 & 0.952 & 0.950 & 0.953 & 0.950 & 0.953 \\
-+
+SF  & 0.0 & 1.002 & 0.983 & 0.997 & 0.994 & 0.991 & 0.997 \\
-+
+    & 0.1 & 1.003 & 0.984 & 0.996 & 0.995 & 0.993 & 0.997 \\
-+
+    & 0.2 & 0.983 & 0.980 & 0.981 & 0.984 & 0.981 & 0.984 \\
-+
+    & 0.3 & 0.950 & 0.952 & 0.950 & 0.953 & 0.950 & 0.953 \\
-+
+\bottomrule
-+
+\end{tabular}
-+
+\label{tab:salt}
-+
+\end{table}
-+
-+
+\begin{table}[htbp]
-+
+\centering
-+
+\caption{VARIANCE RESULTS FROM GAUSSIAN FITS TO ANGULAR
-+
+DISTRIBUTIONS OF THE FORCE VECTORS IN THE CRYSTALLINE SODIUM CHLORIDE
-+
+SYSTEM}
-+
-+
+\footnotesize
-+
+\begin{tabular}{@{} ccrrrrrr @{}}
-+
+\toprule
-+
+\toprule
-+
+& & \multicolumn{3}{c}{Force $\sigma^2$} \\
-+
+\cmidrule(lr){3-5}
-+
+\cmidrule(l){6-8}
-+
+Method & $\alpha$ & 9\AA & 12\AA & 15\AA \\
-+
+\midrule
-+
+PC  &     & 111.945 & 111.824 & 111.866 \\
-+
+SP  & 0.0 & 112.414 & 152.215 & 38.087 \\
-+
+    & 0.1 & 52.361 & 42.574 & 2.819 \\
-+
+    & 0.2 & 10.847 & 9.709 & 9.686 \\
-+
+    & 0.3 & 31.128 & 31.104 & 31.029 \\
-+
+SF  & 0.0 & 10.025 & 3.555 & 1.648 \\
-+
+    & 0.1 & 9.462 & 3.303 & 1.721 \\
-+
+    & 0.2 & 11.454 & 9.813 & 9.701 \\
-+
+    & 0.3 & 31.120 & 31.105 & 31.029 \\
-+
+\bottomrule
-+
+\end{tabular}
-+
+\label{tab:saltAng}
-+
+\end{table}
-+
-+
+The crystalline NaCl system is the most challenging test case for the
-+
+pairwise summation methods, as evidenced by the results in tables
-+
+\ref{tab:salt} and \ref{tab:saltAng}. The undamped and weakly damped
-+
+{\sc sf} methods seem to be the best choices. These methods match well
-+
+with {\sc spme} across the energy gap, force magnitude, and force
-+
+directionality tests.  The {\sc sp} method struggles in all cases,
-+
+with the exception of good dynamics reproduction when using weak
-+
+electrostatic damping with a large cutoff radius.
-+
-+
+The moderate electrostatic damping case is not as good as we would
-+
+expect given the long-time dynamics results observed for this system
-+
+(see section \ref{sec:LongTimeDynamics}). Since the data tabulated in
-+
+tables \ref{tab:salt} and \ref{tab:saltAng} are a test of
-+
+instantaneous dynamics, this indicates that good long-time dynamics
-+
+comes in part at the expense of short-time dynamics.
-+
-+
+\subsection{0.11M NaCl Solution Results}
-+
-+
+In an effort to bridge the charged atomic and neutral molecular
-+
+systems, Na$^+$ and Cl$^-$ ion charge defects were incorporated into
-+
+the liquid water system. This low ionic strength system consists of 4
-+
+ions in the 1000 SPC/E water solvent ($\approx$0.11 M). The results
-+
+for the energy gap comparisons and the force and torque vector
-+
+magnitude comparisons are shown in table \ref{tab:solnWeak}.  The
-+
+force and torque vector directionality results are displayed
-+
+separately in table \ref{tab:solnWeakAng}, where the effect of
-+
+group-based cutoffs and switching functions on the {\sc sp} and {\sc
-+
+sf} potentials are investigated.
-+
-+
+\begin{table}[htbp]
-+
+\centering
-+
+\caption{REGRESSION RESULTS OF THE WEAK SODIUM CHLORIDE SOLUTION
-+
+SYSTEM FOR $\Delta E$ VALUES ({\it upper}), FORCE VECTOR MAGNITUDES
-+
+({\it middle}) AND TORQUE VECTOR MAGNITUDES ({\it lower})}
-+
-+
+\footnotesize
-+
+\begin{tabular}{@{} ccrrrrrr @{}}
-+
+\toprule
-+
+\toprule
-+
+& & \multicolumn{2}{c}{9\AA} & \multicolumn{2}{c}{12\AA} & \multicolumn{2}{c}{15\AA}\\
-+
+\cmidrule(lr){3-4}
-+
+\cmidrule(lr){5-6}
-+
+\cmidrule(l){7-8}
-+
+Method & $\alpha$ & slope & $R^2$ & slope & $R^2$ & slope & $R^2$ \\
-+
+\midrule
-+
+PC  &     & 0.247 & 0.000 & -1.103 & 0.001 & 5.480 & 0.015 \\
-+
+SP  & 0.0 & 0.935 & 0.388 & 0.984 & 0.541 & 1.010 & 0.685 \\
-+
+    & 0.1 & 0.951 & 0.603 & 0.993 & 0.875 & 1.001 & 0.979 \\
-+
+    & 0.2 & 0.969 & 0.968 & 0.996 & 0.997 & 0.994 & 0.997 \\
-+
+    & 0.3 & 0.955 & 0.966 & 0.984 & 0.992 & 0.978 & 0.991 \\
-+
+SF  & 0.0 & 0.963 & 0.971 & 0.989 & 0.996 & 0.991 & 0.998 \\
-+
+    & 0.1 & 0.970 & 0.971 & 0.995 & 0.997 & 0.997 & 0.999 \\
-+
+    & 0.2 & 0.972 & 0.975 & 0.996 & 0.997 & 0.994 & 0.997 \\
-+
+    & 0.3 & 0.955 & 0.966 & 0.984 & 0.992 & 0.978 & 0.991 \\
-+
+GSC &     & 0.964 & 0.731 & 0.984 & 0.704 & 1.005 & 0.770 \\
-+
+RF  &     & 0.968 & 0.605 & 0.974 & 0.541 & 1.014 & 0.614 \\
-+
+\midrule
-+
+PC  &     & 1.354 & 0.000 & -1.190 & 0.000 & -0.314 & 0.000 \\
-+
+SP  & 0.0 & 0.720 & 0.338 & 0.808 & 0.523 & 0.860 & 0.643 \\
-+
+    & 0.1 & 0.839 & 0.583 & 0.955 & 0.882 & 0.992 & 0.978 \\
-+
+    & 0.2 & 0.995 & 0.987 & 0.999 & 1.000 & 0.999 & 1.000 \\
-+
+    & 0.3 & 0.995 & 0.996 & 0.996 & 0.998 & 0.996 & 0.998 \\
-+
+SF  & 0.0 & 0.998 & 0.994 & 1.000 & 0.998 & 1.000 & 0.999 \\
-+
+    & 0.1 & 0.997 & 0.994 & 1.000 & 0.999 & 1.000 & 1.000 \\
-+
+    & 0.2 & 0.999 & 0.998 & 0.999 & 1.000 & 0.999 & 1.000 \\
-+
+    & 0.3 & 0.995 & 0.996 & 0.996 & 0.998 & 0.996 & 0.998 \\
-+
+GSC &     & 0.995 & 0.990 & 0.998 & 0.997 & 0.998 & 0.996 \\
-+
+RF  &     & 0.998 & 0.993 & 0.999 & 0.998 & 0.999 & 0.996 \\
-+
+\midrule
-+
+PC  &     & 2.437 & 0.000 & -1.872 & 0.000 & 2.138 & 0.000 \\
-+
+SP  & 0.0 & 0.838 & 0.525 & 0.901 & 0.686 & 0.932 & 0.779 \\
-+
+    & 0.1 & 0.914 & 0.733 & 0.979 & 0.932 & 0.995 & 0.987 \\
-+
+    & 0.2 & 0.977 & 0.969 & 0.988 & 0.990 & 0.989 & 0.990 \\
-+
+    & 0.3 & 0.952 & 0.950 & 0.964 & 0.971 & 0.965 & 0.970 \\
-+
+SF  & 0.0 & 0.969 & 0.977 & 0.987 & 0.996 & 0.993 & 0.998 \\
-+
+    & 0.1 & 0.975 & 0.978 & 0.993 & 0.996 & 0.997 & 0.998 \\
-+
+    & 0.2 & 0.976 & 0.973 & 0.988 & 0.990 & 0.989 & 0.990 \\
-+
+    & 0.3 & 0.952 & 0.950 & 0.964 & 0.971 & 0.965 & 0.970 \\
-+
+GSC &     & 0.980 & 0.959 & 0.990 & 0.983 & 0.992 & 0.989 \\
-+
+RF  &     & 0.984 & 0.975 & 0.996 & 0.995 & 0.998 & 0.998 \\
-+
+\bottomrule
-+
+\end{tabular}
-+
+\label{tab:solnWeak}
-+
+\end{table}
-+
-+
+\begin{table}[htbp]
-+
+\centering
-+
+\caption{VARIANCE RESULTS FROM GAUSSIAN FITS TO ANGULAR
-+
+DISTRIBUTIONS OF THE FORCE AND TORQUE VECTORS IN THE WEAK SODIUM
-+
+CHLORIDE SOLUTION SYSTEM}
-+
-+
+\footnotesize
-+
+\begin{tabular}{@{} ccrrrrrr @{}}
-+
+\toprule
-+
+\toprule
-+
+& & \multicolumn{3}{c}{Force $\sigma^2$} & \multicolumn{3}{c}{Torque $\sigma^2$} \\
-+
+\cmidrule(lr){3-5}
-+
+\cmidrule(l){6-8}
-+
+Method & $\alpha$ & 9\AA & 12\AA & 15\AA & 9\AA & 12\AA & 15\AA \\
-+
+\midrule
-+
+PC  &     & 882.863 & 510.435 & 344.201 & 277.691 & 154.231 & 100.131 \\
-+
+SP  & 0.0 & 732.569 & 405.704 & 257.756 & 261.445 & 142.245 & 91.497 \\
-+
+    & 0.1 & 329.031 & 70.746 & 12.014 & 118.496 & 25.218 & 4.711 \\
-+
+    & 0.2 & 6.772 & 0.153 & 0.118 & 9.780 & 2.101 & 2.102 \\
-+
+    & 0.3 & 0.951 & 0.774 & 0.784 & 12.108 & 7.673 & 7.851 \\
-+
+SF  & 0.0 & 2.555 & 0.762 & 0.313 & 6.590 & 1.328 & 0.558 \\
-+
+    & 0.1 & 2.561 & 0.560 & 0.123 & 6.464 & 1.162 & 0.457 \\
-+
+    & 0.2 & 0.501 & 0.118 & 0.118 & 5.698 & 2.074 & 2.099 \\
-+
+    & 0.3 & 0.943 & 0.774 & 0.784 & 12.118 & 7.674 & 7.851 \\
-+
+GSC &     & 2.915 & 0.643 & 0.261 & 9.576 & 3.133 & 1.812 \\
-+
+RF  &     & 2.415 & 0.452 & 0.130 & 6.915 & 1.423 & 0.507 \\
-+
+\midrule
-+
+GSSP  & 0.0 & 2.915 & 0.643 & 0.261 & 9.576 & 3.133 & 1.812 \\
-+
+      & 0.1 & 2.251 & 0.324 & 0.064 & 7.628 & 1.639 & 0.497 \\
-+
+      & 0.2 & 0.590 & 0.118 & 0.116 & 6.080 & 2.096 & 2.103 \\
-+
+      & 0.3 & 0.953 & 0.759 & 0.780 & 12.347 & 7.683 & 7.849 \\
-+
+GSSF  & 0.0 & 1.541 & 0.301 & 0.096 & 6.407 & 1.316 & 0.496 \\
-+
+      & 0.1 & 1.541 & 0.237 & 0.050 & 6.356 & 1.202 & 0.457 \\
-+
+      & 0.2 & 0.568 & 0.118 & 0.116 & 6.166 & 2.105 & 2.105 \\
-+
+      & 0.3 & 0.954 & 0.759 & 0.780 & 12.337 & 7.684 & 7.849 \\
-+
+\bottomrule
-+
+\end{tabular}
-+
+\label{tab:solnWeakAng}
-+
+\end{table}
-+
-+
+Because this system is a perturbation of the pure liquid water system,
-+
+comparisons are best drawn between these two sets. The {\sc sp} and
-+
+{\sc sf} methods are not significantly affected by the inclusion of a
-+
+few ions. The aspect of cutoff sphere neutralization aids in the
-+
+smooth incorporation of these ions; thus, all of the observations
-+
+regarding these methods carry over from section
-+
+\ref{sec:WaterResults}. The differences between these systems are more
-+
+visible for the {\sc rf} method. Though good force agreement is still
-+
+maintained, the energy gaps show a significant increase in the scatter
-+
+of the data.
-+
-+
+\subsection{1.1M NaCl Solution Results}
-+
-+
+The bridging of the charged atomic and neutral molecular systems was
-+
+further developed by considering a high ionic strength system
-+
+consisting of 40 ions in the 1000 SPC/E water solvent ($\approx$1.1
-+
+M). The results for the energy gap comparisons and the force and
-+
+torque vector magnitude comparisons are shown in table
-+
+\ref{tab:solnStr}.  The force and torque vector directionality
-+
+results are displayed separately in table \ref{tab:solnStrAng}, where
-+
+the effect of group-based cutoffs and switching functions on the {\sc
-+
+sp} and {\sc sf} potentials are investigated.
-+
-+
+\begin{table}[htbp]
-+
+\centering
-+
+\caption{REGRESSION RESULTS OF THE STRONG SODIUM CHLORIDE SOLUTION
-+
+SYSTEM FOR $\Delta E$ VALUES ({\it upper}), FORCE VECTOR MAGNITUDES
-+
+({\it middle}) AND TORQUE VECTOR MAGNITUDES ({\it lower})}
-+
-+
+\footnotesize
-+
+\begin{tabular}{@{} ccrrrrrr @{}}
-+
+\toprule
-+
+\toprule
-+
+& & \multicolumn{2}{c}{9\AA} & \multicolumn{2}{c}{12\AA} & \multicolumn{2}{c}{15\AA}\\
-+
+\cmidrule(lr){3-4}
-+
+\cmidrule(lr){5-6}
-+
+\cmidrule(l){7-8}
-+
+Method & $\alpha$ & slope & $R^2$ & slope & $R^2$ & slope & $R^2$ \\
-+
+\midrule
-+
+PC  &     & -0.081 & 0.000 & 0.945 & 0.001 & 0.073 & 0.000 \\
-+
+SP  & 0.0 & 0.978 & 0.469 & 0.996 & 0.672 & 0.975 & 0.668 \\
-+
+    & 0.1 & 0.944 & 0.645 & 0.997 & 0.886 & 0.991 & 0.978 \\
-+
+    & 0.2 & 0.873 & 0.896 & 0.985 & 0.993 & 0.980 & 0.993 \\
-+
+    & 0.3 & 0.831 & 0.860 & 0.960 & 0.979 & 0.955 & 0.977 \\
-+
+SF  & 0.0 & 0.858 & 0.905 & 0.985 & 0.970 & 0.990 & 0.998 \\
-+
+    & 0.1 & 0.865 & 0.907 & 0.992 & 0.974 & 0.994 & 0.999 \\
-+
+    & 0.2 & 0.862 & 0.894 & 0.985 & 0.993 & 0.980 & 0.993 \\
-+
+    & 0.3 & 0.831 & 0.859 & 0.960 & 0.979 & 0.955 & 0.977 \\
-+
+GSC &     & 1.985 & 0.152 & 0.760 & 0.031 & 1.106 & 0.062 \\
-+
+RF  &     & 2.414 & 0.116 & 0.813 & 0.017 & 1.434 & 0.047 \\
-+
+\midrule
-+
+PC  &     & -7.028 & 0.000 & -9.364 & 0.000 & 0.925 & 0.865 \\
-+
+SP  & 0.0 & 0.701 & 0.319 & 0.909 & 0.773 & 0.861 & 0.665 \\
-+
+    & 0.1 & 0.824 & 0.565 & 0.970 & 0.930 & 0.990 & 0.979 \\
-+
+    & 0.2 & 0.988 & 0.981 & 0.995 & 0.998 & 0.991 & 0.998 \\
-+
+    & 0.3 & 0.983 & 0.985 & 0.985 & 0.991 & 0.978 & 0.990 \\
-+
+SF  & 0.0 & 0.993 & 0.988 & 0.992 & 0.984 & 0.998 & 0.999 \\
-+
+    & 0.1 & 0.993 & 0.989 & 0.993 & 0.986 & 0.998 & 1.000 \\
-+
+    & 0.2 & 0.993 & 0.992 & 0.995 & 0.998 & 0.991 & 0.998 \\
-+
+    & 0.3 & 0.983 & 0.985 & 0.985 & 0.991 & 0.978 & 0.990 \\
-+
+GSC &     & 0.964 & 0.897 & 0.970 & 0.917 & 0.925 & 0.865 \\
-+
+RF  &     & 0.994 & 0.864 & 0.988 & 0.865 & 0.980 & 0.784 \\
-+
+\midrule
-+
+PC  &     & -2.212 & 0.000 & -0.588 & 0.000 & 0.953 & 0.925 \\
-+
+SP  & 0.0 & 0.800 & 0.479 & 0.930 & 0.804 & 0.924 & 0.759 \\
-+
+    & 0.1 & 0.883 & 0.694 & 0.976 & 0.942 & 0.993 & 0.986 \\
-+
+    & 0.2 & 0.952 & 0.943 & 0.980 & 0.984 & 0.980 & 0.983 \\
-+
+    & 0.3 & 0.914 & 0.909 & 0.943 & 0.948 & 0.944 & 0.946 \\
-+
+SF  & 0.0 & 0.945 & 0.953 & 0.980 & 0.984 & 0.991 & 0.998 \\
-+
+    & 0.1 & 0.951 & 0.954 & 0.987 & 0.986 & 0.995 & 0.998 \\
-+
+    & 0.2 & 0.951 & 0.946 & 0.980 & 0.984 & 0.980 & 0.983 \\
-+
+    & 0.3 & 0.914 & 0.908 & 0.943 & 0.948 & 0.944 & 0.946 \\
-+
+GSC &     & 0.882 & 0.818 & 0.939 & 0.902 & 0.953 & 0.925 \\
-+
+RF  &     & 0.949 & 0.939 & 0.988 & 0.988 & 0.992 & 0.993 \\
-+
+\bottomrule
-+
+\end{tabular}
-+
+\label{tab:solnStr}
-+
+\end{table}
-+
-+
+\begin{table}[htbp]
-+
+\centering
-+
+\caption{VARIANCE RESULTS FROM GAUSSIAN FITS TO ANGULAR DISTRIBUTIONS
-+
+OF THE FORCE AND TORQUE VECTORS IN THE STRONG SODIUM CHLORIDE SOLUTION
-+
+SYSTEM}
-+
-+
+\footnotesize
-+
+\begin{tabular}{@{} ccrrrrrr @{}}
-+
+\toprule
-+
+\toprule
-+
+& & \multicolumn{3}{c}{Force $\sigma^2$} & \multicolumn{3}{c}{Torque $\sigma^2$} \\
-+
+\cmidrule(lr){3-5}
-+
+\cmidrule(l){6-8}
-+
+Method & $\alpha$ & 9\AA & 12\AA & 15\AA & 9\AA & 12\AA & 15\AA \\
-+
+\midrule
-+
+PC  &     & 957.784 & 513.373 & 2.260 & 340.043 & 179.443 & 13.079 \\
-+
+SP  & 0.0 & 786.244 & 139.985 & 259.289 & 311.519 & 90.280 & 105.187 \\
-+
+    & 0.1 & 354.697 & 38.614 & 12.274 & 144.531 & 23.787 & 5.401 \\
-+
+    & 0.2 & 7.674 & 0.363 & 0.215 & 16.655 & 3.601 & 3.634 \\
-+
+    & 0.3 & 1.745 & 1.456 & 1.449 & 23.669 & 14.376 & 14.240 \\
-+
+SF  & 0.0 & 3.282 & 8.567 & 0.369 & 11.904 & 6.589 & 0.717 \\
-+
+    & 0.1 & 3.263 & 7.479 & 0.142 & 11.634 & 5.750 & 0.591 \\
-+
+    & 0.2 & 0.686 & 0.324 & 0.215 & 10.809 & 3.580 & 3.635 \\
-+
+    & 0.3 & 1.749 & 1.456 & 1.449 & 23.635 & 14.375 & 14.240 \\
-+
+GSC &     & 6.181 & 2.904 & 2.263 & 44.349 & 19.442 & 12.873 \\
-+
+RF  &     & 3.891 & 0.847 & 0.323 & 18.628 & 3.995 & 2.072 \\
-+
+\midrule
-+
+GSSP  & 0.0 & 6.197 & 2.929 & 2.290 & 44.441 & 19.442 & 12.873 \\
-+
+      & 0.1 & 4.688 & 1.064 & 0.260 & 31.208 & 6.967 & 2.303 \\
-+
+      & 0.2 & 1.021 & 0.218 & 0.213 & 14.425 & 3.629 & 3.649 \\
-+
+      & 0.3 & 1.752 & 1.454 & 1.451 & 23.540 & 14.390 & 14.245 \\
-+
+GSSF  & 0.0 & 2.494 & 0.546 & 0.217 & 16.391 & 3.230 & 1.613 \\
-+
+      & 0.1 & 2.448 & 0.429 & 0.106 & 16.390 & 2.827 & 1.159 \\
-+
+      & 0.2 & 0.899 & 0.214 & 0.213 & 13.542 & 3.583 & 3.645 \\
-+
+      & 0.3 & 1.752 & 1.454 & 1.451 & 23.587 & 14.390 & 14.245 \\
-+
+\bottomrule
-+
+\end{tabular}
-+
+\label{tab:solnStrAng}
-+
+\end{table}
-+
-+
+The {\sc rf} method struggles with the jump in ionic strength. The
-+
+configuration energy differences degrade to unusable levels while the
-+
+forces and torques show a more modest reduction in the agreement with
-+
+{\sc spme}. The {\sc rf} method was designed for homogeneous systems,
-+
+and this attribute is apparent in these results.
-+
-+
+The {\sc sp} and {\sc sf} methods require larger cutoffs to maintain
-+
+their agreement with {\sc spme}. With these results, we still
-+
+recommend undamped to moderate damping for the {\sc sf} method and
-+
+moderate damping for the {\sc sp} method, both with cutoffs greater
-+
+than 12\AA.
-+
+\subsection{6\AA\ Argon Sphere in SPC/E Water Results}
-<
+\section{Short-Time Dynamics: Velocity Autocorrelation Functions of NaCl Crystals}
->
+The final model system studied was a 6\AA\ sphere of Argon solvated
->
+by SPC/E water. This serves as a test case of a specifically sized
->
+electrostatic defect in a disordered molecular system. The results for
->
+the energy gap comparisons and the force and torque vector magnitude
->
+comparisons are shown in table \ref{tab:argon}.  The force and torque
->
+vector directionality results are displayed separately in table
->
+\ref{tab:argonAng}, where the effect of group-based cutoffs and
->
+switching functions on the {\sc sp} and {\sc sf} potentials are
->
+investigated.
-+
+\begin{table}[htbp]
-+
+\centering
-+
+\caption{REGRESSION RESULTS OF THE 6\AA\ ARGON SPHERE IN LIQUID
-+
+WATER SYSTEM FOR $\Delta E$ VALUES ({\it upper}), FORCE VECTOR
-+
+MAGNITUDES ({\it middle}) AND TORQUE VECTOR MAGNITUDES ({\it lower})}
-+
-+
+\footnotesize
-+
+\begin{tabular}{@{} ccrrrrrr @{}}
-+
+\toprule
-+
+\toprule
-+
+& & \multicolumn{2}{c}{9\AA} & \multicolumn{2}{c}{12\AA} & \multicolumn{2}{c}{15\AA}\\
-+
+\cmidrule(lr){3-4}
-+
+\cmidrule(lr){5-6}
-+
+\cmidrule(l){7-8}
-+
+Method & $\alpha$ & slope & $R^2$ & slope & $R^2$ & slope & $R^2$ \\
-+
+\midrule
-+
+PC  &     & 2.320 & 0.008 & -0.650 & 0.001 & 3.848 & 0.029 \\
-+
+SP  & 0.0 & 1.053 & 0.711 & 0.977 & 0.820 & 0.974 & 0.882 \\
-+
+    & 0.1 & 1.032 & 0.846 & 0.989 & 0.965 & 0.992 & 0.994 \\
-+
+    & 0.2 & 0.993 & 0.995 & 0.982 & 0.998 & 0.986 & 0.998 \\
-+
+    & 0.3 & 0.968 & 0.995 & 0.954 & 0.992 & 0.961 & 0.994 \\
-+
+SF  & 0.0 & 0.982 & 0.996 & 0.992 & 0.999 & 0.993 & 1.000 \\
-+
+    & 0.1 & 0.987 & 0.996 & 0.996 & 0.999 & 0.997 & 1.000 \\
-+
+    & 0.2 & 0.989 & 0.998 & 0.984 & 0.998 & 0.989 & 0.998 \\
-+
+    & 0.3 & 0.971 & 0.995 & 0.957 & 0.992 & 0.965 & 0.994 \\
-+
+GSC &     & 1.002 & 0.983 & 0.992 & 0.973 & 0.996 & 0.971 \\
-+
+RF  &     & 0.998 & 0.995 & 0.999 & 0.998 & 0.998 & 0.998 \\
-+
+\midrule
-+
+PC  &     & -36.559 & 0.002 & -44.917 & 0.004 & -52.945 & 0.006 \\
-+
+SP  & 0.0 & 0.890 & 0.786 & 0.927 & 0.867 & 0.949 & 0.909 \\
-+
+    & 0.1 & 0.942 & 0.895 & 0.984 & 0.974 & 0.997 & 0.995 \\
-+
+    & 0.2 & 0.999 & 0.997 & 1.000 & 1.000 & 1.000 & 1.000 \\
-+
+    & 0.3 & 1.001 & 0.999 & 1.001 & 1.000 & 1.001 & 1.000 \\
-+
+SF  & 0.0 & 1.000 & 0.999 & 1.000 & 1.000 & 1.000 & 1.000 \\
-+
+    & 0.1 & 1.000 & 0.999 & 1.000 & 1.000 & 1.000 & 1.000 \\
-+
+    & 0.2 & 1.000 & 1.000 & 1.000 & 1.000 & 1.000 & 1.000 \\
-+
+    & 0.3 & 1.001 & 0.999 & 1.001 & 1.000 & 1.001 & 1.000 \\
-+
+GSC &     & 0.999 & 0.999 & 1.000 & 1.000 & 1.000 & 1.000 \\
-+
+RF  &     & 0.999 & 0.999 & 1.000 & 1.000 & 1.000 & 1.000 \\
-+
+\midrule
-+
+PC  &     & 1.984 & 0.000 & 0.012 & 0.000 & 1.357 & 0.000 \\
-+
+SP  & 0.0 & 0.850 & 0.552 & 0.907 & 0.703 & 0.938 & 0.793 \\
-+
+    & 0.1 & 0.924 & 0.755 & 0.980 & 0.936 & 0.995 & 0.988 \\
-+
+    & 0.2 & 0.985 & 0.983 & 0.986 & 0.988 & 0.987 & 0.988 \\
-+
+    & 0.3 & 0.961 & 0.966 & 0.959 & 0.964 & 0.960 & 0.966 \\
-+
+SF  & 0.0 & 0.977 & 0.989 & 0.987 & 0.995 & 0.992 & 0.998 \\
-+
+    & 0.1 & 0.982 & 0.989 & 0.992 & 0.996 & 0.997 & 0.998 \\
-+
+    & 0.2 & 0.984 & 0.987 & 0.986 & 0.987 & 0.987 & 0.988 \\
-+
+    & 0.3 & 0.961 & 0.966 & 0.959 & 0.964 & 0.960 & 0.966 \\
-+
+GSC &     & 0.995 & 0.981 & 0.999 & 0.990 & 1.000 & 0.993 \\
-+
+RF  &     & 0.993 & 0.988 & 0.997 & 0.995 & 0.999 & 0.998 \\
-+
+\bottomrule
-+
+\end{tabular}
-+
+\label{tab:argon}
-+
+\end{table}
-+
-+
+\begin{table}[htbp]
-+
+\centering
-+
+\caption{VARIANCE RESULTS FROM GAUSSIAN FITS TO ANGULAR
-+
+DISTRIBUTIONS OF THE FORCE AND TORQUE VECTORS IN THE 6\AA\ SPHERE OF
-+
+ARGON IN LIQUID WATER SYSTEM}
-+
-+
+\footnotesize
-+
+\begin{tabular}{@{} ccrrrrrr @{}}
-+
+\toprule
-+
+\toprule
-+
+& & \multicolumn{3}{c}{Force $\sigma^2$} & \multicolumn{3}{c}{Torque $\sigma^2$} \\
-+
+\cmidrule(lr){3-5}
-+
+\cmidrule(l){6-8}
-+
+Method & $\alpha$ & 9\AA & 12\AA & 15\AA & 9\AA & 12\AA & 15\AA \\
-+
+\midrule
-+
+PC  &     & 568.025 & 265.993 & 195.099 & 246.626 & 138.600 & 91.654 \\
-+
+SP  & 0.0 & 504.578 & 251.694 & 179.932 & 231.568 & 131.444 & 85.119 \\
-+
+    & 0.1 & 224.886 & 49.746 & 9.346 & 104.482 & 23.683 & 4.480 \\
-+
+    & 0.2 & 4.889 & 0.197 & 0.155 & 6.029 & 2.507 & 2.269 \\
-+
+    & 0.3 & 0.817 & 0.833 & 0.812 & 8.286 & 8.436 & 8.135 \\
-+
+SF  & 0.0 & 1.924 & 0.675 & 0.304 & 3.658 & 1.448 & 0.600 \\
-+
+    & 0.1 & 1.937 & 0.515 & 0.143 & 3.565 & 1.308 & 0.546 \\
-+
+    & 0.2 & 0.407 & 0.166 & 0.156 & 3.086 & 2.501 & 2.274 \\
-+
+    & 0.3 & 0.815 & 0.833 & 0.812 & 8.330 & 8.437 & 8.135 \\
-+
+GSC &     & 2.098 & 0.584 & 0.284 & 5.391 & 2.414 & 1.501 \\
-+
+RF  &     & 1.822 & 0.408 & 0.142 & 3.799 & 1.362 & 0.550 \\
-+
+\midrule
-+
+GSSP  & 0.0 & 2.098 & 0.584 & 0.284 & 5.391 & 2.414 & 1.501 \\
-+
+      & 0.1 & 1.652 & 0.309 & 0.087 & 4.197 & 1.401 & 0.590 \\
-+
+      & 0.2 & 0.465 & 0.165 & 0.153 & 3.323 & 2.529 & 2.273 \\
-+
+      & 0.3 & 0.813 & 0.825 & 0.816 & 8.316 & 8.447 & 8.132 \\
-+
+GSSF  & 0.0 & 1.173 & 0.292 & 0.113 & 3.452 & 1.347 & 0.583 \\
-+
+      & 0.1 & 1.166 & 0.240 & 0.076 & 3.381 & 1.281 & 0.575 \\
-+
+      & 0.2 & 0.459 & 0.165 & 0.153 & 3.430 & 2.542 & 2.273 \\
-+
+      & 0.3 & 0.814 & 0.825 & 0.816 & 8.325 & 8.447 & 8.132 \\
-+
+\bottomrule
-+
+\end{tabular}
-+
+\label{tab:argonAng}
-+
+\end{table}
-+
-+
+This system does not appear to show any significant deviations from
-+
+the previously observed results. The {\sc sp} and {\sc sf} methods
-+
+have agreements similar to those observed in section
-+
+\ref{sec:WaterResults}. The only significant difference is the
-+
+improvement in the configuration energy differences for the {\sc rf}
-+
+method. This is surprising in that we are introducing an inhomogeneity
-+
+to the system; however, this inhomogeneity is charge-neutral and does
-+
+not result in charged cutoff spheres. The charge-neutrality of the
-+
+cutoff spheres, which the {\sc sp} and {\sc sf} methods explicitly
-+
+enforce, seems to play a greater role in the stability of the {\sc rf}
-+
+method than the required homogeneity of the environment.
-+
-+
-+
+\section{Short-Time Dynamics: Velocity Autocorrelation Functions of NaCl Crystals}\label{sec:ShortTimeDynamics}
-+
+Zahn {\it et al.} investigated the structure and dynamics of water
+using eqs. (\ref{eq:ZahnPot}) and
+(\ref{eq:WolfForces}).\cite{Zahn02,Kast03} Their results indicated
+\centering
+\includegraphics[width = \linewidth]{./figures/vCorrPlot.pdf}
+\caption{Velocity autocorrelation functions of NaCl crystals at
-<
+K using {\sc spme}, {\sc sf} ($\alpha$ = 0.0, 0.1, \& 0.2), and {\sc
-<
+sp} ($\alpha$ = 0.2). The inset is a magnification of the area around
-<
+the first minimum.  The times to first collision are nearly identical,
-<
+but differences can be seen in the peaks and troughs, where the
-<
+undamped and weakly damped methods are stiffer than the moderately
-<
+damped and {\sc spme} methods.}
->
+K using {\sc spme}, {\sc sf} ($\alpha$ = 0.0, 0.1, \&
->
+.2\AA$^{-1}$), and {\sc sp} ($\alpha$ = 0.2\AA$^{-1}$). The inset is
->
+a magnification of the area around the first minimum.  The times to
->
+first collision are nearly identical, but differences can be seen in
->
+the peaks and troughs, where the undamped and weakly damped methods
->
+are stiffer than the moderately damped and {\sc spme} methods.}
+\label{fig:vCorrPlot}
+\end{figure}
+constructed out of the damped electrostatic interaction are less
+important.
-<
+\section{Collective Motion: Power Spectra of NaCl Crystals}
->
+\section{Collective Motion: Power Spectra of NaCl Crystals}\label{sec:LongTimeDynamics}
+To evaluate how the differences between the methods affect the
+collective long-time motion, we computed power spectra from long-time
+\includegraphics[width = \linewidth]{./figures/spectraSquare.pdf}
+\caption{Power spectra obtained from the velocity auto-correlation
+functions of NaCl crystals at 1000K while using {\sc spme}, {\sc sf}
-<
+($\alpha$ = 0, 0.1, \& 0.2), and {\sc sp} ($\alpha$ = 0.2).  The inset
-<
+shows the frequency region below 100 cm$^{-1}$ to highlight where the
-<
+spectra differ.}
->
+($\alpha$ = 0, 0.1, \& 0.2\AA$^{-1}$), and {\sc sp} ($\alpha$ =
->
+.2\AA$^{-1}$).  The inset shows the frequency region below 100
->
+cm$^{-1}$ to highlight where the spectra differ.}
+\label{fig:methodPS}
+\end{figure}
+the NaCl crystal at 1000K.  The undamped shifted force ({\sc sf})
+method is off by less than 10 cm$^{-1}$, and increasing the
+electrostatic damping to 0.25\AA$^{-1}$ gives quantitative agreement
-<
+with the power spectrum obtained using the Ewald sum.  Overdamping can
->
+with the power spectrum obtained using the Ewald sum.  Over-damping can
+result in underestimates of frequencies of the long-wavelength
+motions.}
+\label{fig:dampInc}
-+
+\end{figure}
-+
-+
+\section{An Application: TIP5P-E Water}\label{sec:t5peApplied}
-+
-+
+The above sections focused on the energetics and dynamics of a variety
-+
+of systems when utilizing the {\sc sp} and {\sc sf} pairwise
-+
+techniques.  A unitary correlation with results obtained using the
-+
+Ewald summation should result in a successful reproduction of both the
-+
+static and dynamic properties of any selected system.  To test this,
-+
+we decided to calculate a series of properties for the TIP5P-E water
-+
+model when using the {\sc sf} technique.
-+
-+
+The TIP5P-E water model is a variant of Mahoney and Jorgensen's
-+
+five-point transferable intermolecular potential (TIP5P) model for
-+
+water.\cite{Mahoney00} TIP5P was developed to reproduce the density
-+
+maximum anomaly present in liquid water near 4$^\circ$C. As with many
-+
+previous point charge water models (such as ST2, TIP3P, TIP4P, SPC,
-+
+and SPC/E), TIP5P was parametrized using a simple cutoff with no
-+
+long-range electrostatic
-+
+correction.\cite{Stillinger74,Jorgensen83,Berendsen81,Berendsen87}
-+
+Without this correction, the pressure term on the central particle
-+
+from the surroundings is missing. Because they expand to compensate
-+
+for this added pressure term when this correction is included, systems
-+
+composed of these particles tend to underpredict the density of water
-+
+under standard conditions. When using any form of long-range
-+
+electrostatic correction, it has become common practice to develop or
-+
+utilize a reparametrized water model that corrects for this
-+
+effect.\cite{vanderSpoel98,Fennell04,Horn04} The TIP5P-E model follows
-+
+this practice and was optimized specifically for use with the Ewald
-+
+summation.\cite{Rick04} In his publication, Rick preserved the
-+
+geometry and point charge magnitudes in TIP5P and focused on altering
-+
+the Lennard-Jones parameters to correct the density at
-+
+K.\cite{Rick04} With the density corrected, he compared common
-+
+water properties for TIP5P-E using the Ewald sum with TIP5P using a
-+
+\AA\ cutoff.
-+
-+
+In the following sections, we compared these same water properties
-+
+calculated from TIP5P-E using the Ewald sum with TIP5P-E using the
-+
+{\sc sf} technique.  In the above evaluation of the pairwise
-+
+techniques, we observed some flexibility in the choice of parameters.
-+
+Because of this, the following comparisons include the {\sc sf}
-+
+technique with a 12\AA\ cutoff and an $\alpha$ = 0.0, 0.1, and
-+
+.2\AA$^{-1}$, as well as a 9\AA\ cutoff with an $\alpha$ =
-+
+.2\AA$^{-1}$.
-+
-+
+\subsection{Density}\label{sec:t5peDensity}
-+
-+
+As stated previously, the property that prompted the development of
-+
+TIP5P-E was the density at 1 atm.  The density depends upon the
-+
+internal pressure of the system in the $NPT$ ensemble, and the
-+
+calculation of the pressure includes a components from both the
-+
+kinetic energy and the virial. More specifically, the instantaneous
-+
+molecular pressure ($p(t)$) is given by
-+
+\begin{equation}
-+
+p(t) =  \frac{1}{\textrm{d}V}\sum_\mu
-+
+        \left[\frac{\mathbf{P}_{\mu}^2}{M_{\mu}}
-+
+        + \mathbf{R}_{\mu}\cdot\sum_i\mathbf{F}_{\mu i}\right],
-+
+\label{eq:MolecularPressure}
-+
+\end{equation}
-+
+where $V$ is the volume, $\mathbf{P}_{\mu}$ is the momentum of
-+
+molecule $\mu$, $\mathbf{R}_\mu$ is the position of the center of mass
-+
+($M_\mu$) of molecule $\mu$, and $\mathbf{F}_{\mu i}$ is the force on
-+
+atom $i$ of molecule $\mu$.\cite{Melchionna93} The virial term (the
-+
+right term in the brackets of eq. \ref{eq:MolecularPressure}) is
-+
+directly dependent on the interatomic forces.  Since the {\sc sp}
-+
+method does not modify the forces (see
-+
+sec. \ref{sec:PairwiseDerivation}), the pressure using {\sc sp} will
-+
+be identical to that obtained without an electrostatic correction.
-+
+The {\sc sf} method does alter the virial component and, by way of the
-+
+modified pressures, should provide densities more in line with those
-+
+obtained using the Ewald summation.
-+
-+
+To compare densities, $NPT$ simulations were performed with the same
-+
+temperatures as those selected by Rick in his Ewald summation
-+
+simulations.\cite{Rick04} In order to improve statistics around the
-+
+density maximum, 3ns trajectories were accumulated at 0, 12.5, and
-+
+$^\circ$C, while 2ns trajectories were obtained at all other
-+
+temperatures. The average densities were calculated from the later
-+
+three-fourths of each trajectory. Similar to Mahoney and Jorgensen's
-+
+method for accumulating statistics, these sequences were spliced into
-+
+segments to calculate the average density and standard deviation
-+
+at each temperature.\cite{Mahoney00}
-+
-+
+\begin{figure}
-+
+\includegraphics[width=\linewidth]{./figures/tip5peDensities.pdf}
-+
+\caption{Density versus temperature for the TIP5P-E water model when
-+
+using the Ewald summation (Ref. \cite{Rick04}) and the {\sc sf} method
-+
+with various parameters. The pressure term from the image-charge shell
-+
+is larger than that provided by the reciprocal-space portion of the
-+
+Ewald summation, leading to slightly lower densities. This effect is
-+
+more visible with the 9\AA\ cutoff, where the image charges exert a
-+
+greater force on the central particle. The error bars for the {\sc sf}
-+
+methods show plus or minus the standard deviation of the density
-+
+measurement at each temperature.}
-+
+\label{fig:t5peDensities}
+\end{figure}
-+
+Figure \ref{fig:t5peDensities} shows the densities calculated for
-+
+TIP5P-E using differing electrostatic corrections overlaid on the
-+
+experimental values.\cite{CRC80} The densities when using the {\sc sf}
-+
+technique are close to, though typically lower than, those calculated
-+
+while using the Ewald summation. These slightly reduced densities
-+
+indicate that the pressure component from the image charges at
-+
+R$_\textrm{c}$ is larger than that exerted by the reciprocal-space
-+
+portion of the Ewald summation. Bringing the image charges closer to
-+
+the central particle by choosing a 9\AA\ R$_\textrm{c}$ (rather than
-+
+the preferred 12\AA\ R$_\textrm{c}$) increases the strength of their
-+
+interactions, resulting in a further reduction of the densities.
-+
+Because the strength of the image charge interactions has a noticable
-+
+effect on the density, we would expect the use of electrostatic
-+
+damping to also play a role in these calculations. Larger values of
-+
+$\alpha$ weaken the pair-interactions; and since electrostatic damping
-+
+is distance-dependent, force components from the image charges will be
-+
+reduced more than those from particles close the the central
-+
+charge. This effect is visible in figure \ref{fig:t5peDensities} with
-+
+the damped {\sc sf} sums showing slightly higher densities; however,
-+
+it is apparent that the choice of cutoff radius plays a much more
-+
+important role in the resulting densities.
-+
-+
+As a final note, all of the above density calculations were performed
-+
+with systems of 512 water molecules. Rick observed a system sized
-+
+dependence of the computed densities when using the Ewald summation,
-+
+most likely due to his tying of the convergence parameter to the box
-+
+dimensions.\cite{Rick04} For systems of 256 water molecules, the
-+
+calculated densities were on average 0.002 to 0.003 g/cm$^3$ lower. A
-+
+system size of 256 molecules would force the use of a shorter
-+
+R$_\textrm{c}$ when using the {\sc sf} technique, and this would also
-+
+lower the densities. Moving to larger systems, as long as the
-+
+R$_\textrm{c}$ remains at a fixed value, we would expect the densities
-+
+to remain constant.
-+
-+
+\subsection{Liquid Structure}\label{sec:t5peLiqStructure}
-+
-+
+A common function considered when developing and comparing water
-+
+models is the oxygen-oxygen radial distribution function
-+
+($g_\textrm{OO}(r)$). The $g_\textrm{OO}(r)$ is the probability of
-+
+finding a pair of oxygen atoms some distance ($r$) apart relative to a
-+
+random distribution at the same density.\cite{Allen87} It is
-+
+calculated via
-+
+\begin{equation}
-+
+g_\textrm{OO}(r) = \frac{V}{N^2}\left\langle\sum_i\sum_{j\ne i}
-+
+        \delta(\mathbf{r}-\mathbf{r}_{ij})\right\rangle,
-+
+\label{eq:GOOofR}
-+
+\end{equation}
-+
+where the double sum is over all $i$ and $j$ pairs of $N$ oxygen
-+
+atoms. The $g_\textrm{OO}(r)$ can be directly compared to X-ray or
-+
+neutron scattering experiments through the oxygen-oxygen structure
-+
+factor ($S_\textrm{OO}(k)$) by the following relationship:
-+
+\begin{equation}
-+
+S_\textrm{OO}(k) = 1 + 4\pi\rho
-+
+        \int_0^\infty r^2\frac{\sin kr}{kr}g_\textrm{OO}(r)\textrm{d}r.
-+
+\label{eq:SOOofK}
-+
+\end{equation}
-+
+Thus, $S_\textrm{OO}(k)$ is related to the Fourier-Laplace transform
-+
+of $g_\textrm{OO}(r)$.
-+
-+
+The expermentally determined $g_\textrm{OO}(r)$ for liquid water has
-+
+been compared in great detail with the various common water models,
-+
+and TIP5P was found to be in better agreement than other rigid,
-+
+non-polarizable models.\cite{Sorenson00} This excellent agreement with
-+
+experiment was maintained when Rick developed TIP5P-E.\cite{Rick04} To
-+
+check whether the choice of using the Ewald summation or the {\sc sf}
-+
+technique alters the liquid structure, the $g_\textrm{OO}(r)$s at 298K
-+
+and 1atm were determined for the systems compared in the previous
-+
+section.
-+
-+
+\begin{figure}
-+
+\includegraphics[width=\linewidth]{./figures/tip5peGofR.pdf}
-+
+\caption{The $g_\textrm{OO}(r)$s calculated for TIP5P-E at 298K and
-+
+atm while using the Ewald summation (Ref. \cite{Rick04}) and the {\sc
-+
+sf} technique with varying parameters. Even with the reduced densities
-+
+using the {\sc sf} technique, the $g_\textrm{OO}(r)$s are essentially
-+
+identical.}
-+
+\label{fig:t5peGofRs}
-+
+\end{figure}
-+
-+
+The $g_\textrm{OO}(r)$s calculated for TIP5P-E while using the {\sc
-+
+sf} technique with a various parameters are overlaid on the
-+
+$g_\textrm{OO}(r)$ while using the Ewald summation. The differences in
-+
+density do not appear to have any effect on the liquid structure as
-+
+the $g_\textrm{OO}(r)$s are indistinquishable. These results indicate
-+
+that the $g_\textrm{OO}(r)$ is insensitive to the choice of
-+
+electrostatic correction.
-+
-+
+\subsection{Thermodynamic Properties}\label{sec:t5peThermo}
-+
-+
+In addition to the density, there are a variety of thermodynamic
-+
+quantities that can be calculated for water and compared directly to
-+
+experimental values. Some of these additional quatities include the
-+
+latent heat of vaporization ($\Delta H_\textrm{vap}$), the constant
-+
+pressure heat capacity ($C_p$), the isothermal compressibility
-+
+($\kappa_T$), the thermal expansivity ($\alpha_p$), and the static
-+
+dielectric constant ($\epsilon$). All of these properties were
-+
+calculated for TIP5P-E with the Ewald summation, so they provide a
-+
+good set for comparisons involving the {\sc sf} technique.
-+
-+
+The $\Delta H_\textrm{vap}$ is the enthalpy change required to
-+
+transform one mol of substance from the liquid phase to the gas
-+
+phase.\cite{Berry00} In molecular simulations, this quantity can be
-+
+determined via
-+
+\begin{equation}
-+
+\begin{split}
-+
+\Delta H_\textrm{vap} &= H_\textrm{gas} - H_\textrm{liq.} \\
-+
+                      &= E_\textrm{gas} - E_\textrm{liq.}
-+
+                        + p(V_\textrm{gas} - V_\textrm{liq.}) \\
-+
+                      &\approx -\frac{\langle U_\textrm{liq.}\rangle}{N}+ RT,
-+
+\end{split}
-+
+\label{eq:DeltaHVap}
-+
+\end{equation}
-+
+where $E$ is the total energy, $U$ is the potential energy, $p$ is the
-+
+pressure, $V$ is the volume, $N$ is the number of molecules, $R$ is
-+
+the gas constant, and $T$ is the temperature.\cite{Jorgensen98b} As
-+
+seen in the last line of equation \ref{eq:DeltaHVap}, we can
-+
+approximate $\Delta H_\textrm{vap}$ by using the ideal gas for the gas
-+
+state. This allows us to cancel the kinetic energy terms, leaving only
-+
+the $U_\textrm{liq.}$ term. Additionally, the $pV$ term for the gas is
-+
+several orders of magnitude larger than that of the liquid, so we can
-+
+neglect the liquid $pV$ term.
-+
-+
+The remaining thermodynamic properties can all be calculated from
-+
+fluctuations of the enthalpy, volume, and system dipole
-+
+moment.\cite{Allen87} $C_p$ can be calculated from fluctuations in the
-+
+enthalpy in constant pressure simulations via
-+
+\begin{equation}
-+
+\begin{split}
-+
+C_p     = \left(\frac{\partial H}{\partial T}\right)_{N,p}
-+
+        = \frac{(\langle H^2\rangle - \langle H\rangle^2)}{Nk_BT^2},
-+
+\end{split}
-+
+\label{eq:Cp}
-+
+\end{equation}
-+
+where $k_B$ is Boltzmann's constant.  $\kappa_T$ can be calculated via
-+
+\begin{equation}
-+
+\begin{split}
-+
+\kappa_T = \frac{1}{V}\left(\frac{\partial V}{\partial p}\right)_{N,T}
-+
+         = \frac{(\langle V^2\rangle_{N,P,T} - \langle V\rangle^2_{N,P,T})}
-+
+                {k_BT\langle V\rangle_{N,P,T}},
-+
+\end{split}
-+
+\label{eq:kappa}
-+
+\end{equation}
-+
+and $\alpha_p$ can be calculated via
-+
+\begin{equation}
-+
+\begin{split}
-+
+\alpha_p        = \frac{1}{V}\left(\frac{\partial V}{\partial T}\right)_{N,P}
-+
+                = \frac{(\langle VH\rangle_{N,P,T}
-+
+                        - \langle V\rangle_{N,P,T}\langle H\rangle_{N,P,T})}
-+
+                        {k_BT^2\langle V\rangle_{N,P,T}}.
-+
+\end{split}
-+
+\label{eq:alpha}
-+
+\end{equation}
-+
+Using the Ewald sum under tin-foil boundary conditions, $\epsilon$ can
-+
+be calculated for systems of non-polarizable substances via
-+
+\begin{equation}
-+
+\epsilon = 1 + \frac{\langle M^2\rangle}{3\epsilon_0k_\textrm{B}TV},
-+
+\label{eq:staticDielectric}
-+
+\end{equation}
-+
+where $\epsilon_0$ is the permittivity of free space and $\langle
-+
+M^2\rangle$ is the fluctuation of the system dipole
-+
+moment.\cite{Allen87} The numerator in the fractional term in equation
-+
+\ref{eq:staticDielectric} is the fluctuation of the simulation-box
-+
+dipole moment, identical to the quantity calculated in the
-+
+finite-system Kirkwood $g$ factor ($G_k$):
-+
+\begin{equation}
-+
+G_k = \frac{\langle M^2\rangle}{N\mu^2},
-+
+\label{eq:KirkwoodFactor}
-+
+\end{equation}
-+
+where $\mu$ is the dipole moment of a single molecule of the
-+
+homogeneous system.\cite{Neumann83,Neumann84,Neumann85} The
-+
+fluctuation term in both equation \ref{eq:staticDielectric} and
-+
+\ref{eq:KirkwoodFactor} is calculated as follows,
-+
+\begin{equation}
-+
+\begin{split}
-+
+\langle M^2\rangle &= \langle\bm{M}\cdot\bm{M}\rangle
-+
+                        - \langle\bm{M}\rangle\cdot\langle\bm{M}\rangle \\
-+
+                   &= \langle M_x^2+M_y^2+M_z^2\rangle
-+
+                        - (\langle M_x\rangle^2 + \langle M_x\rangle^2
-+
+                                + \langle M_x\rangle^2).
-+
+\end{split}
-+
+\label{eq:fluctBoxDipole}
-+
+\end{equation}
-+
+This fluctuation term can be accumulated during the simulation;
-+
+however, it converges rather slowly, thus requiring multi-nanosecond
-+
+simulation times.\cite{Horn04} In the case of tin-foil boundary
-+
+conditions, the dielectric/surface term of equation \ref{eq:EwaldSum}
-+
+is equal to zero. Since the {\sc sf} method also lacks this
-+
+dielectric/surface term, equation \ref{eq:staticDielectric} is still
-+
+valid for determining static dielectric constants.
-+
-+
+All of the above properties were calculated from the same trajectories
-+
+used to determine the densities in section \ref{sec:t5peDensity}
-+
+except for the static dielectric constants. The $\epsilon$ values were
-+
+accumulated from 2ns $NVE$ ensemble trajectories with system densities
-+
+fixed at the average values from the $NPT$ simulations at each of the
-+
+temperatures. The resulting values are displayed in figure
-+
+\ref{fig:t5peThermo}.
-+
+\begin{figure}
-+
+\centering
-+
+\includegraphics[width=5.5in]{./figures/t5peThermo.pdf}
-+
+\caption{Thermodynamic properties for TIP5P-E using the Ewald summation
-+
+and the {\sc sf} techniques along with the experimental values. Units
-+
+for the properties are kcal mol$^{-1}$ for $\Delta H_\textrm{vap}$,
-+
+cal mol$^{-1}$ K$^{-1}$ for $C_p$, 10$^6$ atm$^{-1}$ for $\kappa_T$,
-+
+and 10$^5$ K$^{-1}$ for $\alpha_p$. Ewald summation results are from
-+
+reference \cite{Rick04}. Experimental values for $\Delta
-+
+H_\textrm{vap}$, $\kappa_T$, and $\alpha_p$ are from reference
-+
+\cite{Kell75}. Experimental values for $C_p$ are from reference
-+
+\cite{Wagner02}. Experimental values for $\epsilon$ are from reference
-+
+\cite{Malmberg56}.}
-+
+\label{fig:t5peThermo}
-+
+\end{figure}
-+
-+
+As observed for the density in section \ref{sec:t5peDensity}, the
-+
+property trends with temperature seen when using the Ewald summation
-+
+are reproduced with the {\sc sf} technique. Differences include the
-+
+calculated values of $\Delta H_\textrm{vap}$ underpredicting the Ewald
-+
+values. This is to be expected due to the direct weakening of the
-+
+electrostatic interaction through forced neutralization in {\sc
-+
+sf}. This results in an increase of the intermolecular potential
-+
+producing lower values from equation \ref{eq:DeltaHVap}. The slopes of
-+
+these values with temperature are similar to that seen using the Ewald
-+
+summation; however, they are both steeper than the experimental trend,
-+
+indirectly resulting in the inflated $C_p$ values at all temperatures.
-+
-+
+Above the supercooled regim\'{e}, $C_p$, $\kappa_T$, and $\alpha_p$
-+
+values all overlap within error. As indicated for the $\Delta
-+
+H_\textrm{vap}$ and $C_p$ results discussed in the previous paragraph,
-+
+the deviations between experiment and simulation in this region are
-+
+not the fault of the electrostatic summation methods but are due to
-+
+the TIP5P class model itself. Like most rigid, non-polarizable,
-+
+point-charge water models, the density decreases with temperature at a
-+
+much faster rate than experiment (see figure
-+
+\ref{fig:t5peDensities}). The reduced density leads to the inflated
-+
+compressibility and expansivity values at higher temperatures seen
-+
+here in figure \ref{fig:t5peThermo}. Incorporation of polarizability
-+
+and many-body effects are required in order for simulation to overcome
-+
+these differences with experiment.\cite{Laasonen93,Donchev06}
-+
-+
+At temperatures below the freezing point for experimental water, the
-+
+differences between {\sc sf} and the Ewald summation results are more
-+
+apparent. The larger $C_p$ and lower $\alpha_p$ values in this region
-+
+indicate a more pronounced transition in the supercooled regim\'{e},
-+
+particularly in the case of {\sc sf} without damping. This points to
-+
+the onset of a more frustrated or glassy behavior for TIP5P-E at
-+
+temperatures below 250K in these simulations. Because the systems are
-+
+locked in different regions of phase-space, comparisons between
-+
+properties at these temperatures are not exactly fair. This
-+
+observation is explored in more detail in section
-+
+\ref{sec:t5peDynamics}.
-+
-+
+The final thermodynamic property displayed in figure
-+
+\ref{fig:t5peThermo}, $\epsilon$, shows the greatest discrepancy
-+
+between the Ewald summation and the {\sc sf} technique (and experiment
-+
+for that matter). It is known that the dielectric constant is
-+
+dependent upon and quite sensitive to the imposed boundary
-+
+conditions.\cite{Neumann80,Neumann83} This is readily apparent in the
-+
+converged $\epsilon$ values accumulated for the {\sc sf}
-+
+simulations. Lack of a damping function results in dielectric
-+
+constants significantly smaller than that obtained using the Ewald
-+
+sum. Increasing the damping coefficient to 0.2\AA$^{-1}$ improves the
-+
+agreement considerably. It should be noted that the choice of the
-+
+``Ewald coefficient'' value also has a significant effect on the
-+
+calculated value when using the Ewald summation. In the simulations of
-+
+TIP5P-E with the Ewald sum, this screening parameter was tethered to
-+
+the simulation box size (as was the $R_\textrm{c}$).\cite{Rick04}
-+
+Systems with larger screening parameters reported larger dielectric
-+
+constant values, the same behavior we see here with {\sc sf}. In
-+
+section \ref{sec:dampingDielectric}, this connection is further
-+
+explored as optimal damping coefficients are determined for {\sc
-+
+sf} for capturing the dielectric behavior.
-+
-+
+\subsection{Dynamic Properties}\label{sec:t5peDynamics}
-+
-+
+To look at the dynamic properties of TIP5P-E when using the {\sc sf}
-+
+method, 200ps $NVE$ simulations were performed for each temperature at
-+
+the average density reported by the $NPT$ simulations. The
-+
+self-diffusion constants ($D$) were calculated with the Einstein
-+
+relation using the mean square displacement (MSD),
-+
+\begin{equation}
-+
+D = \frac{\langle\left|\mathbf{r}_i(t)-\mathbf{r}_i(0)\right|^2\rangle}{6t},
-+
+\label{eq:MSD}
-+
+\end{equation}
-+
+where $t$ is time, and $\mathbf{r}_i$ is the position of particle
-+
+$i$.\cite{Allen87} Figure \ref{fig:ExampleMSD} shows an example MSD
-+
+plot. As labeled in the figure, MSD plots consist of three distinct
-+
+regions:
-+
-+
+\begin{enumerate}[itemsep=0pt]
-+
+\item parabolic short-time ballistic motion,
-+
+\item linear diffusive regime, and
-+
+\item poor statistic region at long-time.
-+
+\end{enumerate}
-+
+The slope from the linear region (region 2) is used to calculate $D$.
-+
+\begin{figure}
-+
+\centering
-+
+\includegraphics[width=3.5in]{./figures/ExampleMSD.pdf}
-+
+\caption{Example plot of mean square displacement verses time. The
-+
+left red region is the ballistic motion regime, the middle green
-+
+region is the linear diffusive regime, and the right blue region is
-+
+the region with poor statistics.}
-+
+\label{fig:ExampleMSD}
-+
+\end{figure}
-+
-+
+\begin{figure}
-+
+\centering
-+
+\includegraphics[width=3.5in]{./figures/waterFrame.pdf}
-+
+\caption{Body-fixed coordinate frame for a water molecule. The
-+
+respective molecular principle axes point in the direction of the
-+
+labeled frame axes.}
-+
+\label{fig:waterFrame}
-+
+\end{figure}
-+
+In addition to translational diffusion, reorientational time constants
-+
+were calculated for comparisons with the Ewald simulations and with
-+
+experiments. These values were determined from 25ps $NVE$ trajectories
-+
+through calculation of the orientational time correlation function,
-+
+\begin{equation}
-+
+C_l^\alpha(t) = \left\langle P_l\left[\hat{\mathbf{u}}_i^\alpha(t)
-+
+                \cdot\hat{\mathbf{u}}_i^\alpha(0)\right]\right\rangle,
-+
+\label{eq:OrientCorr}
-+
+\end{equation}
-+
+where $P_l$ is the Legendre polynomial of order $l$ and
-+
+$\hat{\mathbf{u}}_i^\alpha$ is the unit vector of molecule $i$ along
-+
+principle axis $\alpha$. The principle axis frame for these water
-+
+molecules is shown in figure \ref{fig:waterFrame}. As an example,
-+
+$C_l^y$ is calculated from the time evolution of the unit vector
-+
+connecting the two hydrogen atoms.
-+
-+
+\begin{figure}
-+
+\centering
-+
+\includegraphics[width=3.5in]{./figures/ExampleOrientCorr.pdf}
-+
+\caption{Example plots of the orientational autocorrelation functions
-+
+for the first and second Legendre polynomials. These curves show the
-+
+time decay of the unit vector along the $y$ principle axis.}
-+
+\label{fig:OrientCorr}
-+
+\end{figure}
-+
+From the orientation autocorrelation functions, we can obtain time
-+
+constants for rotational relaxation. Figure \ref{fig:OrientCorr} shows
-+
+some example plots of orientational autocorrelation functions for the
-+
+first and second Legendre polynomials. The relatively short time
-+
+portions (between 1 and 3ps for water) of these curves can be fit to
-+
+an exponential decay to obtain these constants, and they are directly
-+
+comparable to water orientational relaxation times from nuclear
-+
+magnetic resonance (NMR). The relaxation constant obtained from
-+
+$C_2^y(t)$ is of particular interest because it is about the principle
-+
+axis with the minimum moment of inertia and should thereby dominate
-+
+the orientational relaxation of the molecule.\cite{Impey82} This means
-+
+that $C_2^y(t)$ should provide the best comparison to the NMR
-+
+relaxation data.
-+
-+
+\begin{figure}
-+
+\centering
-+
+\includegraphics[width=5.5in]{./figures/t5peDynamics.pdf}
-+
+\caption{Diffusion constants ({\it upper}) and reorientational time
-+
+constants ({\it lower}) for TIP5P-E using the Ewald sum and {\sc sf}
-+
+technique compared with experiment. Data at temperatures less that
-+
+$^\circ$C were not included in the $\tau_2^y$ plot to allow for
-+
+easier comparisons in the more relavent temperature regime.}
-+
+\label{fig:t5peDynamics}
-+
+\end{figure}
-+
+Results for the diffusion constants and reorientational time constants
-+
+are shown in figure \ref{fig:t5peDynamics}. From this figure, it is
-+
+apparent that the trends for both $D$ and $\tau_2^y$ of TIP5P-E using
-+
+the Ewald sum are reproduced with the {\sc sf} techinque. The enhanced
-+
+diffusion at high temperatures are again the product of the lower
-+
+densities in comparison with experiment and do not provide any special
-+
+insight into differences between the electrostatic summation
-+
+techniques. With the undamped {\sc sf} technique, TIP5P-E tends to
-+
+diffuse a little faster than with the Ewald sum; however, use of light
-+
+to moderate damping results in indistiguishable $D$ values. Though not
-+
+apparent in this figure, {\sc sf} values at the lowest temperature are
-+
+approximately an order of magnitude lower than with Ewald. These
-+
+values support the observation from section \ref{sec:t5peThermo} that
-+
+there appeared to be a change to a more glassy-like phase with the
-+
+{\sc sf} technique at these lower temperatures.
-+
-+
+The $\tau_2^y$ results in the lower frame of figure
-+
+\ref{fig:t5peDynamics} show a much greater difference between the {\sc
-+
+sf} results and the Ewald results. At all temperatures shown, TIP5P-E
-+
+relaxes faster than experiment with the Ewald sum while tracking
-+
+experiment fairly well when using the {\sc sf} technique, independent
-+
+of the choice of damping constant. Their are several possible reasons
-+
+for this deviation between techniques. The Ewald results were taken
-+
+shorter (10ps) trajectories than the {\sc sf} results (25ps). A quick
-+
+calculation from a 10ps trajectory with {\sc sf} with an $\alpha$ of
-+
+.2\AA$^-1$ at 25$^\circ$C showed a 0.4ps drop in $\tau_2^y$, placing
-+
+the result more in line with that obtained using the Ewald sum. These
-+
+results support this explanation; however, recomputing the results to
-+
+meet a poorer statistical standard is counter-productive. Assuming the
-+
+Ewald results are not the product of poor statistics, differences in
-+
+techniques to integrate the orientational motion could also play a
-+
+role. {\sc shake} is the most commonly used technique for
-+
+approximating rigid-body orientational motion,\cite{Ryckaert77} where
-+
+as in {\sc oopse}, we maintain and integrate the entire rotation
-+
+matrix using the {\sc dlm} method.\cite{Meineke05} Since {\sc shake}
-+
+is an iterative constraint technique, if the convergence tolerances
-+
+are raised for increased performance, error will accumulate in the
-+
+orientational motion. Finally, the Ewald results were calculated using
-+
+the $NVT$ ensemble, while the $NVE$ ensemble was used for {\sc sf}
-+
+calculations. The additional mode of motion due to the thermostat will
-+
+alter the dynamics, resulting in differences between $NVT$ and $NVE$
-+
+results. These differences are increasingly noticable as the
-+
+thermostat time constant decreases.
-+
-+
+\section{Damping of Point Multipoles}\label{sec:dampingMultipoles}
-+
-+
-+
-+
+\section{Damping and Dielectric Constants}\label{sec:dampingDielectric}
-+
-+
+\section{Conclusions}\label{sec:PairwiseConclusions}
-+
-+
+The above investigation of pairwise electrostatic summation techniques
-+
+shows that there are viable and computationally efficient alternatives
-+
+to the Ewald summation.  These methods are derived from the damped and
-+
+cutoff-neutralized Coulombic sum originally proposed by Wolf
-+
+\textit{et al.}\cite{Wolf99} In particular, the {\sc sf}
-+
+method, reformulated above as eqs. (\ref{eq:DSFPot}) and
-+
+(\ref{eq:DSFForces}), shows a remarkable ability to reproduce the
-+
+energetic and dynamic characteristics exhibited by simulations
-+
+employing lattice summation techniques.  The cumulative energy
-+
+difference results showed the undamped {\sc sf} and moderately damped
-+
+{\sc sp} methods produced results nearly identical to {\sc spme}.
-+
+Similarly for the dynamic features, the undamped or moderately damped
-+
+{\sc sf} and moderately damped {\sc sp} methods produce force and
-+
+torque vector magnitude and directions very similar to the expected
-+
+values.  These results translate into long-time dynamic behavior
-+
+equivalent to that produced in simulations using {\sc spme}.
-+
-+
+As in all purely-pairwise cutoff methods, these methods are expected
-+
+to scale approximately {\it linearly} with system size, and they are
-+
+easily parallelizable.  This should result in substantial reductions
-+
+in the computational cost of performing large simulations.
-+
-+
+Aside from the computational cost benefit, these techniques have
-+
+applicability in situations where the use of the Ewald sum can prove
-+
+problematic.  Of greatest interest is their potential use in
-+
+interfacial systems, where the unmodified lattice sum techniques
-+
+artificially accentuate the periodicity of the system in an
-+
+undesirable manner.  There have been alterations to the standard Ewald
-+
+techniques, via corrections and reformulations, to compensate for
-+
+these systems; but the pairwise techniques discussed here require no
-+
+modifications, making them natural tools to tackle these problems.
-+
+Additionally, this transferability gives them benefits over other
-+
+pairwise methods, like reaction field, because estimations of physical
-+
+properties (e.g. the dielectric constant) are unnecessary.
-+
-+
+If a researcher is using Monte Carlo simulations of large chemical
-+
+systems containing point charges, most structural features will be
-+
+accurately captured using the undamped {\sc sf} method or the {\sc sp}
-+
+method with an electrostatic damping of 0.2\AA$^{-1}$.  These methods
-+
+would also be appropriate for molecular dynamics simulations where the
-+
+data of interest is either structural or short-time dynamical
-+
+quantities.  For long-time dynamics and collective motions, the safest
-+
+pairwise method we have evaluated is the {\sc sf} method with an
-+
+electrostatic damping between 0.2 and 0.25\AA$^{-1}$.
-+
-+
+We are not suggesting that there is any flaw with the Ewald sum; in
-+
+fact, it is the standard by which these simple pairwise sums have been
-+
+judged.  However, these results do suggest that in the typical
-+
+simulations performed today, the Ewald summation may no longer be
-+
+required to obtain the level of accuracy most researchers have come to
-+
+expect.
-+
+\chapter{\label{chap:water}SIMPLE MODELS FOR WATER}
+\chapter{\label{chap:ice}PHASE BEHAVIOR OF WATER IN COMPUTER SIMULATIONS}
+\backmatter
+\bibliographystyle{ndthesis}
-<
+\bibliography{dissertation}
->
+\bibliography{dissertation}
+\end{document}

Diff Legend

-–
+Removed lines
-+
+Added lines
-<
+Changed lines
->
+Changed lines

Comparing trunk/fennellDissertation/dissertation.tex (file contents): Revision 2920 by chrisfen, Mon Jul 3 14:51:28 2006 UTC vs. Revision 2971 by chrisfen, Sun Aug 6 22:29:27 2006 UTC

Diff Legend

Comparing trunk/fennellDissertation/dissertation.tex (file contents):
Revision 2920 by chrisfen, Mon Jul 3 14:51:28 2006 UTC vs.
Revision 2971 by chrisfen, Sun Aug 6 22:29:27 2006 UTC