ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/interfacial/interfacial.tex
(Generate patch)

Comparing interfacial/interfacial.tex (file contents):
Revision 3745 by skuang, Thu Jul 21 00:04:26 2011 UTC vs.
Revision 3760 by skuang, Mon Aug 1 17:57:34 2011 UTC

# Line 73 | Line 73 | Interfacial thermal conductance is extensively studied
73   %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
74  
75   \section{Introduction}
76 < Interfacial thermal conductance is extensively studied both
77 < experimentally and computationally\cite{cahill:793}, due to its
78 < importance in nanoscale science and technology. Reliability of
79 < nanoscale devices depends on their thermal transport
80 < properties. Unlike bulk homogeneous materials, nanoscale materials
81 < features significant presence of interfaces, and these interfaces
82 < could dominate the heat transfer behavior of these
83 < materials. Furthermore, these materials are generally heterogeneous,
84 < which challenges traditional research methods for homogeneous
85 < systems.
76 > Due to the importance of heat flow in nanotechnology, interfacial
77 > thermal conductance has been studied extensively both experimentally
78 > and computationally.\cite{cahill:793} Unlike bulk materials, nanoscale
79 > materials have a significant fraction of their atoms at interfaces,
80 > and the chemical details of these interfaces govern the heat transfer
81 > behavior. Furthermore, the interfaces are
82 > heterogeneous (e.g. solid - liquid), which provides a challenge to
83 > traditional methods developed for homogeneous systems.
84  
85 < Heat conductance of molecular and nano-scale interfaces will be
86 < affected by the chemical details of the surface. Experimentally,
87 < various interfaces have been investigated for their thermal
88 < conductance properties. Wang {\it et al.} studied heat transport
89 < through long-chain hydrocarbon monolayers on gold substrate at
90 < individual molecular level\cite{Wang10082007}; Schmidt {\it et al.}
91 < studied the role of CTAB on thermal transport between gold nanorods
92 < and solvent\cite{doi:10.1021/jp8051888}; Juv\'e {\it et al.} studied
85 > Experimentally, various interfaces have been investigated for their
86 > thermal conductance. Cahill and coworkers studied nanoscale thermal
87 > transport from metal nanoparticle/fluid interfaces, to epitaxial
88 > TiN/single crystal oxides interfaces, to hydrophilic and hydrophobic
89 > interfaces between water and solids with different self-assembled
90 > monolayers.\cite{Wilson:2002uq,PhysRevB.67.054302,doi:10.1021/jp048375k,PhysRevLett.96.186101}
91 > Wang {\it et al.} studied heat transport through
92 > long-chain hydrocarbon monolayers on gold substrate at individual
93 > molecular level,\cite{Wang10082007} Schmidt {\it et al.} studied the
94 > role of CTAB on thermal transport between gold nanorods and
95 > solvent,\cite{doi:10.1021/jp8051888} and Juv\'e {\it et al.} studied
96   the cooling dynamics, which is controlled by thermal interface
97   resistence of glass-embedded metal
98 < nanoparticles\cite{PhysRevB.80.195406}. Although interfaces are
99 < commonly barriers for heat transport, Alper {\it et al.} suggested
100 < that specific ligands (capping agents) could completely eliminate this
101 < barrier ($G\rightarrow\infty$)\cite{doi:10.1021/la904855s}.
98 > nanoparticles.\cite{PhysRevB.80.195406} Although interfaces are
99 > normally considered barriers for heat transport, Alper {\it et al.}
100 > suggested that specific ligands (capping agents) could completely
101 > eliminate this barrier
102 > ($G\rightarrow\infty$).\cite{doi:10.1021/la904855s}
103  
104   Theoretical and computational models have also been used to study the
105   interfacial thermal transport in order to gain an understanding of
# Line 105 | Line 107 | atoms)\cite{hase:2010,hase:2011}. However, ensemble av
107   employed Non-Equilibrium Molecular Dynamics (NEMD) simulations to
108   study thermal transport from hot Au(111) substrate to a self-assembled
109   monolayer of alkylthiol with relatively long chain (8-20 carbon
110 < atoms)\cite{hase:2010,hase:2011}. However, ensemble averaged
110 > atoms).\cite{hase:2010,hase:2011} However, ensemble averaged
111   measurements for heat conductance of interfaces between the capping
112 < monolayer on Au and a solvent phase has yet to be studied.
113 < The comparatively low thermal flux through interfaces is
114 < difficult to measure with Equilibrium MD or forward NEMD simulation
115 < methods. Therefore, the Reverse NEMD (RNEMD) methods would have the
116 < advantage of having this difficult to measure flux known when studying
117 < the thermal transport across interfaces, given that the simulation
118 < methods being able to effectively apply an unphysical flux in
119 < non-homogeneous systems.
112 > monolayer on Au and a solvent phase have yet to be studied with their
113 > approach. The comparatively low thermal flux through interfaces is
114 > difficult to measure with Equilibrium
115 > MD\cite{doi:10.1080/0026897031000068578} or forward NEMD simulation
116 > methods. Therefore, the Reverse NEMD (RNEMD)
117 > methods\cite{MullerPlathe:1997xw,kuang:164101} would have the
118 > advantage of applying this difficult to measure flux (while measuring
119 > the resulting gradient), given that the simulation methods being able
120 > to effectively apply an unphysical flux in non-homogeneous systems.
121 > Garde and coworkers\cite{garde:nl2005,garde:PhysRevLett2009} applied
122 > this approach to various liquid interfaces and studied how thermal
123 > conductance (or resistance) is dependent on chemistry details of
124 > interfaces, e.g. hydrophobic and hydrophilic aqueous interfaces.
125  
126 < Recently, we have developed the Non-Isotropic Velocity Scaling (NIVS)
126 > Recently, we have developed a Non-Isotropic Velocity Scaling (NIVS)
127   algorithm for RNEMD simulations\cite{kuang:164101}. This algorithm
128   retains the desirable features of RNEMD (conservation of linear
129   momentum and total energy, compatibility with periodic boundary
# Line 131 | Line 138 | underlying mechanism for this phenomena was investigat
138   properties. Different models were used for both the capping agent and
139   the solvent force field parameters. Using the NIVS algorithm, the
140   thermal transport across these interfaces was studied and the
141 < underlying mechanism for this phenomena was investigated.
141 > underlying mechanism for the phenomena was investigated.
142  
136 [MAY ADD WHY STUDY AU-THIOL SURFACE; CITE SHAOYI JIANG]
137
143   \section{Methodology}
144   \subsection{Imposd-Flux Methods in MD Simulations}
145 < For systems with low interfacial conductivity one must have a method
146 < capable of generating relatively small fluxes, compared to those
147 < required for bulk conductivity. This requirement makes the calculation
148 < even more difficult for those slowly-converging equilibrium
149 < methods\cite{Viscardy:2007lq}.
150 < Forward methods impose gradient, but in interfacail conditions it is
151 < not clear what behavior to impose at the boundary...
152 < Imposed-flux reverse non-equilibrium
153 < methods\cite{MullerPlathe:1997xw} have the flux set {\it a priori} and
154 < the thermal response becomes easier to
155 < measure than the flux. Although M\"{u}ller-Plathe's original momentum
156 < swapping approach can be used for exchanging energy between particles
157 < of different identity, the kinetic energy transfer efficiency is
158 < affected by the mass difference between the particles, which limits
159 < its application on heterogeneous interfacial systems.
145 > Steady state MD simulations have an advantage in that not many
146 > trajectories are needed to study the relationship between thermal flux
147 > and thermal gradients. For systems with low interfacial conductance,
148 > one must have a method capable of generating or measuring relatively
149 > small fluxes, compared to those required for bulk conductivity. This
150 > requirement makes the calculation even more difficult for
151 > slowly-converging equilibrium methods.\cite{Viscardy:2007lq} Forward
152 > NEMD methods impose a gradient (and measure a flux), but at interfaces
153 > it is not clear what behavior should be imposed at the boundaries
154 > between materials.  Imposed-flux reverse non-equilibrium
155 > methods\cite{MullerPlathe:1997xw} set the flux {\it a priori} and
156 > the thermal response becomes an easy-to-measure quantity.  Although
157 > M\"{u}ller-Plathe's original momentum swapping approach can be used
158 > for exchanging energy between particles of different identity, the
159 > kinetic energy transfer efficiency is affected by the mass difference
160 > between the particles, which limits its application on heterogeneous
161 > interfacial systems.
162  
163 < The non-isotropic velocity scaling (NIVS)\cite{kuang:164101} approach to
164 < non-equilibrium MD simulations is able to impose a wide range of
163 > The non-isotropic velocity scaling (NIVS) \cite{kuang:164101} approach
164 > to non-equilibrium MD simulations is able to impose a wide range of
165   kinetic energy fluxes without obvious perturbation to the velocity
166   distributions of the simulated systems. Furthermore, this approach has
167   the advantage in heterogeneous interfaces in that kinetic energy flux
# Line 171 | Line 178 | momenta and energy and does not depend on an external
178   for computing thermal conductivities. The NIVS algorithm conserves
179   momenta and energy and does not depend on an external thermostat.
180  
181 < \subsection{Defining Interfacial Thermal Conductivity $G$}
182 < For interfaces with a relatively low interfacial conductance, the bulk
183 < regions on either side of an interface rapidly come to a state in
184 < which the two phases have relatively homogeneous (but distinct)
185 < temperatures. The interfacial thermal conductivity $G$ can therefore
186 < be approximated as:
181 > \subsection{Defining Interfacial Thermal Conductivity ($G$)}
182 >
183 > For an interface with relatively low interfacial conductance, and a
184 > thermal flux between two distinct bulk regions, the regions on either
185 > side of the interface rapidly come to a state in which the two phases
186 > have relatively homogeneous (but distinct) temperatures. The
187 > interfacial thermal conductivity $G$ can therefore be approximated as:
188   \begin{equation}
189 < G = \frac{E_{total}}{2 t L_x L_y \left( \langle T_\mathrm{hot}\rangle -
189 >  G = \frac{E_{total}}{2 t L_x L_y \left( \langle T_\mathrm{hot}\rangle -
190      \langle T_\mathrm{cold}\rangle \right)}
191   \label{lowG}
192   \end{equation}
193 < where ${E_{total}}$ is the imposed non-physical kinetic energy
194 < transfer and ${\langle T_\mathrm{hot}\rangle}$ and ${\langle
195 <  T_\mathrm{cold}\rangle}$ are the average observed temperature of the
196 < two separated phases.
193 > where ${E_{total}}$ is the total imposed non-physical kinetic energy
194 > transfer during the simulation and ${\langle T_\mathrm{hot}\rangle}$
195 > and ${\langle T_\mathrm{cold}\rangle}$ are the average observed
196 > temperature of the two separated phases.  For an applied flux $J_z$
197 > operating over a simulation time $t$ on a periodically-replicated slab
198 > of dimensions $L_x \times L_y$, $E_{total} = J_z *(t)*(2 L_x L_y)$.
199  
200   When the interfacial conductance is {\it not} small, there are two
201 < ways to define $G$.
202 <
203 < One way is to assume the temperature is discrete on the two sides of
204 < the interface. $G$ can be calculated using the applied thermal flux
205 < $J$ and the maximum temperature difference measured along the thermal
206 < gradient max($\Delta T$), which occurs at the Gibbs deviding surface
207 < (Figure \ref{demoPic}):
201 > ways to define $G$. One common way is to assume the temperature is
202 > discrete on the two sides of the interface. $G$ can be calculated
203 > using the applied thermal flux $J$ and the maximum temperature
204 > difference measured along the thermal gradient max($\Delta T$), which
205 > occurs at the Gibbs deviding surface (Figure \ref{demoPic}). This is
206 > known as the Kapitza conductance, which is the inverse of the Kapitza
207 > resistance.
208   \begin{equation}
209 < G=\frac{J}{\Delta T}
209 >  G=\frac{J}{\Delta T}
210   \label{discreteG}
211   \end{equation}
212  
# Line 216 | Line 226 | the magnitude of thermal conductivity $\lambda$ change
226  
227   The other approach is to assume a continuous temperature profile along
228   the thermal gradient axis (e.g. $z$) and define $G$ at the point where
229 < the magnitude of thermal conductivity $\lambda$ change reach its
229 > the magnitude of thermal conductivity ($\lambda$) change reaches its
230   maximum, given that $\lambda$ is well-defined throughout the space:
231   \begin{equation}
232   G^\prime = \Big|\frac{\partial\lambda}{\partial z}\Big|
# Line 227 | Line 237 | With the temperature profile obtained from simulations
237   \label{derivativeG}
238   \end{equation}
239  
240 < With the temperature profile obtained from simulations, one is able to
240 > With temperature profiles obtained from simulation, one is able to
241   approximate the first and second derivatives of $T$ with finite
242 < difference methods and thus calculate $G^\prime$.
242 > difference methods and calculate $G^\prime$. In what follows, both
243 > definitions have been used, and are compared in the results.
244  
245 < In what follows, both definitions have been used for calculation and
246 < are compared in the results.
247 <
248 < To compare the above definitions ($G$ and $G^\prime$), we have modeled
249 < a metal slab with its (111) surfaces perpendicular to the $z$-axis of
250 < our simulation cells. Both with and without capping agents on the
240 < surfaces, the metal slab is solvated with simple organic solvents, as
241 < illustrated in Figure \ref{demoPic}.
245 > To investigate the interfacial conductivity at metal / solvent
246 > interfaces, we have modeled a metal slab with its (111) surfaces
247 > perpendicular to the $z$-axis of our simulation cells. The metal slab
248 > has been prepared both with and without capping agents on the exposed
249 > surface, and has been solvated with simple organic solvents, as
250 > illustrated in Figure \ref{gradT}.
251  
252   With the simulation cell described above, we are able to equilibrate
253   the system and impose an unphysical thermal flux between the liquid
254   and the metal phase using the NIVS algorithm. By periodically applying
255 < the unphysical flux, we are able to obtain a temperature profile and
256 < its spatial derivatives. These quantities enable the evaluation of the
257 < interfacial thermal conductance of a surface. Figure \ref{gradT} is an
258 < example how those applied thermal fluxes can be used to obtain the 1st
250 < and 2nd derivatives of the temperature profile.
255 > the unphysical flux, we obtained a temperature profile and its spatial
256 > derivatives. Figure \ref{gradT} shows how an applied thermal flux can
257 > be used to obtain the 1st and 2nd derivatives of the temperature
258 > profile.
259  
260   \begin{figure}
261   \includegraphics[width=\linewidth]{gradT}
# Line 261 | Line 269 | OpenMD\cite{Meineke:2005gd,openmd}, and was used for o
269   \section{Computational Details}
270   \subsection{Simulation Protocol}
271   The NIVS algorithm has been implemented in our MD simulation code,
272 < OpenMD\cite{Meineke:2005gd,openmd}, and was used for our
273 < simulations. Different slab thickness (layer numbers of Au) were
274 < simulated. Metal slabs were first equilibrated under atmospheric
275 < pressure (1 atm) and a desired temperature (e.g. 200K). After
276 < equilibration, butanethiol capping agents were placed at three-fold
277 < sites on the Au(111) surfaces. The maximum butanethiol capacity on Au
278 < surface is $1/3$ of the total number of surface Au
279 < atoms\cite{vlugt:cpc2007154}. A series of different coverages was
280 < investigated in order to study the relation between coverage and
281 < interfacial conductance.
272 > OpenMD\cite{Meineke:2005gd,openmd}, and was used for our simulations.
273 > Metal slabs of 6 or 11 layers of Au atoms were first equilibrated
274 > under atmospheric pressure (1 atm) and 200K. After equilibration,
275 > butanethiol capping agents were placed at three-fold hollow sites on
276 > the Au(111) surfaces. These sites are either {\it fcc} or {\it
277 >  hcp} sites, although Hase {\it et al.} found that they are
278 > equivalent in a heat transfer process,\cite{hase:2010} so we did not
279 > distinguish between these sites in our study. The maximum butanethiol
280 > capacity on Au surface is $1/3$ of the total number of surface Au
281 > atoms, and the packing forms a $(\sqrt{3}\times\sqrt{3})R30^\circ$
282 > structure\cite{doi:10.1021/ja00008a001,doi:10.1021/cr9801317}. A
283 > series of lower coverages was also prepared by eliminating
284 > butanethiols from the higher coverage surface in a regular manner. The
285 > lower coverages were prepared in order to study the relation between
286 > coverage and interfacial conductance.
287  
288   The capping agent molecules were allowed to migrate during the
289   simulations. They distributed themselves uniformly and sampled a
290   number of three-fold sites throughout out study. Therefore, the
291 < initial configuration would not noticeably affect the sampling of a
291 > initial configuration does not noticeably affect the sampling of a
292   variety of configurations of the same coverage, and the final
293   conductance measurement would be an average effect of these
294 < configurations explored in the simulations. [MAY NEED FIGURES]
294 > configurations explored in the simulations.
295  
296 < After the modified Au-butanethiol surface systems were equilibrated
297 < under canonical ensemble, organic solvent molecules were packed in the
298 < previously empty part of the simulation cells\cite{packmol}. Two
296 > After the modified Au-butanethiol surface systems were equilibrated in
297 > the canonical (NVT) ensemble, organic solvent molecules were packed in
298 > the previously empty part of the simulation cells.\cite{packmol} Two
299   solvents were investigated, one which has little vibrational overlap
300 < with the alkanethiol and a planar shape (toluene), and one which has
301 < similar vibrational frequencies and chain-like shape ({\it n}-hexane).
300 > with the alkanethiol and which has a planar shape (toluene), and one
301 > which has similar vibrational frequencies to the capping agent and
302 > chain-like shape ({\it n}-hexane).
303  
304 < The space filled by solvent molecules, i.e. the gap between
305 < periodically repeated Au-butanethiol surfaces should be carefully
306 < chosen. A very long length scale for the thermal gradient axis ($z$)
293 < may cause excessively hot or cold temperatures in the middle of the
304 > The simulation cells were not particularly extensive along the
305 > $z$-axis, as a very long length scale for the thermal gradient may
306 > cause excessively hot or cold temperatures in the middle of the
307   solvent region and lead to undesired phenomena such as solvent boiling
308   or freezing when a thermal flux is applied. Conversely, too few
309   solvent molecules would change the normal behavior of the liquid
310   phase. Therefore, our $N_{solvent}$ values were chosen to ensure that
311 < these extreme cases did not happen to our simulations. And the
312 < corresponding spacing is usually $35 \sim 75$\AA.
311 > these extreme cases did not happen to our simulations. The spacing
312 > between periodic images of the gold interfaces is $45 \sim 75$\AA in
313 > our simulations.
314  
315 < The initial configurations generated by Packmol are further
316 < equilibrated with the $x$ and $y$ dimensions fixed, only allowing
317 < length scale change in $z$ dimension. This is to ensure that the
318 < equilibration of liquid phase does not affect the metal crystal
319 < structure in $x$ and $y$ dimensions. Further equilibration are run
320 < under NVT and then NVE ensembles.
315 > The initial configurations generated are further equilibrated with the
316 > $x$ and $y$ dimensions fixed, only allowing the $z$-length scale to
317 > change. This is to ensure that the equilibration of liquid phase does
318 > not affect the metal's crystalline structure. Comparisons were made
319 > with simulations that allowed changes of $L_x$ and $L_y$ during NPT
320 > equilibration. No substantial changes in the box geometry were noticed
321 > in these simulations. After ensuring the liquid phase reaches
322 > equilibrium at atmospheric pressure (1 atm), further equilibration was
323 > carried out under canonical (NVT) and microcanonical (NVE) ensembles.
324  
325 < After the systems reach equilibrium, NIVS is implemented to impose a
326 < periodic unphysical thermal flux between the metal and the liquid
327 < phase. Most of our simulations are under an average temperature of
328 < $\sim$200K. Therefore, this flux usually comes from the metal to the
325 > After the systems reach equilibrium, NIVS was used to impose an
326 > unphysical thermal flux between the metal and the liquid phases. Most
327 > of our simulations were done under an average temperature of
328 > $\sim$200K. Therefore, thermal flux usually came from the metal to the
329   liquid so that the liquid has a higher temperature and would not
330 < freeze due to excessively low temperature. This induced temperature
331 < gradient is stablized and the simulation cell is devided evenly into
332 < N slabs along the $z$-axis and the temperatures of each slab are
333 < recorded. When the slab width $d$ of each slab is the same, the
334 < derivatives of $T$ with respect to slab number $n$ can be directly
335 < used for $G^\prime$ calculations:
336 < \begin{equation}
337 < G^\prime = |J_z|\Big|\frac{\partial^2 T}{\partial z^2}\Big|
330 > freeze due to lowered temperatures. After this induced temperature
331 > gradient had stablized, the temperature profile of the simulation cell
332 > was recorded. To do this, the simulation cell is devided evenly into
333 > $N$ slabs along the $z$-axis. The average temperatures of each slab
334 > are recorded for 1$\sim$2 ns. When the slab width $d$ of each slab is
335 > the same, the derivatives of $T$ with respect to slab number $n$ can
336 > be directly used for $G^\prime$ calculations: \begin{equation}
337 >  G^\prime = |J_z|\Big|\frac{\partial^2 T}{\partial z^2}\Big|
338           \Big/\left(\frac{\partial T}{\partial z}\right)^2
339           = |J_z|\Big|\frac{1}{d^2}\frac{\partial^2 T}{\partial n^2}\Big|
340           \Big/\left(\frac{1}{d}\frac{\partial T}{\partial n}\right)^2
# Line 326 | Line 343 | G^\prime = |J_z|\Big|\frac{\partial^2 T}{\partial z^2}
343   \label{derivativeG2}
344   \end{equation}
345  
346 + All of the above simulation procedures use a time step of 1 fs. Each
347 + equilibration stage took a minimum of 100 ps, although in some cases,
348 + longer equilibration stages were utilized.
349 +
350   \subsection{Force Field Parameters}
351 < Our simulations include various components. Figure \ref{demoMol}
352 < demonstrates the sites defined for both United-Atom and All-Atom
353 < models of the organic solvent and capping agent molecules in our
354 < simulations. Force field parameter descriptions are needed for
351 > Our simulations include a number of chemically distinct components.
352 > Figure \ref{demoMol} demonstrates the sites defined for both
353 > United-Atom and All-Atom models of the organic solvent and capping
354 > agents in our simulations. Force field parameters are needed for
355   interactions both between the same type of particles and between
356   particles of different species.
357  
# Line 340 | Line 361 | particles of different species.
361    these simulations. The chemically-distinct sites (a-e) are expanded
362    in terms of constituent atoms for both United Atom (UA) and All Atom
363    (AA) force fields.  Most parameters are from
364 <  Refs. \protect\cite{TraPPE-UA.alkanes,TraPPE-UA.alkylbenzenes} (UA) and
365 <  \protect\cite{OPLSAA} (AA).  Cross-interactions with the Au atoms are given
366 <  in Table \ref{MnM}.}
364 >  Refs. \protect\cite{TraPPE-UA.alkanes,TraPPE-UA.alkylbenzenes,TraPPE-UA.thiols}
365 >  (UA) and \protect\cite{OPLSAA} (AA). Cross-interactions with the Au
366 >  atoms are given in Table \ref{MnM}.}
367   \label{demoMol}
368   \end{figure}
369  
370   The Au-Au interactions in metal lattice slab is described by the
371 < quantum Sutton-Chen (QSC) formulation\cite{PhysRevB.59.3527}. The QSC
371 > quantum Sutton-Chen (QSC) formulation.\cite{PhysRevB.59.3527} The QSC
372   potentials include zero-point quantum corrections and are
373   reparametrized for accurate surface energies compared to the
374 < Sutton-Chen potentials\cite{Chen90}.
374 > Sutton-Chen potentials.\cite{Chen90}
375  
376 < For both solvent molecules, straight chain {\it n}-hexane and aromatic
377 < toluene, United-Atom (UA) and All-Atom (AA) models are used
378 < respectively. The TraPPE-UA
376 > For the two solvent molecules, {\it n}-hexane and toluene, two
377 > different atomistic models were utilized. Both solvents were modeled
378 > using United-Atom (UA) and All-Atom (AA) models. The TraPPE-UA
379   parameters\cite{TraPPE-UA.alkanes,TraPPE-UA.alkylbenzenes} are used
380   for our UA solvent molecules. In these models, sites are located at
381   the carbon centers for alkyl groups. Bonding interactions, including
382   bond stretches and bends and torsions, were used for intra-molecular
383 < sites not separated by more than 3 bonds. Otherwise, for non-bonded
384 < interactions, Lennard-Jones potentials are used. [MORE CITATION?]
383 > sites closer than 3 bonds. For non-bonded interactions, Lennard-Jones
384 > potentials are used.
385  
386 < By eliminating explicit hydrogen atoms, these models are simple and
387 < computationally efficient, while maintains good accuracy. However, the
388 < TraPPE-UA for alkanes is known to predict a lower boiling point than
389 < experimental values. Considering that after an unphysical thermal flux
390 < is applied to a system, the temperature of ``hot'' area in the liquid
391 < phase would be significantly higher than the average, to prevent over
392 < heating and boiling of the liquid phase, the average temperature in
393 < our simulations should be much lower than the liquid boiling point.
386 > By eliminating explicit hydrogen atoms, the TraPPE-UA models are
387 > simple and computationally efficient, while maintaining good accuracy.
388 > However, the TraPPE-UA model for alkanes is known to predict a slighly
389 > lower boiling point than experimental values. This is one of the
390 > reasons we used a lower average temperature (200K) for our
391 > simulations. If heat is transferred to the liquid phase during the
392 > NIVS simulation, the liquid in the hot slab can actually be
393 > substantially warmer than the mean temperature in the simulation. The
394 > lower mean temperatures therefore prevent solvent boiling.
395  
396 < For UA-toluene model, the non-bonded potentials between
397 < inter-molecular sites have a similar Lennard-Jones formulation. For
398 < intra-molecular interactions, considering the stiffness of the benzene
399 < ring, rigid body constraints are applied for further computational
400 < efficiency. All bonds in the benzene ring and between the ring and the
379 < methyl group remain rigid during the progress of simulations.
396 > For UA-toluene, the non-bonded potentials between intermolecular sites
397 > have a similar Lennard-Jones formulation. The toluene molecules were
398 > treated as a single rigid body, so there was no need for
399 > intramolecular interactions (including bonds, bends, or torsions) in
400 > this solvent model.
401  
402   Besides the TraPPE-UA models, AA models for both organic solvents are
403 < included in our studies as well. For hexane, the OPLS-AA\cite{OPLSAA}
404 < force field is used. Additional explicit hydrogen sites were
403 > included in our studies as well. The OPLS-AA\cite{OPLSAA} force fields
404 > were used. For hexane, additional explicit hydrogen sites were
405   included. Besides bonding and non-bonded site-site interactions,
406   partial charges and the electrostatic interactions were added to each
407 < CT and HC site. For toluene, the United Force Field developed by
408 < Rapp\'{e} {\it et al.}\cite{doi:10.1021/ja00051a040} is
409 < adopted. Without the rigid body constraints, bonding interactions were
389 < included. For the aromatic ring, improper torsions (inversions) were
390 < added as an extra potential for maintaining the planar shape.
391 < [MORE CITATION?]
407 > CT and HC site. For toluene, a flexible model for the toluene molecule
408 > was utilized which included bond, bend, torsion, and inversion
409 > potentials to enforce ring planarity.
410  
411 < The capping agent in our simulations, the butanethiol molecules can
412 < either use UA or AA model. The TraPPE-UA force fields includes
411 > The butanethiol capping agent in our simulations, were also modeled
412 > with both UA and AA model. The TraPPE-UA force field includes
413   parameters for thiol molecules\cite{TraPPE-UA.thiols} and are used for
414   UA butanethiol model in our simulations. The OPLS-AA also provides
415   parameters for alkyl thiols. However, alkyl thiols adsorbed on Au(111)
416 < surfaces do not have the hydrogen atom bonded to sulfur. To adapt this
417 < change and derive suitable parameters for butanethiol adsorbed on
418 < Au(111) surfaces, we adopt the S parameters from Luedtke and
419 < Landman\cite{landman:1998} and modify parameters for its neighbor C
420 < atom for charge balance in the molecule. Note that the model choice
421 < (UA or AA) of capping agent can be different from the
422 < solvent. Regardless of model choice, the force field parameters for
423 < interactions between capping agent and solvent can be derived using
406 < Lorentz-Berthelot Mixing Rule:
416 > surfaces do not have the hydrogen atom bonded to sulfur. To derive
417 > suitable parameters for butanethiol adsorbed on Au(111) surfaces, we
418 > adopt the S parameters from Luedtke and Landman\cite{landman:1998} and
419 > modify the parameters for the CTS atom to maintain charge neutrality
420 > in the molecule.  Note that the model choice (UA or AA) for the capping
421 > agent can be different from the solvent. Regardless of model choice,
422 > the force field parameters for interactions between capping agent and
423 > solvent can be derived using Lorentz-Berthelot Mixing Rule:
424   \begin{eqnarray}
425 < \sigma_{ij} & = & \frac{1}{2} \left(\sigma_{ii} + \sigma_{jj}\right) \\
426 < \epsilon_{ij} & = & \sqrt{\epsilon_{ii}\epsilon_{jj}}
425 >  \sigma_{ij} & = & \frac{1}{2} \left(\sigma_{ii} + \sigma_{jj}\right) \\
426 >  \epsilon_{ij} & = & \sqrt{\epsilon_{ii}\epsilon_{jj}}
427   \end{eqnarray}
428  
429 < To describe the interactions between metal Au and non-metal capping
430 < agent and solvent particles, we refer to an adsorption study of alkyl
431 < thiols on gold surfaces by Vlugt {\it et
432 <  al.}\cite{vlugt:cpc2007154} They fitted an effective Lennard-Jones
433 < form of potential parameters for the interaction between Au and
434 < pseudo-atoms CH$_x$ and S based on a well-established and widely-used
435 < effective potential of Hautman and Klein\cite{hautman:4994} for the
436 < Au(111) surface. As our simulations require the gold lattice slab to
437 < be non-rigid so that it could accommodate kinetic energy for thermal
421 < transport study purpose, the pair-wise form of potentials is
422 < preferred.
429 > To describe the interactions between metal (Au) and non-metal atoms,
430 > we refer to an adsorption study of alkyl thiols on gold surfaces by
431 > Vlugt {\it et al.}\cite{vlugt:cpc2007154} They fitted an effective
432 > Lennard-Jones form of potential parameters for the interaction between
433 > Au and pseudo-atoms CH$_x$ and S based on a well-established and
434 > widely-used effective potential of Hautman and Klein for the Au(111)
435 > surface.\cite{hautman:4994} As our simulations require the gold slab
436 > to be flexible to accommodate thermal excitation, the pair-wise form
437 > of potentials they developed was used for our study.
438  
439 < Besides, the potentials developed from {\it ab initio} calculations by
440 < Leng {\it et al.}\cite{doi:10.1021/jp034405s} are adopted for the
441 < interactions between Au and aromatic C/H atoms in toluene. A set of
442 < pseudo Lennard-Jones parameters were provided for Au in their force
443 < fields. By using the Mixing Rule, this can be used to derive pair-wise
444 < potentials for non-bonded interactions between Au and non-metal sites.
439 > The potentials developed from {\it ab initio} calculations by Leng
440 > {\it et al.}\cite{doi:10.1021/jp034405s} are adopted for the
441 > interactions between Au and aromatic C/H atoms in toluene. However,
442 > the Lennard-Jones parameters between Au and other types of particles,
443 > (e.g. AA alkanes) have not yet been established. For these
444 > interactions, the Lorentz-Berthelot mixing rule can be used to derive
445 > effective single-atom LJ parameters for the metal using the fit values
446 > for toluene. These are then used to construct reasonable mixing
447 > parameters for the interactions between the gold and other atoms.
448 > Table \ref{MnM} summarizes the ``metal/non-metal'' parameters used in
449 > our simulations.
450  
431 However, the Lennard-Jones parameters between Au and other types of
432 particles, such as All-Atom normal alkanes in our simulations are not
433 yet well-established. For these interactions, we attempt to derive
434 their parameters using the Mixing Rule. To do this, Au pseudo
435 Lennard-Jones parameters for ``Metal-non-Metal'' (MnM) interactions
436 were first extracted from the Au-CH$_x$ parameters by applying the
437 Mixing Rule reversely. Table \ref{MnM} summarizes these ``MnM''
438 parameters in our simulations.
439
451   \begin{table*}
452    \begin{minipage}{\linewidth}
453      \begin{center}
# Line 473 | Line 484 | parameters in our simulations.
484   \end{table*}
485  
486  
487 < \section{Results and Discussions}
488 < [MAY HAVE A BRIEF SUMMARY]
489 < \subsection{How Simulation Parameters Affects $G$}
490 < [MAY NOT PUT AT FIRST]
491 < We have varied our protocol or other parameters of the simulations in
492 < order to investigate how these factors would affect the measurement of
493 < $G$'s. It turned out that while some of these parameters would not
494 < affect the results substantially, some other changes to the
484 < simulations would have a significant impact on the measurement
485 < results.
487 > \section{Results}
488 > There are many factors contributing to the measured interfacial
489 > conductance; some of these factors are physically motivated
490 > (e.g. coverage of the surface by the capping agent coverage and
491 > solvent identity), while some are governed by parameters of the
492 > methodology (e.g. applied flux and the formulas used to obtain the
493 > conductance). In this section we discuss the major physical and
494 > calculational effects on the computed conductivity.
495  
496 < In some of our simulations, we allowed $L_x$ and $L_y$ to change
488 < during equilibrating the liquid phase. Due to the stiffness of the
489 < crystalline Au structure, $L_x$ and $L_y$ would not change noticeably
490 < after equilibration. Although $L_z$ could fluctuate $\sim$1\% after a
491 < system is fully equilibrated in the NPT ensemble, this fluctuation, as
492 < well as those of $L_x$ and $L_y$ (which is significantly smaller),
493 < would not be magnified on the calculated $G$'s, as shown in Table
494 < \ref{AuThiolHexaneUA}. This insensivity to $L_i$ fluctuations allows
495 < reliable measurement of $G$'s without the necessity of extremely
496 < cautious equilibration process.
496 > \subsection{Effects due to capping agent coverage}
497  
498 < As stated in our computational details, the spacing filled with
499 < solvent molecules can be chosen within a range. This allows some
500 < change of solvent molecule numbers for the same Au-butanethiol
501 < surfaces. We did this study on our Au-butanethiol/hexane
502 < simulations. Nevertheless, the results obtained from systems of
503 < different $N_{hexane}$ did not indicate that the measurement of $G$ is
504 < susceptible to this parameter. For computational efficiency concern,
505 < smaller system size would be preferable, given that the liquid phase
506 < structure is not affected.
498 > A series of different initial conditions with a range of surface
499 > coverages was prepared and solvated with various with both of the
500 > solvent molecules. These systems were then equilibrated and their
501 > interfacial thermal conductivity was measured with the NIVS
502 > algorithm. Figure \ref{coverage} demonstrates the trend of conductance
503 > with respect to surface coverage.
504  
505 < Our NIVS algorithm allows change of unphysical thermal flux both in
506 < direction and in quantity. This feature extends our investigation of
507 < interfacial thermal conductance. However, the magnitude of this
508 < thermal flux is not arbitary if one aims to obtain a stable and
509 < reliable thermal gradient. A temperature profile would be
510 < substantially affected by noise when $|J_z|$ has a much too low
511 < magnitude; while an excessively large $|J_z|$ that overwhelms the
515 < conductance capacity of the interface would prevent a thermal gradient
516 < to reach a stablized steady state. NIVS has the advantage of allowing
517 < $J$ to vary in a wide range such that the optimal flux range for $G$
518 < measurement can generally be simulated by the algorithm. Within the
519 < optimal range, we were able to study how $G$ would change according to
520 < the thermal flux across the interface. For our simulations, we denote
521 < $J_z$ to be positive when the physical thermal flux is from the liquid
522 < to metal, and negative vice versa. The $G$'s measured under different
523 < $J_z$ is listed in Table \ref{AuThiolHexaneUA} and
524 < \ref{AuThiolToluene}. These results do not suggest that $G$ is
525 < dependent on $J_z$ within this flux range. The linear response of flux
526 < to thermal gradient simplifies our investigations in that we can rely
527 < on $G$ measurement with only a couple $J_z$'s and do not need to test
528 < a large series of fluxes.
505 > \begin{figure}
506 > \includegraphics[width=\linewidth]{coverage}
507 > \caption{Comparison of interfacial thermal conductivity ($G$) values
508 >  for the Au-butanethiol/solvent interface with various UA models and
509 >  different capping agent coverages at $\langle T\rangle\sim$200K.}
510 > \label{coverage}
511 > \end{figure}
512  
513 + In partially covered surfaces, the derivative definition for
514 + $G^\prime$ (Eq. \ref{derivativeG}) becomes difficult to apply, as the
515 + location of maximum change of $\lambda$ becomes washed out.  The
516 + discrete definition (Eq. \ref{discreteG}) is easier to apply, as the
517 + Gibbs dividing surface is still well-defined. Therefore, $G$ (not
518 + $G^\prime$) was used in this section.
519 +
520 + From Figure \ref{coverage}, one can see the significance of the
521 + presence of capping agents. When even a small fraction of the Au(111)
522 + surface sites are covered with butanethiols, the conductivity exhibits
523 + an enhancement by at least a factor of 3.  Cappping agents are clearly
524 + playing a major role in thermal transport at metal / organic solvent
525 + surfaces.
526 +
527 + We note a non-monotonic behavior in the interfacial conductance as a
528 + function of surface coverage. The maximum conductance (largest $G$)
529 + happens when the surfaces are about 75\% covered with butanethiol
530 + caps.  The reason for this behavior is not entirely clear.  One
531 + explanation is that incomplete butanethiol coverage allows small gaps
532 + between butanethiols to form. These gaps can be filled by transient
533 + solvent molecules.  These solvent molecules couple very strongly with
534 + the hot capping agent molecules near the surface, and can then carry
535 + away (diffusively) the excess thermal energy from the surface.
536 +
537 + There appears to be a competition between the conduction of the
538 + thermal energy away from the surface by the capping agents (enhanced
539 + by greater coverage) and the coupling of the capping agents with the
540 + solvent (enhanced by interdigitation at lower coverages).  This
541 + competition would lead to the non-monotonic coverage behavior observed
542 + here.
543 +
544 + Results for rigid body toluene solvent, as well as the UA hexane, are
545 + within the ranges expected from prior experimental
546 + work.\cite{Wilson:2002uq,cahill:793,PhysRevB.80.195406} This suggests
547 + that explicit hydrogen atoms might not be required for modeling
548 + thermal transport in these systems.  C-H vibrational modes do not see
549 + significant excited state population at low temperatures, and are not
550 + likely to carry lower frequency excitations from the solid layer into
551 + the bulk liquid.
552 +
553 + The toluene solvent does not exhibit the same behavior as hexane in
554 + that $G$ remains at approximately the same magnitude when the capping
555 + coverage increases from 25\% to 75\%.  Toluene, as a rigid planar
556 + molecule, cannot occupy the relatively small gaps between the capping
557 + agents as easily as the chain-like {\it n}-hexane.  The effect of
558 + solvent coupling to the capping agent is therefore weaker in toluene
559 + except at the very lowest coverage levels.  This effect counters the
560 + coverage-dependent conduction of heat away from the metal surface,
561 + leading to a much flatter $G$ vs. coverage trend than is observed in
562 + {\it n}-hexane.
563 +
564 + \subsection{Effects due to Solvent \& Solvent Models}
565 + In addition to UA solvent and capping agent models, AA models have
566 + also been included in our simulations.  In most of this work, the same
567 + (UA or AA) model for solvent and capping agent was used, but it is
568 + also possible to utilize different models for different components.
569 + We have also included isotopic substitutions (Hydrogen to Deuterium)
570 + to decrease the explicit vibrational overlap between solvent and
571 + capping agent. Table \ref{modelTest} summarizes the results of these
572 + studies.
573 +
574   \begin{table*}
575    \begin{minipage}{\linewidth}
576      \begin{center}
577 +      
578 +      \caption{Computed interfacial thermal conductance ($G$ and
579 +        $G^\prime$) values for interfaces using various models for
580 +        solvent and capping agent (or without capping agent) at
581 +        $\langle T\rangle\sim$200K. (D stands for deuterated solvent
582 +        or capping agent molecules; ``Avg.'' denotes results that are
583 +        averages of simulations under different applied thermal flux
584 +        values $(J_z)$. Error estimates are indicated in
585 +        parentheses.)}
586 +      
587 +      \begin{tabular}{llccc}
588 +        \hline\hline
589 +        Butanethiol model & Solvent & $J_z$ & $G$ & $G^\prime$ \\
590 +        (or bare surface) & model & (GW/m$^2$) &
591 +        \multicolumn{2}{c}{(MW/m$^2$/K)} \\
592 +        \hline
593 +        UA    & UA hexane    & Avg. & 131(9)    & 87(10)    \\
594 +              & UA hexane(D) & 1.95 & 153(5)    & 136(13)   \\
595 +              & AA hexane    & Avg. & 131(6)    & 122(10)   \\
596 +              & UA toluene   & 1.96 & 187(16)   & 151(11)   \\
597 +              & AA toluene   & 1.89 & 200(36)   & 149(53)   \\
598 +        \hline
599 +        AA    & UA hexane    & 1.94 & 116(9)    & 129(8)    \\
600 +              & AA hexane    & Avg. & 442(14)   & 356(31)   \\
601 +              & AA hexane(D) & 1.93 & 222(12)   & 234(54)   \\
602 +              & UA toluene   & 1.98 & 125(25)   & 97(60)    \\
603 +              & AA toluene   & 3.79 & 487(56)   & 290(42)   \\
604 +        \hline
605 +        AA(D) & UA hexane    & 1.94 & 158(25)   & 172(4)    \\
606 +              & AA hexane    & 1.92 & 243(29)   & 191(11)   \\
607 +              & AA toluene   & 1.93 & 364(36)   & 322(67)   \\
608 +        \hline
609 +        bare  & UA hexane    & Avg. & 46.5(3.2) & 49.4(4.5) \\
610 +              & UA hexane(D) & 0.98 & 43.9(4.6) & 43.0(2.0) \\
611 +              & AA hexane    & 0.96 & 31.0(1.4) & 29.4(1.3) \\
612 +              & UA toluene   & 1.99 & 70.1(1.3) & 65.8(0.5) \\
613 +        \hline\hline
614 +      \end{tabular}
615 +      \label{modelTest}
616 +    \end{center}
617 +  \end{minipage}
618 + \end{table*}
619 +
620 + To facilitate direct comparison between force fields, systems with the
621 + same capping agent and solvent were prepared with the same length
622 + scales for the simulation cells.
623 +
624 + On bare metal / solvent surfaces, different force field models for
625 + hexane yield similar results for both $G$ and $G^\prime$, and these
626 + two definitions agree with each other very well. This is primarily an
627 + indicator of weak interactions between the metal and the solvent, and
628 + is a typical case for acoustic impedance mismatch between these two
629 + phases.  
630 +
631 + For the fully-covered surfaces, the choice of force field for the
632 + capping agent and solvent has a large impact on the calulated values
633 + of $G$ and $G^\prime$.  The AA thiol to AA solvent conductivities are
634 + much larger than their UA to UA counterparts, and these values exceed
635 + the experimental estimates by a large measure.  The AA force field
636 + allows significant energy to go into C-H (or C-D) stretching modes,
637 + and since these modes are high frequency, this non-quantum behavior is
638 + likely responsible for the overestimate of the conductivity.  Compared
639 + to the AA model, the UA model yields more reasonable conductivity
640 + values with much higher computational efficiency.
641 +
642 + \subsubsection{Are electronic excitations in the metal important?}
643 + Because they lack electronic excitations, the QSC and related embedded
644 + atom method (EAM) models for gold are known to predict unreasonably
645 + low values for bulk conductivity
646 + ($\lambda$).\cite{kuang:164101,ISI:000207079300006,Clancy:1992} If the
647 + conductance between the phases ($G$) is governed primarily by phonon
648 + excitation (and not electronic degrees of freedom), one would expect a
649 + classical model to capture most of the interfacial thermal
650 + conductance.  Our results for $G$ and $G^\prime$ indicate that this is
651 + indeed the case, and suggest that the modeling of interfacial thermal
652 + transport depends primarily on the description of the interactions
653 + between the various components at the interface.  When the metal is
654 + chemically capped, the primary barrier to thermal conductivity appears
655 + to be the interface between the capping agent and the surrounding
656 + solvent, so the excitations in the metal have little impact on the
657 + value of $G$.
658 +
659 + \subsection{Effects due to methodology and simulation parameters}
660 +
661 + We have varied the parameters of the simulations in order to
662 + investigate how these factors would affect the computation of $G$.  Of
663 + particular interest are: 1) the length scale for the applied thermal
664 + gradient (modified by increasing the amount of solvent in the system),
665 + 2) the sign and magnitude of the applied thermal flux, 3) the average
666 + temperature of the simulation (which alters the solvent density during
667 + equilibration), and 4) the definition of the interfacial conductance
668 + (Eqs. (\ref{discreteG}) or (\ref{derivativeG})) used in the
669 + calculation.
670 +
671 + Systems of different lengths were prepared by altering the number of
672 + solvent molecules and extending the length of the box along the $z$
673 + axis to accomodate the extra solvent.  Equilibration at the same
674 + temperature and pressure conditions led to nearly identical surface
675 + areas ($L_x$ and $L_y$) available to the metal and capping agent,
676 + while the extra solvent served mainly to lengthen the axis that was
677 + used to apply the thermal flux.  For a given value of the applied
678 + flux, the different $z$ length scale has only a weak effect on the
679 + computed conductivities (Table \ref{AuThiolHexaneUA}).
680 +
681 + \subsubsection{Effects of applied flux}
682 + The NIVS algorithm allows changes in both the sign and magnitude of
683 + the applied flux.  It is possible to reverse the direction of heat
684 + flow simply by changing the sign of the flux, and thermal gradients
685 + which would be difficult to obtain experimentally ($5$ K/\AA) can be
686 + easily simulated.  However, the magnitude of the applied flux is not
687 + arbitary if one aims to obtain a stable and reliable thermal gradient.
688 + A temperature gradient can be lost in the noise if $|J_z|$ is too
689 + small, and excessive $|J_z|$ values can cause phase transitions if the
690 + extremes of the simulation cell become widely separated in
691 + temperature.  Also, if $|J_z|$ is too large for the bulk conductivity
692 + of the materials, the thermal gradient will never reach a stable
693 + state.  
694 +
695 + Within a reasonable range of $J_z$ values, we were able to study how
696 + $G$ changes as a function of this flux.  In what follows, we use
697 + positive $J_z$ values to denote the case where energy is being
698 + transferred by the method from the metal phase and into the liquid.
699 + The resulting gradient therefore has a higher temperature in the
700 + liquid phase.  Negative flux values reverse this transfer, and result
701 + in higher temperature metal phases.  The conductance measured under
702 + different applied $J_z$ values is listed in Tables
703 + \ref{AuThiolHexaneUA} and \ref{AuThiolToluene}. These results do not
704 + indicate that $G$ depends strongly on $J_z$ within this flux
705 + range. The linear response of flux to thermal gradient simplifies our
706 + investigations in that we can rely on $G$ measurement with only a
707 + small number $J_z$ values.  
708 +
709 + \begin{table*}
710 +  \begin{minipage}{\linewidth}
711 +    \begin{center}
712        \caption{Computed interfacial thermal conductivity ($G$ and
713          $G^\prime$) values for the 100\% covered Au-butanethiol/hexane
714          interfaces with UA model and different hexane molecule numbers
# Line 538 | Line 717 | a large series of fluxes.
717        
718        \begin{tabular}{ccccccc}
719          \hline\hline
720 <        $\langle T\rangle$ & $N_{hexane}$ & Fixed & $\rho_{hexane}$ &
720 >        $\langle T\rangle$ & $N_{hexane}$  & $\rho_{hexane}$ &
721          $J_z$ & $G$ & $G^\prime$ \\
722 <        (K) & & $L_x$ \& $L_y$? & (g/cm$^3$) & (GW/m$^2$) &
722 >        (K) &  & (g/cm$^3$) & (GW/m$^2$) &
723          \multicolumn{2}{c}{(MW/m$^2$/K)} \\
724          \hline
725 <        200 & 266 & No  & 0.672 & -0.96 & 102(3)    & 80.0(0.8) \\
726 <            & 200 & Yes & 0.694 &  1.92 & 129(11)   & 87.3(0.3) \\
727 <            &     & Yes & 0.672 &  1.93 & 131(16)   & 78(13)    \\
728 <            &     & No  & 0.688 &  0.96 & 125(16)   & 90.2(15)  \\
729 <            &     &     &       &  1.91 & 139(10)   & 101(10)   \\
730 <            &     &     &       &  2.83 & 141(6)    & 89.9(9.8) \\
552 <            & 166 & Yes & 0.679 &  0.97 & 115(19)   & 69(18)    \\
553 <            &     &     &       &  1.94 & 125(9)    & 87.1(0.2) \\
554 <            &     & No  & 0.681 &  0.97 & 141(30)   & 78(22)    \\
555 <            &     &     &       &  1.92 & 138(4)    & 98.9(9.5) \\
725 >        200 & 266 &  0.672 & -0.96 & 102(3)    & 80.0(0.8) \\
726 >            & 200 &  0.688 &  0.96 & 125(16)   & 90.2(15)  \\
727 >            &     &        &  1.91 & 139(10)   & 101(10)   \\
728 >            &     &        &  2.83 & 141(6)    & 89.9(9.8) \\
729 >            & 166 &  0.681 &  0.97 & 141(30)   & 78(22)    \\
730 >            &     &        &  1.92 & 138(4)    & 98.9(9.5) \\
731          \hline
732 <        250 & 200 & No  & 0.560 &  0.96 & 75(10)    & 61.8(7.3) \\
733 <            &     &     &       & -0.95 & 49.4(0.3) & 45.7(2.1) \\
734 <            & 166 & Yes & 0.570 &  0.98 & 79.0(3.5) & 62.9(3.0) \\
735 <            &     & No  & 0.569 &  0.97 & 80.3(0.6) & 67(11)    \\
736 <            &     &     &       &  1.44 & 76.2(5.0) & 64.8(3.8) \\
737 <            &     &     &       & -0.95 & 56.4(2.5) & 54.4(1.1) \\
563 <            &     &     &       & -1.85 & 47.8(1.1) & 53.5(1.5) \\
732 >        250 & 200 &  0.560 &  0.96 & 75(10)    & 61.8(7.3) \\
733 >            &     &        & -0.95 & 49.4(0.3) & 45.7(2.1) \\
734 >            & 166 &  0.569 &  0.97 & 80.3(0.6) & 67(11)    \\
735 >            &     &        &  1.44 & 76.2(5.0) & 64.8(3.8) \\
736 >            &     &        & -0.95 & 56.4(2.5) & 54.4(1.1) \\
737 >            &     &        & -1.85 & 47.8(1.1) & 53.5(1.5) \\
738          \hline\hline
739        \end{tabular}
740        \label{AuThiolHexaneUA}
# Line 568 | Line 742 | Furthermore, we also attempted to increase system aver
742    \end{minipage}
743   \end{table*}
744  
745 < Furthermore, we also attempted to increase system average temperatures
746 < to above 200K. These simulations are first equilibrated in the NPT
747 < ensemble under normal pressure. As stated above, the TraPPE-UA model
748 < for hexane tends to predict a lower boiling point. In our simulations,
749 < hexane had diffculty to remain in liquid phase when NPT equilibration
750 < temperature is higher than 250K. Additionally, the equilibrated liquid
751 < hexane density under 250K becomes lower than experimental value. This
752 < expanded liquid phase leads to lower contact between hexane and
753 < butanethiol as well.[MAY NEED SLAB DENSITY FIGURE]
580 < And this reduced contact would
581 < probably be accountable for a lower interfacial thermal conductance,
582 < as shown in Table \ref{AuThiolHexaneUA}.
745 > The sign of $J_z$ is a different matter, however, as this can alter
746 > the temperature on the two sides of the interface. The average
747 > temperature values reported are for the entire system, and not for the
748 > liquid phase, so at a given $\langle T \rangle$, the system with
749 > positive $J_z$ has a warmer liquid phase.  This means that if the
750 > liquid carries thermal energy via convective transport, {\it positive}
751 > $J_z$ values will result in increased molecular motion on the liquid
752 > side of the interface, and this will increase the measured
753 > conductivity.
754  
755 < A similar study for TraPPE-UA toluene agrees with the above result as
585 < well. Having a higher boiling point, toluene tends to remain liquid in
586 < our simulations even equilibrated under 300K in NPT
587 < ensembles. Furthermore, the expansion of the toluene liquid phase is
588 < not as significant as that of the hexane. This prevents severe
589 < decrease of liquid-capping agent contact and the results (Table
590 < \ref{AuThiolToluene}) show only a slightly decreased interface
591 < conductance. Therefore, solvent-capping agent contact should play an
592 < important role in the thermal transport process across the interface
593 < in that higher degree of contact could yield increased conductance.
755 > \subsubsection{Effects due to average temperature}
756  
757 + We also studied the effect of average system temperature on the
758 + interfacial conductance.  The simulations are first equilibrated in
759 + the NPT ensemble at 1 atm.  The TraPPE-UA model for hexane tends to
760 + predict a lower boiling point (and liquid state density) than
761 + experiments.  This lower-density liquid phase leads to reduced contact
762 + between the hexane and butanethiol, and this accounts for our
763 + observation of lower conductance at higher temperatures as shown in
764 + Table \ref{AuThiolHexaneUA}.  In raising the average temperature from
765 + 200K to 250K, the density drop of $\sim$20\% in the solvent phase
766 + leads to a $\sim$65\% drop in the conductance. [BUT (125-75)/125 = .4?]
767 +
768 + Similar behavior is observed in the TraPPE-UA model for toluene,
769 + although this model has better agreement with the experimental
770 + densities of toluene.  The expansion of the toluene liquid phase is
771 + not as significant as that of the hexane (8.3\% over 100K), and this
772 + limits the effect to $\sim$20\% drop in thermal conductivity  (Table
773 + \ref{AuThiolToluene}).
774 +
775 + Although we have not mapped out the behavior at a large number of
776 + temperatures, is clear that there will be a strong temperature
777 + dependence in the interfacial conductance when the physical properties
778 + of one side of the interface (notably the density) change rapidly as a
779 + function of temperature.
780 +
781   \begin{table*}
782    \begin{minipage}{\linewidth}
783      \begin{center}
# Line 618 | Line 804 | Besides lower interfacial thermal conductance, surface
804    \end{minipage}
805   \end{table*}
806  
807 < Besides lower interfacial thermal conductance, surfaces in relatively
808 < high temperatures are susceptible to reconstructions, when
809 < butanethiols have a full coverage on the Au(111) surface. These
810 < reconstructions include surface Au atoms migrated outward to the S
811 < atom layer, and butanethiol molecules embedded into the original
812 < surface Au layer. The driving force for this behavior is the strong
813 < Au-S interactions in our simulations. And these reconstructions lead
814 < to higher ratio of Au-S attraction and thus is energetically
815 < favorable. Furthermore, this phenomenon agrees with experimental
816 < results\cite{doi:10.1021/j100035a033,doi:10.1021/la026493y}. Vlugt
817 < {\it et al.} had kept their Au(111) slab rigid so that their
818 < simulations can reach 300K without surface reconstructions. Without
819 < this practice, simulating 100\% thiol covered interfaces under higher
820 < temperatures could hardly avoid surface reconstructions. However, our
821 < measurement is based on assuming homogeneity on $x$ and $y$ dimensions
822 < so that measurement of $T$ at particular $z$ would be an effective
823 < average of the particles of the same type. Since surface
638 < reconstructions could eliminate the original $x$ and $y$ dimensional
639 < homogeneity, measurement of $G$ is more difficult to conduct under
640 < higher temperatures. Therefore, most of our measurements are
641 < undertaken at $\langle T\rangle\sim$200K.
807 > Besides the lower interfacial thermal conductance, surfaces at
808 > relatively high temperatures are susceptible to reconstructions,
809 > particularly when butanethiols fully cover the Au(111) surface. These
810 > reconstructions include surface Au atoms which migrate outward to the
811 > S atom layer, and butanethiol molecules which embed into the surface
812 > Au layer. The driving force for this behavior is the strong Au-S
813 > interactions which are modeled here with a deep Lennard-Jones
814 > potential. This phenomenon agrees with reconstructions that have beeen
815 > experimentally
816 > observed.\cite{doi:10.1021/j100035a033,doi:10.1021/la026493y}.  Vlugt
817 > {\it et al.} kept their Au(111) slab rigid so that their simulations
818 > could reach 300K without surface
819 > reconstructions.\cite{vlugt:cpc2007154} Since surface reconstructions
820 > blur the interface, the measurement of $G$ becomes more difficult to
821 > conduct at higher temperatures.  For this reason, most of our
822 > measurements are undertaken at $\langle T\rangle\sim$200K where
823 > reconstruction is minimized.
824  
825   However, when the surface is not completely covered by butanethiols,
826 < the simulated system is more resistent to the reconstruction
827 < above. Our Au-butanethiol/toluene system had the Au(111) surfaces 90\%
828 < covered by butanethiols, but did not see this above phenomena even at
829 < $\langle T\rangle\sim$300K. The empty three-fold sites not occupied by
830 < capping agents could help prevent surface reconstruction in that they
831 < provide other means of capping agent relaxation. It is observed that
832 < butanethiols can migrate to their neighbor empty sites during a
833 < simulation. Therefore, we were able to obtain $G$'s for these
652 < interfaces even at a relatively high temperature without being
653 < affected by surface reconstructions.
826 > the simulated system appears to be more resistent to the
827 > reconstruction. Our Au / butanethiol / toluene system had the Au(111)
828 > surfaces 90\% covered by butanethiols, but did not see this above
829 > phenomena even at $\langle T\rangle\sim$300K.  That said, we did
830 > observe butanethiols migrating to neighboring three-fold sites during
831 > a simulation.  Since the interface persisted in these simulations,
832 > were able to obtain $G$'s for these interfaces even at a relatively
833 > high temperature without being affected by surface reconstructions.
834  
835 < \subsection{Influence of Capping Agent Coverage on $G$}
656 < To investigate the influence of butanethiol coverage on interfacial
657 < thermal conductance, a series of different coverage Au-butanethiol
658 < surfaces is prepared and solvated with various organic
659 < molecules. These systems are then equilibrated and their interfacial
660 < thermal conductivity are measured with our NIVS algorithm. Figure
661 < \ref{coverage} demonstrates the trend of conductance change with
662 < respect to different coverages of butanethiol. To study the isotope
663 < effect in interfacial thermal conductance, deuterated UA-hexane is
664 < included as well.
835 > \section{Discussion}
836  
837 < It turned out that with partial covered butanethiol on the Au(111)
838 < surface, the derivative definition for $G^\prime$
839 < (Eq. \ref{derivativeG}) was difficult to apply, due to the difficulty
840 < in locating the maximum of change of $\lambda$. Instead, the discrete
841 < definition (Eq. \ref{discreteG}) is easier to apply, as the Gibbs
842 < deviding surface can still be well-defined. Therefore, $G$ (not
672 < $G^\prime$) was used for this section.
837 > The primary result of this work is that the capping agent acts as an
838 > efficient thermal coupler between solid and solvent phases.  One of
839 > the ways the capping agent can carry out this role is to down-shift
840 > between the phonon vibrations in the solid (which carry the heat from
841 > the gold) and the molecular vibrations in the liquid (which carry some
842 > of the heat in the solvent).
843  
844 < From Figure \ref{coverage}, one can see the significance of the
845 < presence of capping agents. Even when a fraction of the Au(111)
846 < surface sites are covered with butanethiols, the conductivity would
847 < see an enhancement by at least a factor of 3. This indicates the
848 < important role cappping agent is playing for thermal transport
849 < phenomena on metal / organic solvent surfaces.
844 > To investigate the mechanism of interfacial thermal conductance, the
845 > vibrational power spectrum was computed. Power spectra were taken for
846 > individual components in different simulations. To obtain these
847 > spectra, simulations were run after equilibration in the
848 > microcanonical (NVE) ensemble and without a thermal
849 > gradient. Snapshots of configurations were collected at a frequency
850 > that is higher than that of the fastest vibrations occuring in the
851 > simulations. With these configurations, the velocity auto-correlation
852 > functions can be computed:
853 > \begin{equation}
854 > C_A (t) = \langle\vec{v}_A (t)\cdot\vec{v}_A (0)\rangle
855 > \label{vCorr}
856 > \end{equation}
857 > The power spectrum is constructed via a Fourier transform of the
858 > symmetrized velocity autocorrelation function,
859 > \begin{equation}
860 >  \hat{f}(\omega) =
861 >  \int_{-\infty}^{\infty} C_A (t) e^{-2\pi it\omega}\,dt
862 > \label{fourier}
863 > \end{equation}
864  
865 < Interestingly, as one could observe from our results, the maximum
866 < conductance enhancement (largest $G$) happens while the surfaces are
867 < about 75\% covered with butanethiols. This again indicates that
868 < solvent-capping agent contact has an important role of the thermal
869 < transport process. Slightly lower butanethiol coverage allows small
870 < gaps between butanethiols to form. And these gaps could be filled with
871 < solvent molecules, which acts like ``heat conductors'' on the
872 < surface. The higher degree of interaction between these solvent
873 < molecules and capping agents increases the enhancement effect and thus
874 < produces a higher $G$ than densely packed butanethiol arrays. However,
875 < once this maximum conductance enhancement is reached, $G$ decreases
692 < when butanethiol coverage continues to decrease. Each capping agent
693 < molecule reaches its maximum capacity for thermal
694 < conductance. Therefore, even higher solvent-capping agent contact
695 < would not offset this effect. Eventually, when butanethiol coverage
696 < continues to decrease, solvent-capping agent contact actually
697 < decreases with the disappearing of butanethiol molecules. In this
698 < case, $G$ decrease could not be offset but instead accelerated. [NEED
699 < SNAPSHOT SHOWING THE PHENOMENA]
865 > \subsection{The role of specific vibrations}
866 > The vibrational spectra for gold slabs in different environments are
867 > shown as in Figure \ref{specAu}. Regardless of the presence of
868 > solvent, the gold surfaces which are covered by butanethiol molecules
869 > exhibit an additional peak observed at a frequency of
870 > $\sim$165cm$^{-1}$.  We attribute this peak to the S-Au bonding
871 > vibration. This vibration enables efficient thermal coupling of the
872 > surface Au layer to the capping agents. Therefore, in our simulations,
873 > the Au / S interfaces do not appear to be the primary barrier to
874 > thermal transport when compared with the butanethiol / solvent
875 > interfaces.
876  
877 < A comparison of the results obtained from differenet organic solvents
878 < can also provide useful information of the interfacial thermal
879 < transport process. The deuterated hexane (UA) results do not appear to
880 < be much different from those of normal hexane (UA), given that
881 < butanethiol (UA) is non-deuterated for both solvents. These UA model
882 < studies, even though eliminating C-H vibration samplings, still have
883 < C-C vibrational frequencies different from each other. However, these
708 < differences in the infrared range do not seem to produce an observable
709 < difference for the results of $G$. [MAY NEED SPECTRA FIGURE]
710 <
711 < Furthermore, results for rigid body toluene solvent, as well as other
712 < UA-hexane solvents, are reasonable within the general experimental
713 < ranges[CITATIONS]. This suggests that explicit hydrogen might not be a
714 < required factor for modeling thermal transport phenomena of systems
715 < such as Au-thiol/organic solvent.
716 <
717 < However, results for Au-butanethiol/toluene do not show an identical
718 < trend with those for Au-butanethiol/hexane in that $G$ remains at
719 < approximately the same magnitue when butanethiol coverage differs from
720 < 25\% to 75\%. This might be rooted in the molecule shape difference
721 < for planar toluene and chain-like {\it n}-hexane. Due to this
722 < difference, toluene molecules have more difficulty in occupying
723 < relatively small gaps among capping agents when their coverage is not
724 < too low. Therefore, the solvent-capping agent contact may keep
725 < increasing until the capping agent coverage reaches a relatively low
726 < level. This becomes an offset for decreasing butanethiol molecules on
727 < its effect to the process of interfacial thermal transport. Thus, one
728 < can see a plateau of $G$ vs. butanethiol coverage in our results.
729 <
730 < \begin{figure}
731 < \includegraphics[width=\linewidth]{coverage}
732 < \caption{Comparison of interfacial thermal conductivity ($G$) values
733 <  for the Au-butanethiol/solvent interface with various UA models and
734 <  different capping agent coverages at $\langle T\rangle\sim$200K
735 <  using certain energy flux respectively.}
736 < \label{coverage}
877 > \begin{figure}
878 > \includegraphics[width=\linewidth]{vibration}
879 > \caption{Vibrational power spectra for gold in different solvent
880 >  environments.  The presence of the butanethiol capping molecules
881 >  adds a vibrational peak at $\sim$165cm$^{-1}$. The butanethiol
882 >  spectra exhibit a corresponding peak.}
883 > \label{specAu}
884   \end{figure}
885  
886 < \subsection{Influence of Chosen Molecule Model on $G$}
887 < [MAY COMBINE W MECHANISM STUDY]
886 > Also in this figure, we show the vibrational power spectrum for the
887 > bound butanethiol molecules, which also exhibits the same
888 > $\sim$165cm$^{-1}$ peak.
889  
890 < In addition to UA solvent/capping agent models, AA models are included
891 < in our simulations as well. Besides simulations of the same (UA or AA)
892 < model for solvent and capping agent, different models can be applied
893 < to different components. Furthermore, regardless of models chosen,
894 < either the solvent or the capping agent can be deuterated, similar to
895 < the previous section. Table \ref{modelTest} summarizes the results of
896 < these studies.
890 > \subsection{Overlap of power spectra}
891 > A comparison of the results obtained from the two different organic
892 > solvents can also provide useful information of the interfacial
893 > thermal transport process.  In particular, the vibrational overlap
894 > between the butanethiol and the organic solvents suggests a highly
895 > efficient thermal exchange between these components.  Very high
896 > thermal conductivity was observed when AA models were used and C-H
897 > vibrations were treated classically.  The presence of extra degrees of
898 > freedom in the AA force field yields higher heat exchange rates
899 > between the two phases and results in a much higher conductivity than
900 > in the UA force field.
901  
902 < \begin{table*}
903 <  \begin{minipage}{\linewidth}
904 <    \begin{center}
905 <      
906 <      \caption{Computed interfacial thermal conductivity ($G$ and
755 <        $G^\prime$) values for interfaces using various models for
756 <        solvent and capping agent (or without capping agent) at
757 <        $\langle T\rangle\sim$200K. (D stands for deuterated solvent
758 <        or capping agent molecules; ``Avg.'' denotes results that are
759 <        averages of simulations under different $J_z$'s. Error
760 <        estimates indicated in parenthesis.)}
761 <      
762 <      \begin{tabular}{llccc}
763 <        \hline\hline
764 <        Butanethiol model & Solvent & $J_z$ & $G$ & $G^\prime$ \\
765 <        (or bare surface) & model & (GW/m$^2$) &
766 <        \multicolumn{2}{c}{(MW/m$^2$/K)} \\
767 <        \hline
768 <        UA    & UA hexane    & Avg. & 131(9)    & 87(10)    \\
769 <              & UA hexane(D) & 1.95 & 153(5)    & 136(13)   \\
770 <              & AA hexane    & Avg. & 131(6)    & 122(10)   \\
771 <              & UA toluene   & 1.96 & 187(16)   & 151(11)   \\
772 <              & AA toluene   & 1.89 & 200(36)   & 149(53)   \\
773 <        \hline
774 <        AA    & UA hexane    & 1.94 & 116(9)    & 129(8)    \\
775 <              & AA hexane    & Avg. & 442(14)   & 356(31)   \\
776 <              & AA hexane(D) & 1.93 & 222(12)   & 234(54)   \\
777 <              & UA toluene   & 1.98 & 125(25)   & 97(60)    \\
778 <              & AA toluene   & 3.79 & 487(56)   & 290(42)   \\
779 <        \hline
780 <        AA(D) & UA hexane    & 1.94 & 158(25)   & 172(4)    \\
781 <              & AA hexane    & 1.92 & 243(29)   & 191(11)   \\
782 <              & AA toluene   & 1.93 & 364(36)   & 322(67)   \\
783 <        \hline
784 <        bare  & UA hexane    & Avg. & 46.5(3.2) & 49.4(4.5) \\
785 <              & UA hexane(D) & 0.98 & 43.9(4.6) & 43.0(2.0) \\
786 <              & AA hexane    & 0.96 & 31.0(1.4) & 29.4(1.3) \\
787 <              & UA toluene   & 1.99 & 70.1(1.3) & 65.8(0.5) \\
788 <        \hline\hline
789 <      \end{tabular}
790 <      \label{modelTest}
791 <    \end{center}
792 <  \end{minipage}
793 < \end{table*}
902 > The similarity in the vibrational modes available to solvent and
903 > capping agent can be reduced by deuterating one of the two components
904 > (Fig. \ref{aahxntln}).  Once either the hexanes or the butanethiols
905 > are deuterated, one can observe a significantly lower $G$ and
906 > $G^\prime$ values (Table \ref{modelTest}).
907  
908 < To facilitate direct comparison, the same system with differnt models
909 < for different components uses the same length scale for their
910 < simulation cells. Without the presence of capping agent, using
911 < different models for hexane yields similar results for both $G$ and
912 < $G^\prime$, and these two definitions agree with eath other very
913 < well. This indicates very weak interaction between the metal and the
914 < solvent, and is a typical case for acoustic impedance mismatch between
915 < these two phases.
908 > \begin{figure}
909 > \includegraphics[width=\linewidth]{aahxntln}
910 > \caption{Spectra obtained for all-atom (AA) Au / butanethiol / solvent
911 >  systems. When butanethiol is deuterated (lower left), its
912 >  vibrational overlap with hexane decreases significantly.  Since
913 >  aromatic molecules and the butanethiol are vibrationally dissimilar,
914 >  the change is not as dramatic when toluene is the solvent (right).}
915 > \label{aahxntln}
916 > \end{figure}
917  
918 < As for Au(111) surfaces completely covered by butanethiols, the choice
805 < of models for capping agent and solvent could impact the measurement
806 < of $G$ and $G^\prime$ quite significantly. For Au-butanethiol/hexane
807 < interfaces, using AA model for both butanethiol and hexane yields
808 < substantially higher conductivity values than using UA model for at
809 < least one component of the solvent and capping agent, which exceeds
810 < the general range of experimental measurement results. This is
811 < probably due to the classically treated C-H vibrations in the AA
812 < model, which should not be appreciably populated at normal
813 < temperatures. In comparison, once either the hexanes or the
814 < butanethiols are deuterated, one can see a significantly lower $G$ and
815 < $G^\prime$. In either of these cases, the C-H(D) vibrational overlap
816 < between the solvent and the capping agent is removed.
817 < [MAY NEED SPECTRA FIGURE] Conclusively, the
818 < improperly treated C-H vibration in the AA model produced
819 < over-predicted results accordingly. Compared to the AA model, the UA
820 < model yields more reasonable results with higher computational
821 < efficiency.
822 <
823 < However, for Au-butanethiol/toluene interfaces, having the AA
918 > For the Au / butanethiol / toluene interfaces, having the AA
919   butanethiol deuterated did not yield a significant change in the
920 < measurement results. Compared to the C-H vibrational overlap between
921 < hexane and butanethiol, both of which have alkyl chains, that overlap
922 < between toluene and butanethiol is not so significant and thus does
923 < not have as much contribution to the ``Intramolecular Vibration
829 < Redistribution''[CITE HASE]. Conversely, extra degrees of freedom such
830 < as the C-H vibrations could yield higher heat exchange rate between
831 < these two phases and result in a much higher conductivity.
920 > measured conductance. Compared to the C-H vibrational overlap between
921 > hexane and butanethiol, both of which have alkyl chains, the overlap
922 > between toluene and butanethiol is not as significant and thus does
923 > not contribute as much to the heat exchange process.
924  
925 < Although the QSC model for Au is known to predict an overly low value
926 < for bulk metal gold conductivity\cite{kuang:164101}, our computational
927 < results for $G$ and $G^\prime$ do not seem to be affected by this
928 < drawback of the model for metal. Instead, our results suggest that the
929 < modeling of interfacial thermal transport behavior relies mainly on
930 < the accuracy of the interaction descriptions between components
931 < occupying the interfaces.
925 > Previous observations of Zhang {\it et al.}\cite{hase:2010} indicate
926 > that the {\it intra}molecular heat transport due to alkylthiols is
927 > highly efficient.  Combining our observations with those of Zhang {\it
928 >  et al.}, it appears that butanethiol acts as a channel to expedite
929 > heat flow from the gold surface and into the alkyl chain.  The
930 > acoustic impedance mismatch between the metal and the liquid phase can
931 > therefore be effectively reduced with the presence of suitable capping
932 > agents.
933  
934 < \subsection{Mechanism of Interfacial Thermal Conductance Enhancement
935 <  by Capping Agent}
936 < [OR: Vibrational Spectrum Study on Conductance Mechanism]
934 > Deuterated models in the UA force field did not decouple the thermal
935 > transport as well as in the AA force field.  The UA models, even
936 > though they have eliminated the high frequency C-H vibrational
937 > overlap, still have significant overlap in the lower-frequency
938 > portions of the infrared spectra (Figure \ref{uahxnua}).  Deuterating
939 > the UA models did not decouple the low frequency region enough to
940 > produce an observable difference for the results of $G$ (Table
941 > \ref{modelTest}).
942  
845 [MAY INTRODUCE PROTOCOL IN METHODOLOGY/COMPUTATIONAL DETAIL, EQN'S]
846
847 To investigate the mechanism of this interfacial thermal conductance,
848 the vibrational spectra of various gold systems were obtained and are
849 shown as in the upper panel of Fig. \ref{vibration}. To obtain these
850 spectra, one first runs a simulation in the NVE ensemble and collects
851 snapshots of configurations; these configurations are used to compute
852 the velocity auto-correlation functions, which is used to construct a
853 power spectrum via a Fourier transform.
854
855 [MAY RELATE TO HASE'S]
856 The gold surfaces covered by butanethiol molecules, compared to bare
857 gold surfaces, exhibit an additional peak observed at the frequency of
858 $\sim$170cm$^{-1}$, which is attributed to the S-Au bonding
859 vibration. This vibration enables efficient thermal transport from
860 surface Au layer to the capping agents.
861 [MAY PUT IN OTHER SECTION] Simultaneously, as shown in
862 the lower panel of Fig. \ref{vibration}, the large overlap of the
863 vibration spectra of butanethiol and hexane in the All-Atom model,
864 including the C-H vibration, also suggests high thermal exchange
865 efficiency. The combination of these two effects produces the drastic
866 interfacial thermal conductance enhancement in the All-Atom model.
867
868 [NEED SEPARATE FIGURE. MAY NEED TO CONVERT TO JPEG]
943   \begin{figure}
944 < \includegraphics[width=\linewidth]{vibration}
945 < \caption{Vibrational spectra obtained for gold in different
946 <  environments.}
947 < \label{vibration}
944 > \includegraphics[width=\linewidth]{uahxnua}
945 > \caption{Vibrational spectra obtained for normal (upper) and
946 >  deuterated (lower) hexane in Au-butanethiol/hexane
947 >  systems. Butanethiol spectra are shown as reference. Both hexane and
948 >  butanethiol were using United-Atom models.}
949 > \label{uahxnua}
950   \end{figure}
951  
876 [MAY ADD COMPARISON OF G AND G', AU SLAB WIDTHS, ETC]
877 % The results show that the two definitions used for $G$ yield
878 % comparable values, though $G^\prime$ tends to be smaller.
879
952   \section{Conclusions}
953 < The NIVS algorithm we developed has been applied to simulations of
954 < Au-butanethiol surfaces with organic solvents. This algorithm allows
955 < effective unphysical thermal flux transferred between the metal and
956 < the liquid phase. With the flux applied, we were able to measure the
957 < corresponding thermal gradient and to obtain interfacial thermal
958 < conductivities. Our simulations have seen significant conductance
959 < enhancement with the presence of capping agent, compared to the bare
960 < gold / liquid interfaces. The acoustic impedance mismatch between the
961 < metal and the liquid phase is effectively eliminated by proper capping
953 > The NIVS algorithm has been applied to simulations of
954 > butanethiol-capped Au(111) surfaces in the presence of organic
955 > solvents. This algorithm allows the application of unphysical thermal
956 > flux to transfer heat between the metal and the liquid phase. With the
957 > flux applied, we were able to measure the corresponding thermal
958 > gradients and to obtain interfacial thermal conductivities. Under
959 > steady states, 2-3 ns trajectory simulations are sufficient for
960 > computation of this quantity.
961 >
962 > Our simulations have seen significant conductance enhancement in the
963 > presence of capping agent, compared with the bare gold / liquid
964 > interfaces. The acoustic impedance mismatch between the metal and the
965 > liquid phase is effectively eliminated by a chemically-bonded capping
966   agent. Furthermore, the coverage precentage of the capping agent plays
967 < an important role in the interfacial thermal transport process.
967 > an important role in the interfacial thermal transport
968 > process. Moderately low coverages allow higher contact between capping
969 > agent and solvent, and thus could further enhance the heat transfer
970 > process, giving a non-monotonic behavior of conductance with
971 > increasing coverage.
972  
973 < Our measurement results, particularly of the UA models, agree with
974 < available experimental data. This indicates that our force field
895 < parameters have a nice description of the interactions between the
896 < particles at the interfaces. AA models tend to overestimate the
973 > Our results, particularly using the UA models, agree well with
974 > available experimental data.  The AA models tend to overestimate the
975   interfacial thermal conductance in that the classically treated C-H
976 < vibration would be overly sampled. Compared to the AA models, the UA
977 < models have higher computational efficiency with satisfactory
978 < accuracy, and thus are preferable in interfacial thermal transport
979 < modelings.
976 > vibrations become too easily populated. Compared to the AA models, the
977 > UA models have higher computational efficiency with satisfactory
978 > accuracy, and thus are preferable in modeling interfacial thermal
979 > transport.
980  
981 < Vlugt {\it et al.} has investigated the surface thiol structures for
982 < nanocrystal gold and pointed out that they differs from those of the
983 < Au(111) surface\cite{vlugt:cpc2007154}. This difference might lead to
984 < change of interfacial thermal transport behavior as well. To
985 < investigate this problem, an effective means to introduce thermal flux
986 < and measure the corresponding thermal gradient is desirable for
987 < simulating structures with spherical symmetry.
981 > Of the two definitions for $G$, the discrete form
982 > (Eq. \ref{discreteG}) was easier to use and gives out relatively
983 > consistent results, while the derivative form (Eq. \ref{derivativeG})
984 > is not as versatile. Although $G^\prime$ gives out comparable results
985 > and follows similar trend with $G$ when measuring close to fully
986 > covered or bare surfaces, the spatial resolution of $T$ profile
987 > required for the use of a derivative form is limited by the number of
988 > bins and the sampling required to obtain thermal gradient information.
989  
990 + Vlugt {\it et al.} have investigated the surface thiol structures for
991 + nanocrystalline gold and pointed out that they differ from those of
992 + the Au(111) surface.\cite{landman:1998,vlugt:cpc2007154} This
993 + difference could also cause differences in the interfacial thermal
994 + transport behavior. To investigate this problem, one would need an
995 + effective method for applying thermal gradients in non-planar
996 + (i.e. spherical) geometries.
997  
998   \section{Acknowledgments}
999   Support for this project was provided by the National Science
1000   Foundation under grant CHE-0848243. Computational time was provided by
1001   the Center for Research Computing (CRC) at the University of Notre
1002 < Dame. \newpage
1002 > Dame.
1003 > \newpage
1004  
1005   \bibliography{interfacial}
1006  

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines