66 |
|
resulting in a pressure of approximately 1 atm at their respective |
67 |
|
temperatures. |
68 |
|
|
69 |
+ |
A single thermodynamic integration involves a sequence of simulations |
70 |
+ |
over which the system of interest is converted into a reference system |
71 |
+ |
for which the free energy is known. This transformation path is then |
72 |
+ |
integrated in order to determine the free energy difference between |
73 |
+ |
the two states: |
74 |
+ |
\begin{equation} |
75 |
+ |
\begin{center} |
76 |
+ |
\Delta A = \int_0^1\left\langle\frac{\partial V(\lambda |
77 |
+ |
)}{\partial\lambda}\right\rangle_\lambda d\lambda, |
78 |
+ |
\end{center} |
79 |
+ |
\end{equation} |
80 |
+ |
where $V$ is the interaction potential and $\lambda$ is the |
81 |
+ |
transformation parameter. Simulations are distributed unevenly along |
82 |
+ |
this path in order to sufficiently sample the regions of greatest |
83 |
+ |
change in the potential. Typical integrations in this study consisted |
84 |
+ |
of $\sim$25 simulations ranging from 300 ps (for the unaltered system) |
85 |
+ |
to 75 ps (near the reference state) in length. |
86 |
+ |
|
87 |
|
For the thermodynamic integration of molecular crystals, the Einstein |
88 |
|
Crystal is chosen as the reference state that the system is converted |
89 |
|
to over the course of the simulation. In an Einstein Crystal, the |
259 |
|
risk of spontaneous crystallization. However, this risk changes when |
260 |
|
applying a longer cutoff. |
261 |
|
|
262 |
+ |
\begin{figure} |
263 |
+ |
\includegraphics[width=\linewidth]{cutoffChange.eps} |
264 |
+ |
\caption{Free energy as a function of cutoff radius for (A) SSD/E, (B) |
265 |
+ |
TIP3P, and (C) SSD/RF. Data points omitted include SSD/E: $I_c$ 12 |
266 |
+ |
\AA\, TIP3P: $I_c$ 12 \AA\ and B 12 \AA\, and SSD/RF: $I_c$ 9 |
267 |
+ |
\AA\. These crystals are unstable at 200 K and rapidly convert into a |
268 |
+ |
liquid. The connecting lines are qualitative visual aid.} |
269 |
+ |
\label{incCutoff} |
270 |
+ |
\end{figure} |
271 |
+ |
|
272 |
|
Increasing the cutoff radius in simulations of the more |
273 |
|
computationally efficient water models was done in order to evaluate |
274 |
|
the trend in free energy values when moving to systems that do not |
306 |
|
\begin{minipage}{\linewidth} |
307 |
|
\renewcommand{\thefootnote}{\thempfootnote} |
308 |
|
\begin{center} |
309 |
< |
\caption{The free energy of the studied ice polymorphs after applying the energy difference attributed to the inclusion of the PME long-range interaction correction. Units are kcal/mol.} |
309 |
> |
\caption{The free energy of the studied ice polymorphs after applying |
310 |
> |
the energy difference attributed to the inclusion of the PME |
311 |
> |
long-range interaction correction. Units are kcal/mol.} |
312 |
|
\begin{tabular}{ l c c c c } |
313 |
|
\hline \\[-7mm] |
314 |
|
\ \ Water Model \ \ & \ \ \ \ \ $I_h$ \ \ & \ \ \ \ \ $I_c$ \ \ & \ \quad \ \ \ \ B \ \ & \ \ \ \ \ Ice-{\it i} \ \ \\ |
323 |
|
|
324 |
|
\section{Conclusions} |
325 |
|
|
326 |
+ |
The free energy for proton ordered variants of hexagonal and cubic ice |
327 |
+ |
$I$, ice B, and recently discovered Ice-{\it i} where calculated under |
328 |
+ |
standard conditions for several common water models via thermodynamic |
329 |
+ |
integration. All the water models studied show Ice-{\it i} to be the |
330 |
+ |
minimum free energy crystal structure in the with a 9 \AA\ switching |
331 |
+ |
function cutoff. Calculated melting and boiling points show |
332 |
+ |
surprisingly good agreement with the experimental values; however, the |
333 |
+ |
solid phase at 1 atm is Ice-{\it i}, not ice $I_h$. The effect of |
334 |
+ |
interaction truncation was investigated through variation of the |
335 |
+ |
cutoff radius, use of a reaction field parameterized model, and |
336 |
+ |
estimation of the results in the presence of the Ewald summation |
337 |
+ |
correction. Interaction truncation has a significant effect on the |
338 |
+ |
computed free energy values, but Ice-{\it i} is still observed to be a |
339 |
+ |
relavent ice polymorph in simulation studies. |
340 |
+ |
|
341 |
|
\section{Acknowledgments} |
342 |
|
Support for this project was provided by the National Science |
343 |
|
Foundation under grant CHE-0134881. Computation time was provided by |
344 |
< |
the Notre Dame Bunch-of-Boxes (B.o.B) computer cluster under NSF grant |
345 |
< |
DMR-0079647. |
344 |
> |
the Notre Dame High Performance Computing Cluster and the Notre Dame |
345 |
> |
Bunch-of-Boxes (B.o.B) computer cluster (NSF grant DMR-0079647). |
346 |
|
|
347 |
|
\newpage |
348 |
|
|