group/interfacial/interfacial.tex

\documentclass[11pt]{article}
\usepackage{amsmath}
\usepackage{amssymb}
\usepackage{setspace}
\usepackage{endfloat}
\usepackage{caption}
%\usepackage{tabularx}
\usepackage{graphicx}
\usepackage{multirow}
%\usepackage{booktabs}
%\usepackage{bibentry}
%\usepackage{mathrsfs}
%\usepackage[ref]{overcite}
\usepackage[square, comma, sort&compress]{natbib}
\usepackage{url}
\pagestyle{plain} \pagenumbering{arabic} \oddsidemargin 0.0cm
\evensidemargin 0.0cm \topmargin -21pt \headsep 10pt \textheight
9.0in \textwidth 6.5in \brokenpenalty=10000

% double space list of tables and figures
\AtBeginDelayedFloats{\renewcommand{\baselinestretch}{1.66}}
\setlength{\abovecaptionskip}{20 pt}
\setlength{\belowcaptionskip}{30 pt}

%\renewcommand\citemid{\ } % no comma in optional referenc note
\bibpunct{[}{]}{,}{s}{}{;}
\bibliographystyle{aip}

\begin{document}

\title{Simulating interfacial thermal conductance at metal-solvent
  interfaces: the role of chemical capping agents}

\author{Shenyu Kuang and J. Daniel
Gezelter\footnote{Corresponding author. \ Electronic mail: gezelter@nd.edu} \\
Department of Chemistry and Biochemistry,\\
University of Notre Dame\\
Notre Dame, Indiana 46556}

\date{\today}

\maketitle 

\begin{doublespace}

\begin{abstract}

We have developed a Non-Isotropic Velocity Scaling algorithm for
setting up and maintaining stable thermal gradients in non-equilibrium
molecular dynamics simulations. This approach effectively imposes
unphysical thermal flux even between particles of different
identities, conserves linear momentum and kinetic energy, and
minimally perturbs the velocity profile of a system when compared with
previous RNEMD methods. We have used this method to simulate thermal
conductance at metal / organic solvent interfaces both with and
without the presence of thiol-based capping agents.  We obtained
values comparable with experimental values, and observed significant
conductance enhancement with the presence of capping agents. Computed
power spectra indicate the acoustic impedance mismatch between metal
and liquid phase is greatly reduced by the capping agents and thus
leads to higher interfacial thermal transfer efficiency.

\end{abstract}

\newpage

%\narrowtext

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%                          BODY OF TEXT
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\section{Introduction}

Interfacial thermal conductance is extensively studied both
experimentally and computationally, and systems with interfaces
present are generally heterogeneous. Although interfaces are commonly
barriers to heat transfer, it has been
reported\cite{doi:10.1021/la904855s} that under specific circustances,
e.g. with certain capping agents present on the surface, interfacial
conductance can be significantly enhanced. However, heat conductance
of molecular and nano-scale interfaces will be affected by the
chemical details of the surface and is challenging to
experimentalist. The lower thermal flux through interfaces is even
more difficult to measure with EMD and forward NEMD simulation
methods. Therefore, developing good simulation methods will be
desirable in order to investigate thermal transport across interfaces.

Recently, we have developed the Non-Isotropic Velocity Scaling (NIVS)
algorithm for RNEMD simulations\cite{kuang:164101}. This algorithm
retains the desirable features of RNEMD (conservation of linear
momentum and total energy, compatibility with periodic boundary
conditions) while establishing true thermal distributions in each of
the two slabs. Furthermore, it allows more effective thermal exchange
between particles of different identities, and thus enables extensive
study of interfacial conductance.

\section{Methodology}
\subsection{Algorithm}
There have been many algorithms for computing thermal conductivity
using molecular dynamics simulations. However, interfacial conductance
is at least an order of magnitude smaller. This would make the
calculation even more difficult for those slowly-converging
equilibrium methods. Imposed-flux non-equilibrium
methods\cite{MullerPlathe:1997xw} have the flux set {\it a priori} and
the response of temperature or momentum gradients are easier to
measure than the flux, if unknown, and thus, is a preferable way to
the forward NEMD methods. Although the momentum swapping approach for
flux-imposing can be used for exchanging energy between particles of
different identity, the kinetic energy transfer efficiency is affected
by the mass difference between the particles, which limits its
application on heterogeneous interfacial systems.

The non-isotropic velocity scaling (NIVS)\cite{kuang:164101} approach in
non-equilibrium MD simulations is able to impose relatively large
kinetic energy flux without obvious perturbation to the velocity
distribution of the simulated systems. Furthermore, this approach has
the advantage in heterogeneous interfaces in that kinetic energy flux
can be applied between regions of particles of arbitary identity, and
the flux quantity is not restricted by particle mass difference.

