| 53 |
|
scaling (NIVS) on the molecules in specific regions of a system, it |
| 54 |
|
is possible to impose momentum or thermal flux between regions of a |
| 55 |
|
simulation while conserving the linear momentum and total energy of |
| 56 |
< |
the system. To test the methods, we have computed the thermal |
| 56 |
> |
the system. To test the method, we have computed the thermal |
| 57 |
|
conductivity of model liquid and solid systems as well as the |
| 58 |
|
interfacial thermal conductivity of a metal-water interface. We |
| 59 |
|
find that the NIVS-RNEMD improves the problematic velocity |
| 108 |
|
methods\cite{ISI:A1988Q205300014,hess:209,Vasquez:2004fk,backer:154503,ISI:000266247600008} |
| 109 |
|
because it imposes what is typically difficult to measure (a flux or |
| 110 |
|
stress) and it is typically much easier to compute the response |
| 111 |
< |
(momentum gradients or strains). For similar reasons, RNEMD is also |
| 111 |
> |
(momentum gradients or strains). For similar reasons, RNEMD is also |
| 112 |
|
preferable to slowly-converging equilibrium methods for measuring |
| 113 |
|
thermal conductivity and shear viscosity (using Green-Kubo |
| 114 |
|
relations\cite{daivis:541,mondello:9327} or the Helfand moment |
| 115 |
|
approach of Viscardy {\it et |
| 116 |
|
al.}\cite{Viscardy:2007bh,Viscardy:2007lq}) because these rely on |
| 117 |
< |
computing difficult to measure quantities. |
| 117 |
> |
computing and integrating long-time correlation functions that are |
| 118 |
> |
subject to noise issues. |
| 119 |
|
|
| 120 |
|
Another attractive feature of RNEMD is that it conserves both total |
| 121 |
|
linear momentum and total energy during the swaps (as long as the two |
| 122 |
|
molecules have the same identity), so the swapped configurations are |
| 123 |
|
typically samples from the same manifold of states in the |
| 124 |
< |
microcanonical ensemble. |
| 124 |
> |
microcanonical ensemble. Furthermore, the method is applicable with |
| 125 |
> |
different ensembles, unlike the heat-exchange algorithm proposed by |
| 126 |
> |
Hafskjold {\it et al.} \cite{HeX:1994,HeX:1993}, which is incompatible |
| 127 |
> |
with non-microcanonical ensemble. |
| 128 |
|
|
| 129 |
|
Recently, Tenney and Maginn\cite{Maginn:2010} have discovered some |
| 130 |
|
problems with the original RNEMD swap technique. Notably, large |
| 293 |
|
isotropic fluid as possible. With this in mind, we would like the |
| 294 |
|
kinetic energies in the three different directions could become as |
| 295 |
|
close as each other as possible after each scaling. Simultaneously, |
| 296 |
< |
one would also like each scaling as gentle as possible, i.e. ${x,y,z |
| 297 |
< |
\rightarrow 1}$, in order to avoid large perturbation to the system. |
| 298 |
< |
To do this, we pick the intersection point which maintains the three |
| 299 |
< |
scaling variables ${x, y, z}$ as well as the ratio of kinetic energies |
| 300 |
< |
${K_c^z/K_c^x, K_c^z/K_c^y}$ as close as possible to 1. |
| 296 |
> |
one would also like each scaling to be as gentle as possible, i.e. |
| 297 |
> |
${x,y,z \rightarrow 1}$, in order to avoid large perturbations to the |
| 298 |
> |
system. To do this, we pick the intersection point which maintains |
| 299 |
> |
the three scaling variables ${x, y, z}$ as well as the ratio of |
| 300 |
> |
kinetic energies ${K_c^z/K_c^x, K_c^z/K_c^y}$ as close as possible to |
| 301 |
> |
1. |
| 302 |
|
|
| 303 |
|
After the valid scaling parameters are arrived at by solving geometric |
| 304 |
|
intersection problems in $x, y, z$ space in order to obtain cold slab |
| 341 |
|
solids, using both the embedded atom method |
| 342 |
|
(EAM)~\cite{PhysRevB.33.7983} and quantum Sutton-Chen |
| 343 |
|
(QSC)~\cite{PhysRevB.59.3527} models for Gold, and heterogeneous |
| 344 |
< |
interfaces (QSC gold - SPC/E water). The last of these systems would |
| 345 |
< |
have been difficult to study using previous RNEMD methods, but the |
| 344 |
> |
interfaces (QSC gold - SPC/E water). Even though previous RNEMD |
| 345 |
> |
methods might remain usable for the last of these systems, energy |
| 346 |
> |
transfer from imaginary elastic collisions would be less effective |
| 347 |
> |
when the two particles involved have larger mass difference, and thus |
| 348 |
> |
affect the actuall implementation of these methods. However, our |
| 349 |
|
current method can easily provide estimates of the interfacial thermal |
| 350 |
|
conductivity ($G$). |
| 351 |
|
|
| 391 |
|
NIVS-RNEMD simulation was carried out using a target momentum flux set |
| 392 |
|
to produce the same flux experienced in the swapping simulation. |
| 393 |
|
|
| 394 |
< |
To test the temperature homogeneity, directional momentum and |
| 395 |
< |
temperature distributions were accumulated for molecules in each of |
| 396 |
< |
the slabs. Transport coefficients were computed using the temperature |
| 397 |
< |
(and momentum) gradients across the $z$-axis as well as the total |
| 398 |
< |
momentum flux the system experienced during data collection portion of |
| 399 |
< |
the simulation. |
| 394 |
> |
To test the temperature homogeneity, momentum and temperature |
| 395 |
> |
distributions (for all three dimensions) were accumulated for |
| 396 |
> |
molecules in each of the slabs. Transport coefficients were computed |
| 397 |
> |
using the temperature (and momentum) gradients across the $z$-axis as |
| 398 |
> |
well as the total momentum flux the system experienced during data |
| 399 |
> |
collection portion of the simulation. |
| 400 |
|
|
| 401 |
|
\subsection{Shear viscosities} |
| 402 |
|
|
| 437 |
|
be approximated as: |
| 438 |
|
|
| 439 |
|
\begin{equation} |
| 440 |
< |
G = \frac{E_{total}}{2 t L_x L_y \left( \langle T_{gold}\rangle - |
| 441 |
< |
\langle T_{water}\rangle \right)} |
| 440 |
> |
G = \frac{E_{total}}{2 t L_x L_y \left( \langle T_\mathrm{hot}\rangle - |
| 441 |
> |
\langle T_\mathrm{cold}\rangle \right)} |
| 442 |
|
\label{interfaceCalc} |
| 443 |
|
\end{equation} |
| 444 |
|
where ${E_{total}}$ is the imposed non-physical kinetic energy |
| 445 |
< |
transfer and ${\langle T_{gold}\rangle}$ and ${\langle |
| 446 |
< |
T_{water}\rangle}$ are the average observed temperature of gold and |
| 447 |
< |
water phases respectively. If the interfacial conductance is {\it |
| 448 |
< |
not} small, it is also be possible to compute the interfacial |
| 449 |
< |
thermal conductivity using this method utilizing the change in the |
| 450 |
< |
slope of the thermal gradient ($\partial^2 \langle T \rangle / \partial |
| 451 |
< |
z^2$) at the interface. |
| 445 |
> |
transfer and ${\langle T_\mathrm{hot}\rangle}$ and ${\langle |
| 446 |
> |
T_\mathrm{cold}\rangle}$ are the average observed temperature of the |
| 447 |
> |
two separated phases. If the interfacial conductance is {\it not} |
| 448 |
> |
small, it would also be possible to compute the interfacial thermal |
| 449 |
> |
conductivity using this method by computing the change in the slope of |
| 450 |
> |
the thermal gradient ($\partial^2 \langle T \rangle / |
| 451 |
> |
\partial z^2$) at the interface. |
| 452 |
|
|
| 453 |
|
\section{Results} |
| 454 |
|
|
| 472 |
|
Our thermal conductivity calculations show that the NIVS method agrees |
| 473 |
|
well with the swapping method. Five different swap intervals were |
| 474 |
|
tested (Table \ref{LJ}). Similar thermal gradients were observed with |
| 475 |
< |
similar thermal flux under the two different methods (Figure |
| 475 |
> |
similar thermal flux under the two different methods (Fig. |
| 476 |
|
\ref{thermalGrad}). Furthermore, the 1-d temperature profiles showed |
| 477 |
< |
no observable differences between the $x$, $y$ and $z$ axes (Figure |
| 478 |
< |
\ref{thermalGrad} c), so even though we are using a non-isotropic |
| 479 |
< |
scaling method, none of the three directions are experience |
| 480 |
< |
disproportionate heating due to the velocity scaling. |
| 477 |
> |
no observable differences between the $x$, $y$ and $z$ axes (lower |
| 478 |
> |
panel of Fig. \ref{thermalGrad}), so even though we are using a |
| 479 |
> |
non-isotropic scaling method, none of the three directions are |
| 480 |
> |
experience disproportionate heating due to the velocity scaling. |
| 481 |
|
|
| 482 |
|
\begin{table*} |
| 483 |
|
\begin{minipage}{\linewidth} |
| 496 |
|
\multicolumn{2}{|c}{Swapping RNEMD} & |
| 497 |
|
\multicolumn{2}{|c|}{NIVS-RNEMD} \\ |
| 498 |
|
\hline |
| 499 |
< |
\multirow{2}{*}{Swap Interval (fs)} & Equivalent $J_z^*$ or & |
| 500 |
< |
\multirow{2}{*}{$\lambda^*_{swap}$} & |
| 499 |
> |
\multirow{2}{*}{Swap Interval} & Equivalent $J_z^*$ or & |
| 500 |
> |
\multirow{2}{*}{$\lambda^*_{swap}$} & |
| 501 |
|
\multirow{2}{*}{$\eta^*_{swap}$} & |
| 502 |
|
\multirow{2}{*}{$\lambda^*_{scale}$} & |
| 503 |
|
\multirow{2}{*}{$\eta^*_{scale}$} \\ |
| 504 |
< |
& $j_z^*(p_x)$ (reduced units) & & & & \\ |
| 504 |
> |
(timesteps) & $j_z^*(p_x)$ (reduced units) & & & & \\ |
| 505 |
|
\hline |
| 506 |
|
250 & 0.16 & 7.03(0.34) & 3.57(0.06) & 7.30(0.10) & 3.54(0.04)\\ |
| 507 |
|
500 & 0.09 & 7.03(0.14) & 3.64(0.05) & 6.95(0.09) & 3.76(0.09)\\ |
| 527 |
|
|
| 528 |
|
\subsubsection*{Velocity Distributions} |
| 529 |
|
|
| 530 |
< |
During these simulations, velocities were recorded every 1000 steps |
| 531 |
< |
and were used to produce distributions of both velocity and speed in |
| 532 |
< |
each of the slabs. From these distributions, we observed that under |
| 533 |
< |
relatively high non-physical kinetic energy flux, the speed of |
| 534 |
< |
molecules in slabs where swapping occured could deviate from the |
| 530 |
> |
To test the effects on the velocity distributions, we accumulated |
| 531 |
> |
velocities every 100 steps and produced distributions of both velocity |
| 532 |
> |
and speed in each of the slabs. From these distributions, we observed |
| 533 |
> |
that under high non-physical kinetic energy flux, the speed of |
| 534 |
> |
molecules in slabs where {\it swapping} occured could deviate from the |
| 535 |
|
Maxwell-Boltzmann distribution. This behavior was also noted by Tenney |
| 536 |
|
and Maginn.\cite{Maginn:2010} Figure \ref{thermalHist} shows how these |
| 537 |
|
distributions deviate from an ideal distribution. In the ``hot'' slab, |
| 538 |
|
the probability density is notched at low speeds and has a substantial |
| 539 |
< |
shoulder at higher speeds relative to the ideal MB distribution. In |
| 540 |
< |
the cold slab, the opposite notching and shouldering occurs. This |
| 541 |
< |
phenomenon is more obvious at higher swapping rates. |
| 539 |
> |
shoulder at higher speeds relative to the ideal distribution. In the |
| 540 |
> |
cold slab, the opposite notching and shouldering occurs. This |
| 541 |
> |
phenomenon is more obvious at high swapping rates. |
| 542 |
|
|
| 543 |
|
The peak of the velocity distribution is substantially flattened in |
| 544 |
|
the hot slab, and is overly sharp (with truncated wings) in the cold |
| 545 |
|
slab. This problem is rooted in the mechanism of the swapping method. |
| 546 |
|
Continually depleting low (high) speed particles in the high (low) |
| 547 |
< |
temperature slab is not complemented by diffusions of low (high) speed |
| 548 |
< |
particles from neighboring slabs, unless the swapping rate is |
| 547 |
> |
temperature slab is not complemented by diffusion of low (high) speed |
| 548 |
> |
particles from neighboring slabs unless the swapping rate is |
| 549 |
|
sufficiently small. Simutaneously, surplus low speed particles in the |
| 550 |
< |
low temperature slab do not have sufficient time to diffuse to |
| 551 |
< |
neighboring slabs. Since the thermal exchange rate must reach a |
| 552 |
< |
minimum level to produce an observable thermal gradient, the |
| 553 |
< |
swapping-method RNEMD has a relatively narrow choice of exchange times |
| 546 |
< |
that can be utilized. |
| 550 |
> |
cold slab do not have sufficient time to diffuse to neighboring slabs. |
| 551 |
> |
Since the thermal exchange rate must reach a minimum level to produce |
| 552 |
> |
an observable thermal gradient, the swapping-method RNEMD has a |
| 553 |
> |
relatively narrow choice of exchange times that can be utilized. |
| 554 |
|
|
| 555 |
|
For comparison, NIVS-RNEMD produces a speed distribution closer to the |
| 556 |
< |
Maxwell-Boltzmann curve (Figure \ref{thermalHist}). The reason for |
| 557 |
< |
this is simple; upon velocity scaling, a Gaussian distribution remains |
| 556 |
> |
Maxwell-Boltzmann curve (Fig. \ref{thermalHist}). The reason for this |
| 557 |
> |
is simple; upon velocity scaling, a Gaussian distribution remains |
| 558 |
|
Gaussian. Although a single scaling move is non-isotropic in three |
| 559 |
|
dimensions, our criteria for choosing a set of scaling coefficients |
| 560 |
|
helps maintain the distributions as close to isotropic as possible. |
| 561 |
< |
Consequently, NIVS-RNEMD is able to impose an unphysical thermal flux |
| 562 |
< |
as the previous RNEMD methods but without large perturbations to the |
| 563 |
< |
velocity distributions in the two slabs. |
| 561 |
> |
Consequently, NIVS-RNEMD is able to impose a non-physical thermal flux |
| 562 |
> |
without large perturbations to the velocity distributions in the two |
| 563 |
> |
slabs. |
| 564 |
|
|
| 565 |
|
\begin{figure} |
| 566 |
|
\includegraphics[width=\linewidth]{thermalHist} |
| 579 |
|
|
| 580 |
|
\subsubsection*{Shear Viscosity} |
| 581 |
|
Our calculations (Table \ref{LJ}) show that velocity-scaling RNEMD |
| 582 |
< |
predicted comparable shear viscosities to swap RNEMD method. The |
| 583 |
< |
average molecular momentum gradients of these samples are shown in |
| 584 |
< |
Figure \ref{shear} (a) and (b). |
| 582 |
> |
predicted similar values for shear viscosities to the swapping RNEMD |
| 583 |
> |
method. The average molecular momentum gradients of these samples are |
| 584 |
> |
shown in the upper two panels of Fig. \ref{shear}. |
| 585 |
|
|
| 586 |
|
\begin{figure} |
| 587 |
|
\includegraphics[width=\linewidth]{shear} |
| 595 |
|
|
| 596 |
|
Observations of the three one-dimensional temperatures in each of the |
| 597 |
|
slabs shows that NIVS-RNEMD does produce some thermal anisotropy, |
| 598 |
< |
particularly in the hot and cold slabs. Figure \ref{shear} (c) |
| 599 |
< |
indicates that with a relatively large imposed momentum flux, |
| 600 |
< |
$j_z(p_x)$, the $x$ direction approaches a different temperature from |
| 601 |
< |
the $y$ and $z$ directions in both the hot and cold bins. This is an |
| 602 |
< |
artifact of the scaling constraints. After the momentum gradient has |
| 603 |
< |
been established, $P_c^x < 0$. Consequently, the scaling factor $x$ |
| 604 |
< |
is nearly always greater than one in order to satisfy the constraints. |
| 605 |
< |
This will continually increase the kinetic energy in $x$-dimension, |
| 606 |
< |
$K_c^x$. If there is not enough time for the kinetic energy to |
| 607 |
< |
exchange among different directions and different slabs, the system |
| 608 |
< |
will exhibit the observed thermal anisotropy in the hot and cold bins. |
| 598 |
> |
particularly in the hot and cold slabs. Note that these temperature |
| 599 |
> |
measurements have been taken into account of the kinetic energy |
| 600 |
> |
contributed by the slab field velocity. However, this contribution has |
| 601 |
> |
only a comparable order of magnitude to the errors of data, and does |
| 602 |
> |
not significantly affect our observation. The lower panel of Fig. |
| 603 |
> |
\ref{shear} indicates that with a relatively large imposed momentum |
| 604 |
> |
flux, $j_z(p_x)$, the $x$ direction approaches a different temperature |
| 605 |
> |
from the $y$ and $z$ directions in both the hot and cold bins. This |
| 606 |
> |
is an artifact of the scaling constraints. After a momentum gradient |
| 607 |
> |
has been established, $P_c^x$ is always less than zero. Consequently, |
| 608 |
> |
the scaling factor $x$ is always greater than one in order to satisfy |
| 609 |
> |
the constraints. This will continually increase the kinetic energy in |
| 610 |
> |
$x$-dimension, $K_c^x$. If there is not enough time for the kinetic |
| 611 |
> |
energy to exchange among different directions and different slabs, the |
| 612 |
> |
system will exhibit the observed thermal anisotropy in the hot and |
| 613 |
> |
cold bins. |
| 614 |
|
|
| 615 |
|
Although results between scaling and swapping methods are comparable, |
| 616 |
|
the inherent temperature anisotropy does make NIVS-RNEMD method less |
| 617 |
|
attractive than swapping RNEMD for shear viscosity calculations. We |
| 618 |
< |
note that this problem appears only when momentum flux is applied, and |
| 619 |
< |
does not appear in thermal flux calculations. |
| 618 |
> |
note that this problem appears only when a large {\it linear} momentum |
| 619 |
> |
flux is applied, and does not appear in {\it thermal} flux |
| 620 |
> |
calculations. |
| 621 |
|
|
| 622 |
|
\subsection{Bulk SPC/E water} |
| 623 |
|
|
| 624 |
|
We compared the thermal conductivity of SPC/E water using NIVS-RNEMD |
| 625 |
|
to the work of Bedrov {\it et al.}\cite{Bedrov:2000}, which employed |
| 626 |
|
the original swapping RNEMD method. Bedrov {\it et |
| 627 |
< |
al.}\cite{Bedrov:2000} argued that exchange of the molecule |
| 628 |
< |
center-of-mass velocities instead of single atom velocities in a |
| 629 |
< |
molecule conserves the total kinetic energy and linear momentum. This |
| 630 |
< |
principle is also adopted Fin our simulations. Scaling was applied to |
| 631 |
< |
the center-of-mass velocities of the rigid bodies of SPC/E model water |
| 619 |
< |
molecules. |
| 627 |
> |
al.}\cite{Bedrov:2000} argued that exchange of the molecular |
| 628 |
> |
center-of-mass velocities instead of single atom velocities conserves |
| 629 |
> |
the total kinetic energy and linear momentum. This principle is also |
| 630 |
> |
adopted in our simulations. Scaling was applied to the center-of-mass |
| 631 |
> |
velocities of SPC/E water molecules. |
| 632 |
|
|
| 633 |
|
To construct the simulations, a simulation box consisting of 1000 |
| 634 |
|
molecules were first equilibrated under ambient pressure and |
| 639 |
|
100~ps period under constant-NVE conditions without any momentum flux. |
| 640 |
|
Another 100~ps was allowed to stabilize the system with an imposed |
| 641 |
|
momentum transfer to create a gradient, and 1~ns was allotted for data |
| 642 |
< |
collection under RNEMD. |
| 642 |
> |
collection under RNEMD. Total system energy is recorded to ensure that |
| 643 |
> |
it is not drifted noticeably without a thermostat although |
| 644 |
> |
electrostatic interactions are involved. |
| 645 |
|
|
| 646 |
|
In our simulations, the established temperature gradients were similar |
| 647 |
|
to the previous work. Our simulation results at 318K are in good |
| 731 |
|
in volume. Following adjustment of the box volume, equilibrations in |
| 732 |
|
both the canonical and microcanonical ensembles were carried out. With |
| 733 |
|
the simulation cell divided evenly into 10 slabs, different thermal |
| 734 |
< |
gradients were established by applying a set of target thermal |
| 721 |
< |
transfer fluxes. |
| 734 |
> |
gradients were established by applying a set of target thermal fluxes. |
| 735 |
|
|
| 736 |
|
The results for the thermal conductivity of gold are shown in Table |
| 737 |
|
\ref{AuThermal}. In these calculations, the end and middle slabs were |
| 738 |
< |
excluded in thermal gradient linear regession. EAM predicts slightly |
| 739 |
< |
larger thermal conductivities than QSC. However, both values are |
| 740 |
< |
smaller than experimental value by a factor of more than 200. This |
| 738 |
> |
excluded from the thermal gradient linear regession. EAM predicts |
| 739 |
> |
slightly larger thermal conductivities than QSC. However, both values |
| 740 |
> |
are smaller than experimental value by a factor of more than 200. This |
| 741 |
|
behavior has been observed previously by Richardson and Clancy, and |
| 742 |
|
has been attributed to the lack of electronic contribution in these |
| 743 |
|
force fields.\cite{Clancy:1992} It should be noted that the density of |
| 798 |
|
relax under ambient temperature and pressure. A separate (but |
| 799 |
|
identically sized) box of SPC/E water was also equilibrated at ambient |
| 800 |
|
conditions. The two boxes were combined by removing all water |
| 801 |
< |
molecules within 3 \AA radius of any gold atom. The final |
| 801 |
> |
molecules within 3 \AA~ radius of any gold atom. The final |
| 802 |
|
configuration contained 1862 SPC/E water molecules. |
| 803 |
|
|
| 804 |
|
The Spohr potential was adopted in depicting the interaction between |
| 812 |
|
phase.\cite{Cahill:793,RevModPhys.61.605} Experiments on the thermal |
| 813 |
|
conductivity of gold nanoparticles and nanorods in solvent and in |
| 814 |
|
glass cages have predicted values for $G$ between 100 and 350 |
| 815 |
< |
(MW/m$^2$/K). The experiments typically have multiple gold surfaces |
| 816 |
< |
that have been protected by a capping agent (citrate or CTAB) or which |
| 817 |
< |
are in direct contact with various glassy solids. In these cases, the |
| 818 |
< |
acoustic impedance mismatch would be substantially reduced, leading to |
| 819 |
< |
much higher interfacial conductances. It is also possible, however, |
| 820 |
< |
that the lack of electronic effects that gives rise to the low thermal |
| 821 |
< |
conductivity of EAM gold is also causing a low reading for this |
| 822 |
< |
particular interface. |
| 815 |
> |
(MW/m$^2$/K), two orders of magnitude larger than the value reported |
| 816 |
> |
here. The experiments typically have multiple surfaces that have been |
| 817 |
> |
protected by ionic surfactants, either |
| 818 |
> |
citrate\cite{Wilson:2002uq,plech:195423} or cetyltrimethylammonium |
| 819 |
> |
bromide (CTAB), or which are in direct contact with various glassy |
| 820 |
> |
solids. In these cases, the acoustic impedance mismatch would be |
| 821 |
> |
substantially reduced, leading to much higher interfacial |
| 822 |
> |
conductances. It is also possible, however, that the lack of |
| 823 |
> |
electronic effects that gives rise to the low thermal conductivity of |
| 824 |
> |
EAM gold is also causing a low reading for this particular interface. |
| 825 |
|
|
| 826 |
< |
Under this low thermal conductance, both gold and water phase have |
| 827 |
< |
sufficient time to eliminate temperature difference inside |
| 828 |
< |
respectively (Figure \ref{interface} b). With indistinguishable |
| 829 |
< |
temperature difference within respective phase, it is valid to assume |
| 830 |
< |
that the temperature difference between gold and water on surface |
| 831 |
< |
would be approximately the same as the difference between the gold and |
| 832 |
< |
water phase. This assumption enables convenient calculation of $G$ |
| 818 |
< |
using Eq. \ref{interfaceCalc} instead of measuring temperatures of |
| 819 |
< |
thin layer of water and gold close enough to surface, which would have |
| 820 |
< |
greater fluctuation and lower accuracy. Reported results (Table |
| 821 |
< |
\ref{interfaceRes}) are of two orders of magnitude smaller than our |
| 822 |
< |
calculations on homogeneous systems, and thus have larger relative |
| 823 |
< |
errors than our calculation results on homogeneous systems. |
| 826 |
> |
Under this low thermal conductance, both gold and water phases have |
| 827 |
> |
sufficient time to eliminate local temperature differences (Fig. |
| 828 |
> |
\ref{interface}). With flat thermal profiles within each phase, it is |
| 829 |
> |
valid to assume that the temperature difference between gold and water |
| 830 |
> |
surfaces would be approximately the same as the difference between the |
| 831 |
> |
gold and water bulk regions. This assumption enables convenient |
| 832 |
> |
calculation of $G$ using Eq. \ref{interfaceCalc}. |
| 833 |
|
|
| 834 |
|
\begin{figure} |
| 835 |
|
\includegraphics[width=\linewidth]{interface} |
| 869 |
|
|
| 870 |
|
|
| 871 |
|
\section{Conclusions} |
| 872 |
< |
NIVS-RNEMD simulation method is developed and tested on various |
| 873 |
< |
systems. Simulation results demonstrate its validity in thermal |
| 874 |
< |
conductivity calculations, from Lennard-Jones fluid to multi-atom |
| 875 |
< |
molecule like water and metal crystals. NIVS-RNEMD improves |
| 876 |
< |
non-Boltzmann-Maxwell distributions, which exist inb previous RNEMD |
| 877 |
< |
methods. Furthermore, it develops a valid means for unphysical thermal |
| 872 |
> |
|
| 873 |
> |
Our simulations demonstrate that validity of non-isotropic velocity |
| 874 |
> |
scaling (NIVS) in RNEMD calculations of thermal conductivity in atomic |
| 875 |
> |
and molecular liquids and solids. NIVS-RNEMD improves the problematic |
| 876 |
> |
velocity distributions which can develop in other RNEMD methods. |
| 877 |
> |
Furthermore, it provides a means for carrying out non-physical thermal |
| 878 |
|
transfer between different species of molecules, and thus extends its |
| 879 |
< |
applicability to interfacial systems. Our calculation of gold/water |
| 880 |
< |
interfacial thermal conductivity demonstrates this advantage over |
| 881 |
< |
previous RNEMD methods. NIVS-RNEMD has also limited application on |
| 882 |
< |
shear viscosity calculations, but could cause temperature difference |
| 883 |
< |
among different dimensions under high momentum flux. Modification is |
| 884 |
< |
necessary to extend the applicability of NIVS-RNEMD in shear viscosity |
| 885 |
< |
calculations. |
| 879 |
> |
applicability to interfacial systems. Our calculation of the gold / |
| 880 |
> |
water interfacial thermal conductivity demonstrates this advantage |
| 881 |
> |
over previous RNEMD methods. NIVS-RNEMD also has limited applications |
| 882 |
> |
for shear viscosity calculations, but has the potential to cause |
| 883 |
> |
temperature anisotropy under high momentum fluxes. Further work will |
| 884 |
> |
be necessary to eliminate the one-dimensional heating if shear |
| 885 |
> |
viscosities are required. |
| 886 |
|
|
| 887 |
|
\section{Acknowledgments} |
| 888 |
|
The authors would like to thank Craig Tenney and Ed Maginn for many |
