| 100 |
|
\label{thermalDemo} |
| 101 |
|
\end{figure} |
| 102 |
|
|
| 103 |
< |
RNEMD is preferable in many ways to the forward NEMD methods |
| 104 |
< |
[CITATIONS NEEDED] because it imposes what is typically difficult to measure |
| 105 |
< |
(a flux or stress) and it is typically much easier to compute momentum |
| 106 |
< |
gradients or strains (the response). For similar reasons, RNEMD is |
| 107 |
< |
also preferable to slowly-converging equilibrium methods for measuring |
| 108 |
< |
thermal conductivity and shear viscosity (using Green-Kubo relations |
| 109 |
< |
[CITATIONS NEEDED] or the Helfand moment approach of Viscardy {\it et |
| 103 |
> |
RNEMD is preferable in many ways to the forward NEMD |
| 104 |
> |
methods\cite{ISI:A1988Q205300014,hess:209,Vasquez:2004fk,backer:154503,ISI:000266247600008} |
| 105 |
> |
because it imposes what is typically difficult to measure (a flux or |
| 106 |
> |
stress) and it is typically much easier to compute momentum gradients |
| 107 |
> |
or strains (the response). For similar reasons, RNEMD is also |
| 108 |
> |
preferable to slowly-converging equilibrium methods for measuring |
| 109 |
> |
thermal conductivity and shear viscosity (using Green-Kubo |
| 110 |
> |
relations\cite{daivis:541,mondello:9327} or the Helfand moment |
| 111 |
> |
approach of Viscardy {\it et |
| 112 |
|
al.}\cite{Viscardy:2007bh,Viscardy:2007lq}) because these rely on |
| 113 |
|
computing difficult to measure quantities. |
| 114 |
|
|
| 118 |
|
typically samples from the same manifold of states in the |
| 119 |
|
microcanonical ensemble. |
| 120 |
|
|
| 121 |
< |
Recently, Tenney and Maginn\cite{ISI:000273472300004} have discovered |
| 121 |
> |
Recently, Tenney and Maginn\cite{Maginn:2010} have discovered |
| 122 |
|
some problems with the original RNEMD swap technique. Notably, large |
| 123 |
|
momentum fluxes (equivalent to frequent momentum swaps between the |
| 124 |
|
slabs) can result in ``notched'', ``peaked'' and generally non-thermal |
| 246 |
|
|
| 247 |
|
\begin{figure} |
| 248 |
|
\includegraphics[width=\linewidth]{ellipsoids} |
| 249 |
< |
\caption{Scaling points which maintain both constant energy and |
| 250 |
< |
constant linear momentum of the system lie on the surface of the |
| 251 |
< |
{\it constraint ellipsoid} while points which generate the target |
| 252 |
< |
momentum flux lie on the surface of the {\it flux ellipsoid}. The |
| 253 |
< |
velocity distributions in the cold bin are scaled by only those |
| 254 |
< |
points which lie on both ellipsoids.} |
| 249 |
> |
\caption{Illustration from the perspective of a space having cold |
| 250 |
> |
slab scaling coefficients as its coordinates. Scaling points which |
| 251 |
> |
maintain both constant energy and constant linear momentum of the |
| 252 |
> |
system lie on the surface of the {\it constraint ellipsoid} while |
| 253 |
> |
points which generate the target momentum flux lie on the surface of |
| 254 |
> |
the {\it flux ellipsoid}. The velocity distributions in the cold bin |
| 255 |
> |
are scaled by only those points which lie on both ellipsoids.} |
| 256 |
|
\label{ellipsoids} |
| 257 |
|
\end{figure} |
| 258 |
|
|
| 336 |
|
periodic boundary condition, and devided into ${N = 20}$ slabs. In each swap, |
| 337 |
|
the top slab ${(n = 1)}$ exchange the most negative $x$ momentum with the |
| 338 |
|
most positive $x$ momentum in the center slab ${(n = N/2 + 1)}$. Referred |
| 339 |
< |
to Tenney {\it et al.}\cite{ISI:000273472300004}, a series of swapping |
| 339 |
> |
to Tenney {\it et al.}\cite{Maginn:2010}, a series of swapping |
| 340 |
|
frequency were chosen. According to each result from swapping |
| 341 |
|
RNEMD, scaling RNEMD simulations were run with the target momentum |
| 342 |
|
flux set to produce a similar momentum flux, and consequently shear |
| 482 |
|
|
| 483 |
|
\begin{figure} |
| 484 |
|
\includegraphics[width=\linewidth]{thermalGrad} |
| 485 |
< |
\caption{Temperature gradients under various kinetic energy flux of |
| 486 |
< |
thermal conductivity simulations} |
| 485 |
> |
\caption{NIVS-RNEMD method introduced similar temperature gradients |
| 486 |
> |
compared to ``swapping'' method under various kinetic energy flux in |
| 487 |
> |
thermal conductivity simulations.} |
| 488 |
|
\label{thermalGrad} |
| 489 |
|
\end{figure} |
| 490 |
|
|
| 495 |
|
that under relatively high unphysical kinetic energy flux, speed and |
| 496 |
|
velocity distribution of molecules in slabs where swapping occured |
| 497 |
|
could deviate from Maxwell-Boltzmann distribution. Figure |
| 498 |
< |
\ref{histSwap} illustrates how these distributions deviate from an |
| 498 |
> |
\ref{thermalHist} a) illustrates how these distributions deviate from an |
| 499 |
|
ideal distribution. In high temperature slab, probability density in |
| 500 |
|
low speed is confidently smaller than ideal curve fit; in low |
| 501 |
|
temperature slab, probability density in high speed is smaller than |
| 515 |
|
a minimum level to produce an observable thermal gradient under noise |
| 516 |
|
interference. Consequently, swapping RNEMD has a relatively narrow |
| 517 |
|
choice of swapping rate to satisfy these above restrictions. |
| 514 |
– |
|
| 515 |
– |
\begin{figure} |
| 516 |
– |
\includegraphics[width=\linewidth]{histSwap} |
| 517 |
– |
\caption{Speed distribution for thermal conductivity using swapping |
| 518 |
– |
RNEMD. Shown is from the simulation with 250 fs exchange interval.} |
| 519 |
– |
\label{histSwap} |
| 520 |
– |
\end{figure} |
| 518 |
|
|
| 519 |
|
Comparatively, NIVS-RNEMD has a speed distribution closer to the ideal |
| 520 |
< |
curve fit (Figure \ref{histScale}). Essentially, after scaling, a |
| 520 |
> |
curve fit (Figure \ref{thermalHist} b). Essentially, after scaling, a |
| 521 |
|
Gaussian distribution function would remain Gaussian. Although a |
| 522 |
|
single scaling is non-isotropic in all three dimensions, our scaling |
| 523 |
|
coefficient criteria could help maintian the scaling region as |
| 530 |
|
to the distribution of velocity and speed in the exchange regions. |
| 531 |
|
|
| 532 |
|
\begin{figure} |
| 533 |
< |
\includegraphics[width=\linewidth]{histScale} |
| 534 |
< |
\caption{Speed distribution for thermal conductivity using scaling |
| 535 |
< |
RNEMD. Shown is from the simulation with an equilvalent thermal flux |
| 536 |
< |
as an 250 fs exchange interval swapping simulation.} |
| 537 |
< |
\label{histScale} |
| 533 |
> |
\includegraphics[width=\linewidth]{thermalHist} |
| 534 |
> |
\caption{Speed distribution for thermal conductivity using a) |
| 535 |
> |
``swapping'' and b) NIVS- RNEMD methods. Shown is from the |
| 536 |
> |
simulations with an exchange or equilvalent exchange interval of 250 |
| 537 |
> |
fs. In circled areas, distributions from ``swapping'' RNEMD |
| 538 |
> |
simulation have deviation from ideal Maxwell-Boltzmann distribution |
| 539 |
> |
(curves fit for each distribution).} |
| 540 |
> |
\label{thermalHist} |
| 541 |
|
\end{figure} |
| 542 |
|
|
| 543 |
|
\subsubsection{SPC/E Water} |
| 544 |
|
Our results of SPC/E water thermal conductivity are comparable to |
| 545 |
< |
Bedrov {\it et al.}\cite{ISI:000090151400044}, which employed the |
| 545 |
> |
Bedrov {\it et al.}\cite{Bedrov:2000}, which employed the |
| 546 |
|
previous swapping RNEMD method for their calculation. Bedrov {\it et |
| 547 |
< |
al.}\cite{ISI:000090151400044} argued that exchange of the molecule |
| 547 |
> |
al.}\cite{Bedrov:2000} argued that exchange of the molecule |
| 548 |
|
center-of-mass velocities instead of single atom velocities in a |
| 549 |
|
molecule conserves the total kinetic energy and linear momentum. This |
| 550 |
|
principle is adopted in our simulations. Scaling is applied to the |
| 560 |
|
|
| 561 |
|
\begin{figure} |
| 562 |
|
\includegraphics[width=\linewidth]{spceGrad} |
| 563 |
< |
\caption{Temperature gradients for SPC/E water thermal conductivity.} |
| 563 |
> |
\caption{Temperature gradients in SPC/E water thermal conductivity |
| 564 |
> |
simulations.} |
| 565 |
|
\label{spceGrad} |
| 566 |
|
\end{figure} |
| 567 |
|
|
| 576 |
|
\begin{tabular}{cccc} |
| 577 |
|
\hline |
| 578 |
|
$\langle dT/dz\rangle$(K/\AA) & & $\lambda$(W/m/K) & \\ |
| 579 |
< |
& This work & Previous simulations\cite{ISI:000090151400044} & |
| 579 |
> |
& This work & Previous simulations\cite{Bedrov:2000} & |
| 580 |
|
Experiment$^a$\\ |
| 581 |
|
\hline |
| 582 |
|
0.38 & 0.816(0.044) & & 0.64\\ |
| 671 |
|
NIVS-RNEMD method were proved valid, calculation of gold/water |
| 672 |
|
interfacial thermal conductivity was followed. It is found out that |
| 673 |
|
the low interfacial conductance is probably due to the hydrophobic |
| 674 |
< |
surface in our system. Figure \ref{interfaceDensity} demonstrates mass |
| 674 |
> |
surface in our system. Figure \ref{interface} (a) demonstrates mass |
| 675 |
|
density change along $z$-axis, which is perpendicular to the |
| 676 |
|
gold/water interface. It is observed that water density significantly |
| 677 |
|
decreases when approaching the surface. Under this low thermal |
| 678 |
|
conductance, both gold and water phase have sufficient time to |
| 679 |
|
eliminate temperature difference inside respectively (Figure |
| 680 |
< |
\ref{interfaceGrad}). With indistinguishable temperature difference |
| 680 |
> |
\ref{interface} b). With indistinguishable temperature difference |
| 681 |
|
within respective phase, it is valid to assume that the temperature |
| 682 |
|
difference between gold and water on surface would be approximately |
| 683 |
|
the same as the difference between the gold and water phase. This |
| 690 |
|
errors than our calculation results on homogeneous systems. |
| 691 |
|
|
| 692 |
|
\begin{figure} |
| 693 |
< |
\includegraphics[width=\linewidth]{interfaceDensity} |
| 694 |
< |
\caption{Density profile for interfacial thermal conductivity |
| 695 |
< |
simulation box. Significant water density decrease is observed on |
| 696 |
< |
gold surface.} |
| 697 |
< |
\label{interfaceDensity} |
| 698 |
< |
\end{figure} |
| 699 |
< |
|
| 699 |
< |
\begin{figure} |
| 700 |
< |
\includegraphics[width=\linewidth]{interfaceGrad} |
| 701 |
< |
\caption{Temperature profiles for interfacial thermal conductivity |
| 702 |
< |
simulation box. Temperatures of different slabs in the same phase |
| 703 |
< |
show no significant difference.} |
| 704 |
< |
\label{interfaceGrad} |
| 693 |
> |
\includegraphics[width=\linewidth]{interface} |
| 694 |
> |
\caption{Simulation results for Gold/Water interfacial thermal |
| 695 |
> |
conductivity: (a) Significant water density decrease is observed on |
| 696 |
> |
crystalline gold surface. (b) Temperature profiles for a series of |
| 697 |
> |
simulations. Temperatures of different slabs in the same phase show |
| 698 |
> |
no significant differences.} |
| 699 |
> |
\label{interface} |
| 700 |
|
\end{figure} |
| 701 |
|
|
| 702 |
|
\begin{table*} |
| 732 |
|
momentum flux would theoretically result in swap method. All the scale |
| 733 |
|
method results were from simulations that had a scaling interval of 10 |
| 734 |
|
time steps. The average molecular momentum gradients of these samples |
| 735 |
< |
are shown in Figure \ref{shearGrad}. |
| 735 |
> |
are shown in Figure \ref{shear} (a) and (b). |
| 736 |
|
|
| 737 |
|
\begin{table*} |
| 738 |
|
\begin{minipage}{\linewidth} |
| 759 |
|
\end{table*} |
| 760 |
|
|
| 761 |
|
\begin{figure} |
| 762 |
< |
\includegraphics[width=\linewidth]{shearGrad} |
| 763 |
< |
\caption{Average momentum gradients of shear viscosity simulations} |
| 764 |
< |
\label{shearGrad} |
| 762 |
> |
\includegraphics[width=\linewidth]{shear} |
| 763 |
> |
\caption{Average momentum gradients in shear viscosity simulations, |
| 764 |
> |
using (a) ``swapping'' method and (b) NIVS-RNEMD method |
| 765 |
> |
respectively. (c) Temperature difference among x and y, z dimensions |
| 766 |
> |
observed when using NIVS-RNEMD with equivalent exchange interval of |
| 767 |
> |
500 fs.} |
| 768 |
> |
\label{shear} |
| 769 |
|
\end{figure} |
| 770 |
|
|
| 772 |
– |
\begin{figure} |
| 773 |
– |
\includegraphics[width=\linewidth]{shearTempScale} |
| 774 |
– |
\caption{Temperature profile for scaling RNEMD simulation.} |
| 775 |
– |
\label{shearTempScale} |
| 776 |
– |
\end{figure} |
| 771 |
|
However, observations of temperatures along three dimensions show that |
| 772 |
|
inhomogeneity occurs in scaling RNEMD simulations, particularly in the |
| 773 |
< |
two slabs which were scaled. Figure \ref{shearTempScale} indicate that with |
| 773 |
> |
two slabs which were scaled. Figure \ref{shear} (c) indicate that with |
| 774 |
|
relatively large imposed momentum flux, the temperature difference among $x$ |
| 775 |
|
and the other two dimensions was significant. This would result from the |
| 776 |
|
algorithm of scaling method. From Eq. \ref{eq:fluxPlane}, after |