The NIVS algorithm scales the velocity vectors in two separate regions
of a simulation system with respective diagonal scaling matricies. To
determine these scaling factors in the matricies, a set of equations
including linear momentum conservation and kinetic energy conservation
constraints and target momentum/energy flux satisfaction is
solved. With the scaling operation applied to the system in a set
frequency, corresponding momentum/temperature gradients can be built,
which can be used for computing transportation properties and other
applications related to momentum/temperature gradients. The NIVS
algorithm conserves momenta and energy and does not depend on an
external thermostat. 

(wondering how much detail of algorithm should be put here...)

\subsection{Force Field Parameters}
Our simulation systems consists of metal gold lattice slab solvated by
organic solvents. In order to study the role of capping agents in
interfacial thermal conductance, butanethiol is chosen to cover gold
surfaces in comparison to no capping agent present.

The Au-Au interactions in metal lattice slab is described by the
quantum Sutton-Chen (QSC) formulation.\cite{PhysRevB.59.3527} The QSC
potentials include zero-point quantum corrections and are
reparametrized for accurate surface energies compared to the
Sutton-Chen potentials\cite{Chen90}. 

Straight chain {\it n}-hexane and aromatic toluene are respectively
used as solvents. For hexane, both United-Atom\cite{TraPPE-UA.alkanes}
and All-Atom\cite{OPLSAA} force fields are used for comparison; for
toluene, United-Atom\cite{TraPPE-UA.alkylbenzenes} force fields are
used with rigid body constraints applied. (maybe needs more details
about rigid body)

Buatnethiol molecules are used as capping agent for some of our
simulations. United-Atom\cite{TraPPE-UA.thiols} and All-Atom models
are respectively used corresponding to the force field type of
solvent.

To describe the interactions between metal Au and non-metal capping
agent and solvent, we refer to Vlugt\cite{vlugt:cpc2007154} and derive
other interactions which are not parametrized in their work. (can add
hautman and klein's paper here and more discussion; need to put
aromatic-metal interaction approximation here)

[TABULATED FORCE FIELD PARAMETERS NEEDED]

\section{Computational Details}
\subsection{System Geometry}
Our simulation systems consists of a lattice Au slab with the (111)
surface perpendicular to the $z$-axis, and a solvent layer between the
periodic Au slabs along the $z$-axis. To set up the interfacial
system, the Au slab is first equilibrated without solvent under room
pressure and a desired temperature. After the metal slab is
equilibrated, United-Atom or All-Atom butanethiols are replicated on
the Au surface, each occupying the (??) among three Au atoms, and is
equilibrated under NVT ensemble. According to (CITATION), the maximal
thiol capacity on Au surface is $1/3$ of the total number of surface
Au atoms. 

\cite{packmol}

\subsection{Simulation Parameters}

When the interfacial conductance is {\it not} small, there are two
ways to define $G$. If we assume the temperature is discretely
different on two sides of the interface, $G$ can be calculated with
the thermal flux applied $J$ and the temperature difference measured
$\Delta T$ as:
\begin{equation}
G=\frac{J}{\Delta T}
\label{discreteG}
\end{equation}
We can as well assume a continuous temperature profile along the
thermal gradient axis $z$ and define $G$ as the change of bulk thermal
conductivity $\lambda$ at a defined interfacial point:
\begin{equation}
G^\prime = \Big|\frac{\partial\lambda}{\partial z}\Big|
         = \Big|\frac{\partial}{\partial z}\left(-J_z\Big/
           \left(\frac{\partial T}{\partial z}\right)\right)\Big|
         = J_z\Big|\frac{\partial^2 T}{\partial z^2}\Big|
         \Big/\left(\frac{\partial T}{\partial z}\right)^2
\label{derivativeG}
\end{equation}
With the temperature profile obtained from simulations, one is able to
approximate the first and second derivatives of $T$ with finite
difference method and thus calculate $G^\prime$.

In what follows, both definitions are used for calculation and comparison.

\section{Results}
\subsection{Toluene Solvent}

The simulations follow a protocol similar to the previous gold/water
interfacial systems. The results (Table \ref{AuThiolToluene}) show a
significant conductance enhancement compared to the gold/water
interface without capping agent and agree with available experimental
data. This indicates that the metal-metal potential, though not
predicting an accurate bulk metal thermal conductivity, does not
greatly interfere with the simulation of the thermal conductance
behavior across a non-metal interface. The solvent model is not
particularly volatile, so the simulation cell does not expand
significantly under higher temperature. We did not observe a
significant conductance decrease when the temperature was increased to
300K. The results show that the two definitions used for $G$ yield
comparable values, though $G^\prime$ tends to be smaller.

