| 263 |
|
each of the three coordinates, finding the the intersection points of |
| 264 |
|
two ellipsoids is isomorphic to finding the roots a polynomial of |
| 265 |
|
degree 16. There are a number of polynomial root-finding methods in |
| 266 |
< |
the literature, [CITATIONS NEEDED] but numerically finding the roots |
| 266 |
> |
the literature,\cite{384119} [CITATIONS NEEDED] but numerically finding the roots |
| 267 |
|
of high-degree polynomials is generally an ill-conditioned |
| 268 |
< |
problem.[CITATION NEEDED] One simplification is to maintain velocity |
| 268 |
> |
problem.\cite{Hoffman:2001sf}[P157] One simplification is to maintain velocity |
| 269 |
|
scalings that are {\it as isotropic as possible}. To do this, we |
| 270 |
|
impose $x=y$, and to treat both the constraint and flux ellipsoids as |
| 271 |
|
2-dimensional ellipses. In reduced dimensionality, the |
| 326 |
|
|
| 327 |
|
We have implemented this methodology in our molecular dynamics code, |
| 328 |
|
OpenMD,\cite{Meineke:2005gd,openmd} performing the NIVS scaling moves |
| 329 |
< |
after each MD step. We have tested the method in a variety of |
| 330 |
< |
different systems, including homogeneous fluids (Lennard-Jones and |
| 331 |
< |
SPC/E water), crystalline solids ({\sc eam}~\cite{PhysRevB.33.7983} and |
| 332 |
< |
quantum Sutton-Chen ({\sc q-sc})~\cite{PhysRevB.59.3527} |
| 333 |
< |
models for Gold), and heterogeneous interfaces (QSC gold - SPC/E |
| 334 |
< |
water). The last of these systems would have been difficult to study |
| 335 |
< |
using previous RNEMD methods, but using velocity scaling moves, we can |
| 336 |
< |
even obtain estimates of the interfacial thermal conductivities ($G$). |
| 329 |
> |
after an MD step with a variable frequency. We have tested the method |
| 330 |
> |
in a variety of different systems, including homogeneous fluids |
| 331 |
> |
(Lennard-Jones and SPC/E water), crystalline solids ({\sc |
| 332 |
> |
eam}~\cite{PhysRevB.33.7983} and quantum Sutton-Chen ({\sc |
| 333 |
> |
q-sc})~\cite{PhysRevB.59.3527} models for Gold), and heterogeneous |
| 334 |
> |
interfaces (QSC gold - SPC/E water). The last of these systems would |
| 335 |
> |
have been difficult to study using previous RNEMD methods, but using |
| 336 |
> |
velocity scaling moves, we can even obtain estimates of the |
| 337 |
> |
interfacial thermal conductivities ($G$). |
| 338 |
|
|
| 339 |
|
\subsection{Simulation Cells} |
| 340 |
|
|
| 342 |
|
rectangular simulation cell using periodic boundary conditions in all |
| 343 |
|
three dimensions. The cells were longer along the $z$ axis and the |
| 344 |
|
space was divided into $N$ slabs along this axis (typically $N=20$). |
| 345 |
< |
The top slab ($n=1$) was designated the ``cold'' slab, while the |
| 346 |
< |
central slab ($n= N/2 + 1$) was designated the ``hot'' slab. In all |
| 345 |
> |
The top slab ($n=1$) was designated the ``hot'' slab, while the |
| 346 |
> |
central slab ($n= N/2 + 1$) was designated the ``cold'' slab. In all |
| 347 |
|
cases, simulations were first thermalized in canonical ensemble (NVT) |
| 348 |
|
using a Nos\'{e}-Hoover thermostat.\cite{Hoover85} then equilibrated in |
| 349 |
|
microcanonical ensemble (NVE) before introducing any non-equilibrium |
| 437 |
|
|
| 438 |
|
Our thermal conductivity calculations show that the NIVS method agrees |
| 439 |
|
well with the swapping method. Four different swap intervals were |
| 440 |
< |
tested (Table \ref{LJ}). With a fixed 10 fs [WHY NOT REDUCED |
| 441 |
< |
UNITS???] scaling interval, the target exchange kinetic energy |
| 442 |
< |
produced equivalent kinetic energy flux as in the swapping method. |
| 443 |
< |
Similar thermal gradients were observed with similar thermal flux |
| 444 |
< |
under the two different methods (Figure \ref{thermalGrad}). |
| 440 |
> |
tested (Table \ref{LJ}). With a fixed scaling interval of 10 time steps, |
| 441 |
> |
the target exchange kinetic energy produced equivalent kinetic energy |
| 442 |
> |
flux as in the swapping method. Similar thermal gradients were |
| 443 |
> |
observed with similar thermal flux under the two different methods |
| 444 |
> |
(Figure \ref{thermalGrad}). |
| 445 |
|
|
| 446 |
|
\begin{table*} |
| 447 |
|
\begin{minipage}{\linewidth} |
| 493 |
|
During these simulations, velocities were recorded every 1000 steps |
| 494 |
|
and was used to produce distributions of both velocity and speed in |
| 495 |
|
each of the slabs. From these distributions, we observed that under |
| 496 |
< |
relatively high non-physical kinetic energy flux, the spee of |
| 496 |
> |
relatively high non-physical kinetic energy flux, the speed of |
| 497 |
|
molecules in slabs where swapping occured could deviate from the |
| 498 |
|
Maxwell-Boltzmann distribution. This behavior was also noted by Tenney |
| 499 |
|
and Maginn.\cite{Maginn:2010} Figure \ref{thermalHist} shows how these |
| 596 |
|
ensemble.\cite{melchionna93} A fixed volume was chosen to match the |
| 597 |
|
average volume observed in the NPT simulations, and this was followed |
| 598 |
|
by equilibration, first in the canonical (NVT) ensemble, followed by a |
| 599 |
< |
[XXX ps] period under constant-NVE conditions without any momentum |
| 600 |
< |
flux. [YYY ps] was allowed to stabilize the system with an imposed |
| 601 |
< |
momentum transfer to create a gradient, and [ZZZ ps] was alotted for |
| 599 |
> |
100ps period under constant-NVE conditions without any momentum |
| 600 |
> |
flux. 100ps was allowed to stabilize the system with an imposed |
| 601 |
> |
momentum transfer to create a gradient, and 1ns was alotted for |
| 602 |
|
data collection under RNEMD. |
| 603 |
|
|
| 604 |
|
As shown in Figure \ref{spceGrad}, temperature gradients were |
| 684 |
|
potentials to test the behavior of scaling RNEMD. |
| 685 |
|
|
| 686 |
|
A face-centered-cubic (FCC) lattice was prepared containing 2880 Au |
| 687 |
< |
atoms. [LxMxN UNIT CELLS]. Simulations were run both with and |
| 688 |
< |
without isobaric-isothermal (NPT)~\cite{melchionna93} |
| 687 |
> |
atoms (i.e. ${6\times 6\times 20}$ unit cells). Simulations were run |
| 688 |
> |
both with and without isobaric-isothermal (NPT)~\cite{melchionna93} |
| 689 |
|
pre-equilibration at a target pressure of 1 atm. When equilibrated |
| 690 |
|
under NPT conditions, our simulation box expanded by approximately 1\% |
| 691 |
< |
[IN VOLUME OR LINEAR DIMENSIONS ?]. Following adjustment of the box |
| 692 |
< |
volume, equilibrations in both the canonical and microcanonical |
| 693 |
< |
ensembles were carried out. With the simulation cell divided evenly |
| 694 |
< |
into 10 slabs, different thermal gradients were established by |
| 695 |
< |
applying a set of target thermal transfer fluxes. |
| 691 |
> |
in volume. Following adjustment of the box volume, equilibrations in |
| 692 |
> |
both the canonical and microcanonical ensembles were carried out. With |
| 693 |
> |
the simulation cell divided evenly into 10 slabs, different thermal |
| 694 |
> |
gradients were established by applying a set of target thermal |
| 695 |
> |
transfer fluxes. |
| 696 |
|
|
| 697 |
|
The results for the thermal conductivity of gold are shown in Table |
| 698 |
|
\ref{AuThermal}. In these calculations, the end and middle slabs were |
| 702 |
|
200. This behavior has been observed previously by Richardson and |
| 703 |
|
Clancy, and has been attributed to the lack of electronic effects in |
| 704 |
|
these force fields.\cite{Clancy:1992} The non-equilibrium MD method |
| 705 |
< |
employed in their simulations produced comparable results to ours. It |
| 705 |
> |
employed in their simulations gave an thermal conductance estimation |
| 706 |
> |
of [FORCE FIELD] gold as [RESULT IN REF], which is comparable to ours. It |
| 707 |
|
should be noted that the density of the metal being simulated also |
| 708 |
|
greatly affects the thermal conductivity. With an expanded lattice, |
| 709 |
|
lower thermal conductance is expected (and observed). We also observed |
| 719 |
|
using two related force fields. Calculations were done at both |
| 720 |
|
experimental and equilibrated densities and at a range of |
| 721 |
|
temperatures and thermal flux rates. Uncertainties are |
| 722 |
< |
indicated in parentheses. [CLANCY COMPARISON? SWAPPING |
| 721 |
< |
COMPARISON?]} |
| 722 |
> |
indicated in parentheses. [SWAPPING COMPARISON?]} |
| 723 |
|
|
| 724 |
|
\begin{tabular}{|c|c|c|cc|} |
| 725 |
|
\hline |
| 759 |
|
relax under ambient temperature and pressure. A separate (but |
| 760 |
|
identically sized) box of SPC/E water was also equilibrated at ambient |
| 761 |
|
conditions. The two boxes were combined by removing all water |
| 762 |
< |
molecules withing 3 \AA radius of any gold atom. The final |
| 762 |
> |
molecules within 3 \AA radius of any gold atom. The final |
| 763 |
|
configuration contained 1862 SPC/E water molecules. |
| 764 |
|
|
| 765 |
|
After simulations of bulk water and crystal gold, a mixture system was |