\begin{table*}
  \begin{minipage}{\linewidth}
    \begin{center}
      \caption{Computed interfacial thermal conductivity ($G$ and
        $G^\prime$) values for the Au/butanethiol/toluene interface at
        different temperatures using a range of energy fluxes.}
      
      \begin{tabular}{cccc}
        \hline\hline
        $\langle T\rangle$ & $J_z$ & $G$ & $G^\prime$ \\
        (K) & (GW/m$^2$) & \multicolumn{2}{c}{(MW/m$^2$/K)} \\
        \hline
        200 & 1.86 & 180 & 135 \\
            & 2.15 & 204 & 113 \\
            & 3.93 & 175 & 114 \\
        300 & 1.91 & 143 & 125 \\
            & 4.19 & 134 & 113 \\
        \hline\hline
      \end{tabular}
      \label{AuThiolToluene}
    \end{center}
  \end{minipage}
\end{table*}

\subsection{Hexane Solvent}

Using the united-atom model, different coverages of capping agent,
temperatures of simulations and numbers of solvent molecules were all
investigated and Table \ref{AuThiolHexaneUA} shows the results of
these computations. The number of hexane molecules in our simulations
does not affect the calculations significantly. However, a very long
length scale for the thermal gradient axis ($z$) may cause excessively
hot or cold temperatures in the middle of the solvent region and lead
to undesired phenomena such as solvent boiling or freezing, while too
few solvent molecules would change the normal behavior of the liquid
phase. Our $N_{hexane}$ values were chosen to ensure that these
extreme cases did not happen to our simulations.

Table \ref{AuThiolHexaneUA} enables direct comparison between
different coverages of capping agent, when other system parameters are
held constant. With high coverage of butanethiol on the gold surface,
the interfacial thermal conductance is enhanced
significantly. Interestingly, a slightly lower butanethiol coverage
leads to a moderately higher conductivity. This is probably due to
more solvent/capping agent contact when butanethiol molecules are
not densely packed, which enhances the interactions between the two
phases and lowers the thermal transfer barrier of this interface.
% [COMPARE TO AU/WATER IN PAPER]

It is also noted that the overall simulation temperature is another
factor that affects the interfacial thermal conductance. One
possibility of this effect may be rooted in the decrease in density of
the liquid phase. We observed that when the average temperature
increases from 200K to 250K, the bulk hexane density becomes lower
than experimental value, as the system is equilibrated under NPT
ensemble. This leads to lower contact between solvent and capping
agent, and thus lower conductivity.

Conductivity values are more difficult to obtain under higher
temperatures. This is because the Au surface tends to undergo
reconstructions in relatively high temperatures. Surface Au atoms can
migrate outward to reach higher Au-S contact; and capping agent
molecules can be embedded into the surface Au layer due to the same
driving force. This phenomenon agrees with experimental
results\cite{doi:10.1021/j100035a033,doi:10.1021/la026493y}. A surface
fully covered in capping agent is more susceptible to reconstruction,
possibly because fully coverage prevents other means of capping agent
relaxation, such as migration to an empty neighbor three-fold site.

%MAY ADD MORE DATA TO TABLE
\begin{table*}
  \begin{minipage}{\linewidth}
    \begin{center}
      \caption{Computed interfacial thermal conductivity ($G$ and
        $G^\prime$) values for the Au/butanethiol/hexane interface
        with united-atom model and different capping agent coverage
        and solvent molecule numbers at different temperatures using a
        range of energy fluxes.}
      
      \begin{tabular}{cccccc}
        \hline\hline
        Thiol & $\langle T\rangle$ & & $J_z$ & $G$ & $G^\prime$ \\
        coverage (\%) & (K) & $N_{hexane}$ & (GW/m$^2$) &
        \multicolumn{2}{c}{(MW/m$^2$/K)} \\
        \hline
        0.0   & 200 & 200 & 0.96 & 43.3 & 42.7 \\
              &     &     & 1.91 & 45.7 & 42.9 \\
              &     & 166 & 0.96 & 43.1 & 53.4 \\
        88.9  & 200 & 166 & 1.94 & 172  & 108  \\
        100.0 & 250 & 200 & 0.96 & 81.8 & 67.0 \\
              &     & 166 & 0.98 & 79.0 & 62.9 \\
              &     &     & 1.44 & 76.2 & 64.8 \\
              & 200 & 200 & 1.92 & 129  & 87.3 \\
              &     &     & 1.93 & 131  & 77.5 \\
              &     & 166 & 0.97 & 115  & 69.3 \\
              &     &     & 1.94 & 125  & 87.1 \\
        \hline\hline
      \end{tabular}
      \label{AuThiolHexaneUA}
    \end{center}
  \end{minipage}
\end{table*}

For the all-atom model, the liquid hexane phase was not stable under NPT
conditions. Therefore, the simulation length scale parameters are
adopted from previous equilibration results of the united-atom model
at 200K. Table \ref{AuThiolHexaneAA} shows the results of these
simulations. The conductivity values calculated with full capping
agent coverage are substantially larger than observed in the
united-atom model, and is even higher than predicted by
experiments. It is possible that our parameters for metal-non-metal
particle interactions lead to an overestimate of the interfacial
thermal conductivity, although the active C-H vibrations in the
all-atom model (which should not be appreciably populated at normal
temperatures) could also account for this high conductivity. The major
thermal transfer barrier of Au/butanethiol/hexane interface is between
the liquid phase and the capping agent, so extra degrees of freedom
such as the C-H vibrations could enhance heat exchange between these
two phases and result in a much higher conductivity.

\begin{table*}
  \begin{minipage}{\linewidth}
    \begin{center}
      
      \caption{Computed interfacial thermal conductivity ($G$ and
        $G^\prime$) values for the Au/butanethiol/hexane interface
        with all-atom model and different capping agent coverage at
        200K using a range of energy fluxes.}
      
      \begin{tabular}{cccc}
        \hline\hline
        Thiol & $J_z$ & $G$ & $G^\prime$ \\
        coverage (\%) & (GW/m$^2$) & \multicolumn{2}{c}{(MW/m$^2$/K)} \\
        \hline
        0.0   & 0.95 & 28.5 & 27.2 \\
              & 1.88 & 30.3 & 28.9 \\
        100.0 & 2.87 & 551  & 294  \\
              & 3.81 & 494  & 193  \\
        \hline\hline
      \end{tabular}
      \label{AuThiolHexaneAA}
    \end{center}
  \end{minipage}
\end{table*}

%subsubsection{Vibrational spectrum study on conductance mechanism}
To investigate the mechanism of this interfacial thermal conductance,
the vibrational spectra of various gold systems were obtained and are
shown as in the upper panel of Fig. \ref{vibration}. To obtain these
spectra, one first runs a simulation in the NVE ensemble and collects
snapshots of configurations; these configurations are used to compute
the velocity auto-correlation functions, which is used to construct a
power spectrum via a Fourier transform. The gold surfaces covered by
butanethiol molecules exhibit an additional peak observed at a
frequency of $\sim$170cm$^{-1}$, which is attributed to the vibration
of the S-Au bond. This vibration enables efficient thermal transport
from surface Au atoms to the capping agents. Simultaneously, as shown
in the lower panel of Fig. \ref{vibration}, the large overlap of the
vibration spectra of butanethiol and hexane in the all-atom model,
including the C-H vibration, also suggests high thermal exchange
efficiency. The combination of these two effects produces the drastic
interfacial thermal conductance enhancement in the all-atom model.

\begin{figure}
\includegraphics[width=\linewidth]{vibration}
\caption{Vibrational spectra obtained for gold in different
  environments (upper panel) and for Au/thiol/hexane simulation in
  all-atom model (lower panel).}
\label{vibration}
\end{figure}
% 600dpi, letter size. too large?


\section{Acknowledgments}
Support for this project was provided by the National Science
Foundation under grant CHE-0848243. Computational time was provided by
the Center for Research Computing (CRC) at the University of Notre
Dame.  \newpage

\bibliography{interfacial}

\end{doublespace}
\end{document}

Revision:	3725
Committed:	Mon May 9 19:08:08 2011 UTC (13 years, 1 month ago) by skuang
Content type:	application/x-tex
File size:	18396 byte(s)
Log Message:	integrated candidacy document into paper
#	User	Rev	Content
1	gezelter	3717	\documentclass[11pt]{article}
2			\usepackage{amsmath}
3			\usepackage{amssymb}
4			\usepackage{setspace}
5			\usepackage{endfloat}
6			\usepackage{caption}
7			%\usepackage{tabularx}
8			\usepackage{graphicx}
9			\usepackage{multirow}
10			%\usepackage{booktabs}
11			%\usepackage{bibentry}
12			%\usepackage{mathrsfs}
13			%\usepackage[ref]{overcite}
14			\usepackage[square, comma, sort&compress]{natbib}
15			\usepackage{url}
16			\pagestyle{plain} \pagenumbering{arabic} \oddsidemargin 0.0cm
17			\evensidemargin 0.0cm \topmargin -21pt \headsep 10pt \textheight
18			9.0in \textwidth 6.5in \brokenpenalty=10000
19
20			% double space list of tables and figures
21			\AtBeginDelayedFloats{\renewcommand{\baselinestretch}{1.66}}
22			\setlength{\abovecaptionskip}{20 pt}
23			\setlength{\belowcaptionskip}{30 pt}
24
25			%\renewcommand\citemid{\ } % no comma in optional referenc note
26			\bibpunct{[}{]}{,}{s}{}{;}
27			\bibliographystyle{aip}
28
29			\begin{document}
30
31			\title{Simulating interfacial thermal conductance at metal-solvent
32			interfaces: the role of chemical capping agents}
33
34			\author{Shenyu Kuang and J. Daniel
35			Gezelter\footnote{Corresponding author. \ Electronic mail: gezelter@nd.edu} \\
36			Department of Chemistry and Biochemistry,\\
37			University of Notre Dame\\
38			Notre Dame, Indiana 46556}
39
40			\date{\today}
41
42			\maketitle
43
44			\begin{doublespace}
45
46			\begin{abstract}
47	skuang	3725
48			We have developed a Non-Isotropic Velocity Scaling algorithm for
49			setting up and maintaining stable thermal gradients in non-equilibrium
50			molecular dynamics simulations. This approach effectively imposes
51			unphysical thermal flux even between particles of different
52			identities, conserves linear momentum and kinetic energy, and
53			minimally perturbs the velocity profile of a system when compared with
54			previous RNEMD methods. We have used this method to simulate thermal
55			conductance at metal / organic solvent interfaces both with and
56			without the presence of thiol-based capping agents. We obtained
57			values comparable with experimental values, and observed significant
58			conductance enhancement with the presence of capping agents. Computed
59			power spectra indicate the acoustic impedance mismatch between metal
60			and liquid phase is greatly reduced by the capping agents and thus
61			leads to higher interfacial thermal transfer efficiency.
62
63	gezelter	3717	\end{abstract}
64
65			\newpage
66
67			%\narrowtext
68
69			%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
70			% BODY OF TEXT
71			%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
72
73			\section{Introduction}
74
75	skuang	3725	Interfacial thermal conductance is extensively studied both
76			experimentally and computationally, and systems with interfaces
77			present are generally heterogeneous. Although interfaces are commonly
78			barriers to heat transfer, it has been
79			reported\cite{doi:10.1021/la904855s} that under specific circustances,
80			e.g. with certain capping agents present on the surface, interfacial
81			conductance can be significantly enhanced. However, heat conductance
82			of molecular and nano-scale interfaces will be affected by the
83			chemical details of the surface and is challenging to
84			experimentalist. The lower thermal flux through interfaces is even
85			more difficult to measure with EMD and forward NEMD simulation
86			methods. Therefore, developing good simulation methods will be
87			desirable in order to investigate thermal transport across interfaces.
88	gezelter	3717
89	skuang	3725	Recently, we have developed the Non-Isotropic Velocity Scaling (NIVS)
90			algorithm for RNEMD simulations\cite{kuang:164101}. This algorithm
91			retains the desirable features of RNEMD (conservation of linear
92			momentum and total energy, compatibility with periodic boundary
93			conditions) while establishing true thermal distributions in each of
94			the two slabs. Furthermore, it allows more effective thermal exchange
95			between particles of different identities, and thus enables extensive
96			study of interfacial conductance.
97
98	skuang	3721	\section{Methodology}
99			\subsection{Algorithm}
100			There have been many algorithms for computing thermal conductivity
101			using molecular dynamics simulations. However, interfacial conductance
102			is at least an order of magnitude smaller. This would make the
103			calculation even more difficult for those slowly-converging
104			equilibrium methods. Imposed-flux non-equilibrium
105			methods\cite{MullerPlathe:1997xw} have the flux set {\it a priori} and
106			the response of temperature or momentum gradients are easier to
107			measure than the flux, if unknown, and thus, is a preferable way to
108			the forward NEMD methods. Although the momentum swapping approach for
109			flux-imposing can be used for exchanging energy between particles of
110			different identity, the kinetic energy transfer efficiency is affected
111			by the mass difference between the particles, which limits its
112			application on heterogeneous interfacial systems.
113
114			The non-isotropic velocity scaling (NIVS)\cite{kuang:164101} approach in
115			non-equilibrium MD simulations is able to impose relatively large
116			kinetic energy flux without obvious perturbation to the velocity
117			distribution of the simulated systems. Furthermore, this approach has
118			the advantage in heterogeneous interfaces in that kinetic energy flux
119			can be applied between regions of particles of arbitary identity, and
120			the flux quantity is not restricted by particle mass difference.
121
122			The NIVS algorithm scales the velocity vectors in two separate regions
123			of a simulation system with respective diagonal scaling matricies. To
124			determine these scaling factors in the matricies, a set of equations
125			including linear momentum conservation and kinetic energy conservation
126			constraints and target momentum/energy flux satisfaction is
127			solved. With the scaling operation applied to the system in a set
128			frequency, corresponding momentum/temperature gradients can be built,
129			which can be used for computing transportation properties and other
130			applications related to momentum/temperature gradients. The NIVS
131			algorithm conserves momenta and energy and does not depend on an
132			external thermostat.
133
134			(wondering how much detail of algorithm should be put here...)
135
136	skuang	3725	\subsection{Force Field Parameters}
137	skuang	3721	Our simulation systems consists of metal gold lattice slab solvated by
138			organic solvents. In order to study the role of capping agents in
139			interfacial thermal conductance, butanethiol is chosen to cover gold
140			surfaces in comparison to no capping agent present.
141
142			The Au-Au interactions in metal lattice slab is described by the
143			quantum Sutton-Chen (QSC) formulation.\cite{PhysRevB.59.3527} The QSC
144			potentials include zero-point quantum corrections and are
145			reparametrized for accurate surface energies compared to the
146			Sutton-Chen potentials\cite{Chen90}.
147
148			Straight chain {\it n}-hexane and aromatic toluene are respectively
149			used as solvents. For hexane, both United-Atom\cite{TraPPE-UA.alkanes}
150			and All-Atom\cite{OPLSAA} force fields are used for comparison; for
151			toluene, United-Atom\cite{TraPPE-UA.alkylbenzenes} force fields are
152			used with rigid body constraints applied. (maybe needs more details
153			about rigid body)
154
155			Buatnethiol molecules are used as capping agent for some of our
156			simulations. United-Atom\cite{TraPPE-UA.thiols} and All-Atom models
157			are respectively used corresponding to the force field type of
158			solvent.
159
160			To describe the interactions between metal Au and non-metal capping
161			agent and solvent, we refer to Vlugt\cite{vlugt:cpc2007154} and derive
162			other interactions which are not parametrized in their work. (can add
163			hautman and klein's paper here and more discussion; need to put
164			aromatic-metal interaction approximation here)
165
166	skuang	3725	[TABULATED FORCE FIELD PARAMETERS NEEDED]
167
168	skuang	3721	\section{Computational Details}
169	skuang	3725	\subsection{System Geometry}
170	skuang	3721	Our simulation systems consists of a lattice Au slab with the (111)
171			surface perpendicular to the $z$-axis, and a solvent layer between the
172			periodic Au slabs along the $z$-axis. To set up the interfacial
173			system, the Au slab is first equilibrated without solvent under room
174			pressure and a desired temperature. After the metal slab is
175			equilibrated, United-Atom or All-Atom butanethiols are replicated on
176			the Au surface, each occupying the (??) among three Au atoms, and is
177			equilibrated under NVT ensemble. According to (CITATION), the maximal
178			thiol capacity on Au surface is $1/3$ of the total number of surface
179			Au atoms.
180
181			\cite{packmol}
182
183	skuang	3725	\subsection{Simulation Parameters}
184
185			When the interfacial conductance is {\it not} small, there are two
186			ways to define $G$. If we assume the temperature is discretely
187			different on two sides of the interface, $G$ can be calculated with
188			the thermal flux applied $J$ and the temperature difference measured
189			$\Delta T$ as:
190			\begin{equation}
191			G=\frac{J}{\Delta T}
192			\label{discreteG}
193			\end{equation}
194			We can as well assume a continuous temperature profile along the
195			thermal gradient axis $z$ and define $G$ as the change of bulk thermal
196			conductivity $\lambda$ at a defined interfacial point:
197			\begin{equation}
198			G^\prime = \Big\|\frac{\partial\lambda}{\partial z}\Big\|
199			= \Big\|\frac{\partial}{\partial z}\left(-J_z\Big/
200			\left(\frac{\partial T}{\partial z}\right)\right)\Big\|
201			= J_z\Big\|\frac{\partial^2 T}{\partial z^2}\Big\|
202			\Big/\left(\frac{\partial T}{\partial z}\right)^2
203			\label{derivativeG}
204			\end{equation}
205			With the temperature profile obtained from simulations, one is able to
206			approximate the first and second derivatives of $T$ with finite
207			difference method and thus calculate $G^\prime$.
208
209			In what follows, both definitions are used for calculation and comparison.
210
211			\section{Results}
212			\subsection{Toluene Solvent}
213
214			The simulations follow a protocol similar to the previous gold/water
215			interfacial systems. The results (Table \ref{AuThiolToluene}) show a
216			significant conductance enhancement compared to the gold/water
217			interface without capping agent and agree with available experimental
218			data. This indicates that the metal-metal potential, though not
219			predicting an accurate bulk metal thermal conductivity, does not
220			greatly interfere with the simulation of the thermal conductance
221			behavior across a non-metal interface. The solvent model is not
222			particularly volatile, so the simulation cell does not expand
223			significantly under higher temperature. We did not observe a
224			significant conductance decrease when the temperature was increased to
225			300K. The results show that the two definitions used for $G$ yield
226			comparable values, though $G^\prime$ tends to be smaller.
227
228			\begin{table*}
229			\begin{minipage}{\linewidth}
230			\begin{center}
231			\caption{Computed interfacial thermal conductivity ($G$ and
232			$G^\prime$) values for the Au/butanethiol/toluene interface at
233			different temperatures using a range of energy fluxes.}
234
235			\begin{tabular}{cccc}
236			\hline\hline
237			$\langle T\rangle$ & $J_z$ & $G$ & $G^\prime$ \\
238			(K) & (GW/m$^2$) & \multicolumn{2}{c}{(MW/m$^2$/K)} \\
239			\hline
240			200 & 1.86 & 180 & 135 \\
241			& 2.15 & 204 & 113 \\
242			& 3.93 & 175 & 114 \\
243			300 & 1.91 & 143 & 125 \\
244			& 4.19 & 134 & 113 \\
245			\hline\hline
246			\end{tabular}
247			\label{AuThiolToluene}
248			\end{center}
249			\end{minipage}
250			\end{table*}
251
252			\subsection{Hexane Solvent}
253
254			Using the united-atom model, different coverages of capping agent,
255			temperatures of simulations and numbers of solvent molecules were all
256			investigated and Table \ref{AuThiolHexaneUA} shows the results of
257			these computations. The number of hexane molecules in our simulations
258			does not affect the calculations significantly. However, a very long
259			length scale for the thermal gradient axis ($z$) may cause excessively
260			hot or cold temperatures in the middle of the solvent region and lead
261			to undesired phenomena such as solvent boiling or freezing, while too
262			few solvent molecules would change the normal behavior of the liquid
263			phase. Our $N_{hexane}$ values were chosen to ensure that these
264			extreme cases did not happen to our simulations.
265
266			Table \ref{AuThiolHexaneUA} enables direct comparison between
267			different coverages of capping agent, when other system parameters are
268			held constant. With high coverage of butanethiol on the gold surface,
269			the interfacial thermal conductance is enhanced
270			significantly. Interestingly, a slightly lower butanethiol coverage
271			leads to a moderately higher conductivity. This is probably due to
272			more solvent/capping agent contact when butanethiol molecules are
273			not densely packed, which enhances the interactions between the two
274			phases and lowers the thermal transfer barrier of this interface.
275			% [COMPARE TO AU/WATER IN PAPER]
276
277			It is also noted that the overall simulation temperature is another
278			factor that affects the interfacial thermal conductance. One
279			possibility of this effect may be rooted in the decrease in density of
280			the liquid phase. We observed that when the average temperature
281			increases from 200K to 250K, the bulk hexane density becomes lower
282			than experimental value, as the system is equilibrated under NPT
283			ensemble. This leads to lower contact between solvent and capping
284			agent, and thus lower conductivity.
285
286			Conductivity values are more difficult to obtain under higher
287			temperatures. This is because the Au surface tends to undergo
288			reconstructions in relatively high temperatures. Surface Au atoms can
289			migrate outward to reach higher Au-S contact; and capping agent
290			molecules can be embedded into the surface Au layer due to the same
291			driving force. This phenomenon agrees with experimental
292			results\cite{doi:10.1021/j100035a033,doi:10.1021/la026493y}. A surface
293			fully covered in capping agent is more susceptible to reconstruction,
294			possibly because fully coverage prevents other means of capping agent
295			relaxation, such as migration to an empty neighbor three-fold site.
296
297			%MAY ADD MORE DATA TO TABLE
298			\begin{table*}
299			\begin{minipage}{\linewidth}
300			\begin{center}
301			\caption{Computed interfacial thermal conductivity ($G$ and
302			$G^\prime$) values for the Au/butanethiol/hexane interface
303			with united-atom model and different capping agent coverage
304			and solvent molecule numbers at different temperatures using a
305			range of energy fluxes.}
306
307			\begin{tabular}{cccccc}
308			\hline\hline
309			Thiol & $\langle T\rangle$ & & $J_z$ & $G$ & $G^\prime$ \\
310			coverage (\%) & (K) & $N_{hexane}$ & (GW/m$^2$) &
311			\multicolumn{2}{c}{(MW/m$^2$/K)} \\
312			\hline
313			0.0 & 200 & 200 & 0.96 & 43.3 & 42.7 \\
314			& & & 1.91 & 45.7 & 42.9 \\
315			& & 166 & 0.96 & 43.1 & 53.4 \\
316			88.9 & 200 & 166 & 1.94 & 172 & 108 \\
317			100.0 & 250 & 200 & 0.96 & 81.8 & 67.0 \\
318			& & 166 & 0.98 & 79.0 & 62.9 \\
319			& & & 1.44 & 76.2 & 64.8 \\
320			& 200 & 200 & 1.92 & 129 & 87.3 \\
321			& & & 1.93 & 131 & 77.5 \\
322			& & 166 & 0.97 & 115 & 69.3 \\
323			& & & 1.94 & 125 & 87.1 \\
324			\hline\hline
325			\end{tabular}
326			\label{AuThiolHexaneUA}
327			\end{center}
328			\end{minipage}
329			\end{table*}
330
331			For the all-atom model, the liquid hexane phase was not stable under NPT
332			conditions. Therefore, the simulation length scale parameters are
333			adopted from previous equilibration results of the united-atom model
334			at 200K. Table \ref{AuThiolHexaneAA} shows the results of these
335			simulations. The conductivity values calculated with full capping
336			agent coverage are substantially larger than observed in the
337			united-atom model, and is even higher than predicted by
338			experiments. It is possible that our parameters for metal-non-metal
339			particle interactions lead to an overestimate of the interfacial
340			thermal conductivity, although the active C-H vibrations in the
341			all-atom model (which should not be appreciably populated at normal
342			temperatures) could also account for this high conductivity. The major
343			thermal transfer barrier of Au/butanethiol/hexane interface is between
344			the liquid phase and the capping agent, so extra degrees of freedom
345			such as the C-H vibrations could enhance heat exchange between these
346			two phases and result in a much higher conductivity.
347
348			\begin{table*}
349			\begin{minipage}{\linewidth}
350			\begin{center}
351
352			\caption{Computed interfacial thermal conductivity ($G$ and
353			$G^\prime$) values for the Au/butanethiol/hexane interface
354			with all-atom model and different capping agent coverage at
355			200K using a range of energy fluxes.}
356
357			\begin{tabular}{cccc}
358			\hline\hline
359			Thiol & $J_z$ & $G$ & $G^\prime$ \\
360			coverage (\%) & (GW/m$^2$) & \multicolumn{2}{c}{(MW/m$^2$/K)} \\
361			\hline
362			0.0 & 0.95 & 28.5 & 27.2 \\
363			& 1.88 & 30.3 & 28.9 \\
364			100.0 & 2.87 & 551 & 294 \\
365			& 3.81 & 494 & 193 \\
366			\hline\hline
367			\end{tabular}
368			\label{AuThiolHexaneAA}
369			\end{center}
370			\end{minipage}
371			\end{table*}
372
373			%subsubsection{Vibrational spectrum study on conductance mechanism}
374			To investigate the mechanism of this interfacial thermal conductance,
375			the vibrational spectra of various gold systems were obtained and are
376			shown as in the upper panel of Fig. \ref{vibration}. To obtain these
377			spectra, one first runs a simulation in the NVE ensemble and collects
378			snapshots of configurations; these configurations are used to compute
379			the velocity auto-correlation functions, which is used to construct a
380			power spectrum via a Fourier transform. The gold surfaces covered by
381			butanethiol molecules exhibit an additional peak observed at a
382			frequency of $\sim$170cm$^{-1}$, which is attributed to the vibration
383			of the S-Au bond. This vibration enables efficient thermal transport
384			from surface Au atoms to the capping agents. Simultaneously, as shown
385			in the lower panel of Fig. \ref{vibration}, the large overlap of the
386			vibration spectra of butanethiol and hexane in the all-atom model,
387			including the C-H vibration, also suggests high thermal exchange
388			efficiency. The combination of these two effects produces the drastic
389			interfacial thermal conductance enhancement in the all-atom model.
390
391			\begin{figure}
392			\includegraphics[width=\linewidth]{vibration}
393			\caption{Vibrational spectra obtained for gold in different
394			environments (upper panel) and for Au/thiol/hexane simulation in
395			all-atom model (lower panel).}
396			\label{vibration}
397			\end{figure}
398			% 600dpi, letter size. too large?
399
400
401	gezelter	3717	\section{Acknowledgments}
402			Support for this project was provided by the National Science
403			Foundation under grant CHE-0848243. Computational time was provided by
404			the Center for Research Computing (CRC) at the University of Notre
405			Dame. \newpage
406
407			\bibliography{interfacial}
408
409			\end{doublespace}
410			\end{document}
411