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 3755 by skuang, Fri Jul 29 15:45:14 2011 UTC vs.
Revision 3767 by gezelter, Fri Sep 30 19:37:13 2011 UTC

# Line 28 | Line 28
28  
29   \begin{document}
30  
31 < \title{Simulating interfacial thermal conductance at metal-solvent
32 <  interfaces: the role of chemical capping agents}
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} \\
# Line 44 | Line 44 | Notre Dame, Indiana 46556}
44   \begin{doublespace}
45  
46   \begin{abstract}
47 +  With the Non-Isotropic Velocity Scaling (NIVS) approach to Reverse
48 +  Non-Equilibrium Molecular Dynamics (RNEMD) it is possible to impose
49 +  an unphysical thermal flux between different regions of
50 +  inhomogeneous systems such as solid / liquid interfaces.  We have
51 +  applied NIVS to compute the interfacial thermal conductance at a
52 +  metal / organic solvent interface that has been chemically capped by
53 +  butanethiol molecules.  Our calculations suggest that vibrational
54 +  coupling between the metal and liquid phases is enhanced by the
55 +  capping agents, leading to a greatly enhanced conductivity at the
56 +  interface.  Specifically, the chemical bond between the metal and
57 +  the capping agent introduces a vibrational overlap that is not
58 +  present without the capping agent, and the overlap between the
59 +  vibrational spectra (metal to cap, cap to solvent) provides a
60 +  mechanism for rapid thermal transport across the interface. Our
61 +  calculations also suggest that this is a non-monotonic function of
62 +  the fractional coverage of the surface, as moderate coverages allow
63 +  diffusive heat transport of solvent molecules that have been in
64 +  close contact with the capping agent.
65  
66 < With the Non-Isotropic Velocity Scaling algorithm (NIVS) we have
67 < developed, an unphysical thermal flux can be effectively set up even
50 < for non-homogeneous systems like interfaces in non-equilibrium
51 < molecular dynamics simulations. In this work, this algorithm is
52 < applied for simulating thermal conductance at metal / organic solvent
53 < interfaces with various coverages of butanethiol capping
54 < agents. Different solvents and force field models were tested. Our
55 < results suggest that the United-Atom models are able to provide an
56 < estimate of the interfacial thermal conductivity comparable to
57 < experiments in our simulations with satisfactory computational
58 < efficiency. From our results, the acoustic impedance mismatch between
59 < metal and liquid phase is effectively reduced by the capping
60 < agents, and thus leads to interfacial thermal conductance
61 < enhancement. Furthermore, this effect is closely related to the
62 < capping agent coverage on the metal surfaces and the type of solvent
63 < molecules, and is affected by the models used in the simulations.
64 <
66 > Keywords: non-equilibrium, molecular dynamics, vibrational overlap,
67 > coverage dependent.
68   \end{abstract}
69  
70   \newpage
# Line 73 | Line 76 | Due to the importance of heat flow in nanotechnology,
76   %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
77  
78   \section{Introduction}
79 < Due to the importance of heat flow in nanotechnology, interfacial
80 < thermal conductance has been studied extensively both experimentally
81 < and computationally.\cite{cahill:793} Unlike bulk materials, nanoscale
82 < materials have a significant fraction of their atoms at interfaces,
83 < and the chemical details of these interfaces govern the heat transfer
84 < behavior. Furthermore, the interfaces are
79 > Due to the importance of heat flow (and heat removal) in
80 > nanotechnology, interfacial thermal conductance has been studied
81 > extensively both experimentally and computationally.\cite{cahill:793}
82 > Nanoscale materials have a significant fraction of their atoms at
83 > interfaces, and the chemical details of these interfaces govern the
84 > thermal transport properties.  Furthermore, the interfaces are often
85   heterogeneous (e.g. solid - liquid), which provides a challenge to
86 < traditional methods developed for homogeneous systems.
86 > computational methods which have been developed for homogeneous or
87 > bulk systems.
88  
89 < Experimentally, various interfaces have been investigated for their
90 < thermal conductance. Cahill and coworkers studied nanoscale thermal
89 > Experimentally, the thermal properties of a number of interfaces have
90 > been investigated.  Cahill and coworkers studied nanoscale thermal
91   transport from metal nanoparticle/fluid interfaces, to epitaxial
92 < TiN/single crystal oxides interfaces, to hydrophilic and hydrophobic
92 > TiN/single crystal oxides interfaces, and hydrophilic and hydrophobic
93   interfaces between water and solids with different self-assembled
94   monolayers.\cite{Wilson:2002uq,PhysRevB.67.054302,doi:10.1021/jp048375k,PhysRevLett.96.186101}
95 < Wang {\it et al.} studied heat transport through
96 < long-chain hydrocarbon monolayers on gold substrate at individual
97 < molecular level,\cite{Wang10082007} Schmidt {\it et al.} studied the
98 < role of CTAB on thermal transport between gold nanorods and
99 < solvent,\cite{doi:10.1021/jp8051888} and Juv\'e {\it et al.} studied
100 < the cooling dynamics, which is controlled by thermal interface
101 < resistence of glass-embedded metal
95 > Wang {\it et al.} studied heat transport through long-chain
96 > hydrocarbon monolayers on gold substrate at individual molecular
97 > level,\cite{Wang10082007} Schmidt {\it et al.} studied the role of
98 > cetyltrimethylammonium bromide (CTAB) on the thermal transport between
99 > gold nanorods and solvent,\cite{doi:10.1021/jp8051888} and Juv\'e {\it
100 >  et al.} studied the cooling dynamics, which is controlled by thermal
101 > interface resistance of glass-embedded metal
102   nanoparticles.\cite{PhysRevB.80.195406} Although interfaces are
103   normally considered barriers for heat transport, Alper {\it et al.}
104   suggested that specific ligands (capping agents) could completely
105   eliminate this barrier
106   ($G\rightarrow\infty$).\cite{doi:10.1021/la904855s}
107  
108 < Theoretical and computational models have also been used to study the
108 > The acoustic mismatch model for interfacial conductance utilizes the
109 > acoustic impedance ($Z_a = \rho_a v^s_a$) on both sides of the
110 > interface.\cite{schwartz} Here, $\rho_a$ and $v^s_a$ are the density
111 > and speed of sound in material $a$.  The phonon transmission
112 > probability at the $a-b$ interface is
113 > \begin{equation}
114 > t_{ab} = \frac{4 Z_a Z_b}{\left(Z_a + Z_b \right)^2},
115 > \end{equation}
116 > and the interfacial conductance can then be approximated as
117 > \begin{equation}
118 > G_{ab} \approx \frac{1}{4} C_D v_D t_{ab}
119 > \end{equation}
120 > where $C_D$ is the Debye heat capacity of the hot side, and $v_D$ is
121 > the Debye phonon velocity ($1/v_D^3 = 1/3v_L^3 + 2/3 v_T^3$) where
122 > $v_L$ and $v_T$ are the longitudinal and transverse speeds of sound,
123 > respectively.  For the Au/hexane and Au/toluene interfaces, the
124 > acoustic mismatch model predicts room-temperature $G \approx 87 \mbox{
125 >  and } 129$ MW m$^{-2}$ K$^{-1}$, respectively.  However, it is not
126 > clear how one might apply the acoustic mismatch model to a
127 > chemically-modified surface, particularly when the acoustic properties
128 > of a monolayer film may not be well characterized.
129 >
130 > More precise computational models have also been used to study the
131   interfacial thermal transport in order to gain an understanding of
132   this phenomena at the molecular level. Recently, Hase and coworkers
133   employed Non-Equilibrium Molecular Dynamics (NEMD) simulations to
# Line 114 | Line 140 | methods\cite{MullerPlathe:1997xw,kuang:164101} would h
140   difficult to measure with Equilibrium
141   MD\cite{doi:10.1080/0026897031000068578} or forward NEMD simulation
142   methods. Therefore, the Reverse NEMD (RNEMD)
143 < methods\cite{MullerPlathe:1997xw,kuang:164101} would have the
144 < advantage of applying this difficult to measure flux (while measuring
145 < the resulting gradient), given that the simulation methods being able
146 < to effectively apply an unphysical flux in non-homogeneous systems.
143 > methods\cite{MullerPlathe:1997xw,kuang:164101} would be advantageous
144 > in that they {\it apply} the difficult to measure quantity (flux),
145 > while {\it measuring} the easily-computed quantity (the thermal
146 > gradient).  This is particularly true for inhomogeneous interfaces
147 > where it would not be clear how to apply a gradient {\it a priori}.
148   Garde and coworkers\cite{garde:nl2005,garde:PhysRevLett2009} applied
149   this approach to various liquid interfaces and studied how thermal
150 < conductance (or resistance) is dependent on chemistry details of
151 < interfaces, e.g. hydrophobic and hydrophilic aqueous interfaces.
150 > conductance (or resistance) is dependent on chemical details of a
151 > number of hydrophobic and hydrophilic aqueous interfaces. And
152 > recently, Luo {\it et al.} studied the thermal conductance of
153 > Au-SAM-Au junctions using the same approach, comparing to a constant
154 > temperature difference method.\cite{Luo20101} While this latter
155 > approach establishes more ideal Maxwell-Boltzmann distributions than
156 > previous RNEMD methods, it does not guarantee momentum or kinetic
157 > energy conservation.
158  
159   Recently, we have developed a Non-Isotropic Velocity Scaling (NIVS)
160   algorithm for RNEMD simulations\cite{kuang:164101}. This algorithm
# Line 141 | Line 174 | underlying mechanism for the phenomena was investigate
174   underlying mechanism for the phenomena was investigated.
175  
176   \section{Methodology}
177 < \subsection{Imposd-Flux Methods in MD Simulations}
177 > \subsection{Imposed-Flux Methods in MD Simulations}
178   Steady state MD simulations have an advantage in that not many
179   trajectories are needed to study the relationship between thermal flux
180   and thermal gradients. For systems with low interfacial conductance,
# Line 165 | Line 198 | can be applied between regions of particles of arbitar
198   kinetic energy fluxes without obvious perturbation to the velocity
199   distributions of the simulated systems. Furthermore, this approach has
200   the advantage in heterogeneous interfaces in that kinetic energy flux
201 < can be applied between regions of particles of arbitary identity, and
201 > can be applied between regions of particles of arbitrary identity, and
202   the flux will not be restricted by difference in particle mass.
203  
204   The NIVS algorithm scales the velocity vectors in two separate regions
205 < of a simulation system with respective diagonal scaling matricies. To
206 < determine these scaling factors in the matricies, a set of equations
205 > of a simulation system with respective diagonal scaling matrices. To
206 > determine these scaling factors in the matrices, a set of equations
207   including linear momentum conservation and kinetic energy conservation
208   constraints and target energy flux satisfaction is solved. With the
209   scaling operation applied to the system in a set frequency, bulk
# Line 193 | Line 226 | temperature of the two separated phases.
226   where ${E_{total}}$ is the total imposed non-physical kinetic energy
227   transfer during the simulation and ${\langle T_\mathrm{hot}\rangle}$
228   and ${\langle T_\mathrm{cold}\rangle}$ are the average observed
229 < temperature of the two separated phases.
229 > temperature of the two separated phases.  For an applied flux $J_z$
230 > operating over a simulation time $t$ on a periodically-replicated slab
231 > of dimensions $L_x \times L_y$, $E_{total} = J_z *(t)*(2 L_x L_y)$.
232  
233   When the interfacial conductance is {\it not} small, there are two
234   ways to define $G$. One common way is to assume the temperature is
235   discrete on the two sides of the interface. $G$ can be calculated
236   using the applied thermal flux $J$ and the maximum temperature
237   difference measured along the thermal gradient max($\Delta T$), which
238 < occurs at the Gibbs deviding surface (Figure \ref{demoPic}). This is
238 > occurs at the Gibbs dividing surface (Figure \ref{demoPic}). This is
239   known as the Kapitza conductance, which is the inverse of the Kapitza
240   resistance.
241   \begin{equation}
# Line 213 | Line 248 | resistance.
248   \caption{Interfacial conductance can be calculated by applying an
249    (unphysical) kinetic energy flux between two slabs, one located
250    within the metal and another on the edge of the periodic box.  The
251 <  system responds by forming a thermal response or a gradient.  In
252 <  bulk liquids, this gradient typically has a single slope, but in
253 <  interfacial systems, there are distinct thermal conductivity
254 <  domains.  The interfacial conductance, $G$ is found by measuring the
255 <  temperature gap at the Gibbs dividing surface, or by using second
256 <  derivatives of the thermal profile.}
251 >  system responds by forming a thermal gradient.  In bulk liquids,
252 >  this gradient typically has a single slope, but in interfacial
253 >  systems, there are distinct thermal conductivity domains.  The
254 >  interfacial conductance, $G$ is found by measuring the temperature
255 >  gap at the Gibbs dividing surface, or by using second derivatives of
256 >  the thermal profile.}
257   \label{demoPic}
258   \end{figure}
259  
260 < The other approach is to assume a continuous temperature profile along
261 < the thermal gradient axis (e.g. $z$) and define $G$ at the point where
262 < the magnitude of thermal conductivity ($\lambda$) change reaches its
263 < maximum, given that $\lambda$ is well-defined throughout the space:
264 < \begin{equation}
265 < G^\prime = \Big|\frac{\partial\lambda}{\partial z}\Big|
266 <         = \Big|\frac{\partial}{\partial z}\left(-J_z\Big/
267 <           \left(\frac{\partial T}{\partial z}\right)\right)\Big|
268 <         = |J_z|\Big|\frac{\partial^2 T}{\partial z^2}\Big|
269 <         \Big/\left(\frac{\partial T}{\partial z}\right)^2
270 < \label{derivativeG}
271 < \end{equation}
260 > Another approach is to assume that the temperature is continuous and
261 > differentiable throughout the space. Given that $\lambda$ is also
262 > differentiable, $G$ can be defined as its gradient ($\nabla\lambda$)
263 > projected along a vector normal to the interface ($\mathbf{\hat{n}}$)
264 > and evaluated at the interface location ($z_0$). This quantity,
265 > \begin{align}
266 > G^\prime &= \left(\nabla\lambda \cdot \mathbf{\hat{n}}\right)_{z_0} \\
267 >         &= \frac{\partial}{\partial z}\left(-\frac{J_z}{
268 >           \left(\frac{\partial T}{\partial z}\right)}\right)_{z_0} \\
269 >         &= J_z\left(\frac{\partial^2 T}{\partial z^2}\right)_{z_0} \Big/
270 >         \left(\frac{\partial T}{\partial z}\right)_{z_0}^2 \label{derivativeG}
271 > \end{align}
272 > has the same units as the common definition for $G$, and the maximum
273 > of its magnitude denotes where thermal conductivity has the largest
274 > change, i.e. the interface.  In the geometry used in this study, the
275 > vector normal to the interface points along the $z$ axis, as do
276 > $\vec{J}$ and the thermal gradient.  This yields the simplified
277 > expressions in Eq. \ref{derivativeG}.
278  
279   With temperature profiles obtained from simulation, one is able to
280   approximate the first and second derivatives of $T$ with finite
# Line 257 | Line 298 | profile.
298  
299   \begin{figure}
300   \includegraphics[width=\linewidth]{gradT}
301 < \caption{A sample of Au-butanethiol/hexane interfacial system and the
302 <  temperature profile after a kinetic energy flux is imposed to
303 <  it. The 1st and 2nd derivatives of the temperature profile can be
304 <  obtained with finite difference approximation (lower panel).}
301 > \caption{A sample of Au (111) / butanethiol / hexane interfacial
302 >  system with the temperature profile after a kinetic energy flux has
303 >  been imposed.  Note that the largest temperature jump in the thermal
304 >  profile (corresponding to the lowest interfacial conductance) is at
305 >  the interface between the butanethiol molecules (blue) and the
306 >  solvent (grey).  First and second derivatives of the temperature
307 >  profile are obtained using a finite difference approximation (lower
308 >  panel).}
309   \label{gradT}
310   \end{figure}
311  
# Line 307 | Line 352 | between periodic images of the gold interfaces is $45
352   solvent molecules would change the normal behavior of the liquid
353   phase. Therefore, our $N_{solvent}$ values were chosen to ensure that
354   these extreme cases did not happen to our simulations. The spacing
355 < between periodic images of the gold interfaces is $45 \sim 75$\AA.
355 > between periodic images of the gold interfaces is $45 \sim 75$\AA in
356 > our simulations.
357  
358   The initial configurations generated are further equilibrated with the
359   $x$ and $y$ dimensions fixed, only allowing the $z$-length scale to
# Line 325 | Line 371 | gradient had stablized, the temperature profile of the
371   $\sim$200K. Therefore, thermal flux usually came from the metal to the
372   liquid so that the liquid has a higher temperature and would not
373   freeze due to lowered temperatures. After this induced temperature
374 < gradient had stablized, the temperature profile of the simulation cell
375 < was recorded. To do this, the simulation cell is devided evenly into
374 > gradient had stabilized, the temperature profile of the simulation cell
375 > was recorded. To do this, the simulation cell is divided evenly into
376   $N$ slabs along the $z$-axis. The average temperatures of each slab
377   are recorded for 1$\sim$2 ns. When the slab width $d$ of each slab is
378   the same, the derivatives of $T$ with respect to slab number $n$ can
# Line 339 | Line 385 | be directly used for $G^\prime$ calculations: \begin{e
385           \Big/\left(\frac{\partial T}{\partial n}\right)^2
386   \label{derivativeG2}
387   \end{equation}
388 + The absolute values in Eq. \ref{derivativeG2} appear because the
389 + direction of the flux $\vec{J}$ is in an opposing direction on either
390 + side of the metal slab.
391  
392   All of the above simulation procedures use a time step of 1 fs. Each
393   equilibration stage took a minimum of 100 ps, although in some cases,
# Line 357 | Line 406 | particles of different species.
406   \caption{Structures of the capping agent and solvents utilized in
407    these simulations. The chemically-distinct sites (a-e) are expanded
408    in terms of constituent atoms for both United Atom (UA) and All Atom
409 <  (AA) force fields.  Most parameters are from
410 <  Refs. \protect\cite{TraPPE-UA.alkanes,TraPPE-UA.alkylbenzenes,TraPPE-UA.thiols}
409 >  (AA) force fields.  Most parameters are from References
410 >  \protect\cite{TraPPE-UA.alkanes,TraPPE-UA.alkylbenzenes,TraPPE-UA.thiols}
411    (UA) and \protect\cite{OPLSAA} (AA). Cross-interactions with the Au
412    atoms are given in Table \ref{MnM}.}
413   \label{demoMol}
# Line 382 | Line 431 | However, the TraPPE-UA model for alkanes is known to p
431  
432   By eliminating explicit hydrogen atoms, the TraPPE-UA models are
433   simple and computationally efficient, while maintaining good accuracy.
434 < However, the TraPPE-UA model for alkanes is known to predict a slighly
434 > However, the TraPPE-UA model for alkanes is known to predict a slightly
435   lower boiling point than experimental values. This is one of the
436   reasons we used a lower average temperature (200K) for our
437   simulations. If heat is transferred to the liquid phase during the
# Line 501 | Line 550 | with respect to surface coverage.
550  
551   \begin{figure}
552   \includegraphics[width=\linewidth]{coverage}
553 < \caption{Comparison of interfacial thermal conductivity ($G$) values
554 <  for the Au-butanethiol/solvent interface with various UA models and
555 <  different capping agent coverages at $\langle T\rangle\sim$200K.}
553 > \caption{The interfacial thermal conductivity ($G$) has a
554 >  non-monotonic dependence on the degree of surface capping.  This
555 >  data is for the Au(111) / butanethiol / solvent interface with
556 >  various UA force fields at $\langle T\rangle \sim $200K.}
557   \label{coverage}
558   \end{figure}
559  
560 < In partially covered surfaces, the derivative definition for $G^\prime$ (Eq. \ref{derivativeG}) becomes difficult to apply, as the location of maximum change of $\lambda$ becomes washed out.  The discrete definition (Eq. \ref{discreteG}) is easier to apply, as the Gibbs dividing surface is still well-defined. Therefore, $G$ (not $G^\prime$) was used in this section.
560 > In partially covered surfaces, the derivative definition for
561 > $G^\prime$ (Eq. \ref{derivativeG}) becomes difficult to apply, as the
562 > location of maximum change of $\lambda$ becomes washed out.  The
563 > discrete definition (Eq. \ref{discreteG}) is easier to apply, as the
564 > Gibbs dividing surface is still well-defined. Therefore, $G$ (not
565 > $G^\prime$) was used in this section.
566  
567 < From Figure \ref{coverage}, one can see the significance of the presence of capping agents. When even a small fraction of the Au(111) surface sites are covered with butanethiols, the conductivity exhibits an enhancement by at least a factor of 3.  Cappping agents are clearly playing a major role in thermal transport at metal / organic solvent surfaces.
568 <
569 < We note a non-monotonic behavior in the interfacial conductance as a function of surface coverage. The maximum conductance (largest $G$) happens when the surfaces are about 75\% covered with butanethiol caps.  The reason for this behavior is not entirely clear.  One explanation is that incomplete butanethiol coverage allows small gaps between butanethiols to form. These gaps can be filled by transient solvent molecules.  These solvent molecules couple very strongly with the hot capping agent molecules near the surface, and can then carry away (diffusively) the excess thermal energy from the surface.
570 <
571 < There appears to be a competition between the conduction of the thermal energy away from the surface by the capping agents (enhanced by greater coverage) and the coupling of the capping agents with the solvent (enhanced by interdigitation at lower coverages).  This competition would lead to the non-monotonic coverage behavior observed here.
567 > From Figure \ref{coverage}, one can see the significance of the
568 > presence of capping agents. When even a small fraction of the Au(111)
569 > surface sites are covered with butanethiols, the conductivity exhibits
570 > an enhancement by at least a factor of 3.  Capping agents are clearly
571 > playing a major role in thermal transport at metal / organic solvent
572 > surfaces.
573  
574 < Results for rigid body toluene solvent, as well as the UA hexane, are within the ranges expected from prior experimental work.\cite{Wilson:2002uq,cahill:793,PhysRevB.80.195406} This suggests that explicit hydrogen atoms might not be  required for modeling thermal transport in these systems.  C-H vibrational modes do not see significant excited state population at low temperatures, and are not likely to carry lower frequency excitations from the solid layer into the bulk liquid.
574 > We note a non-monotonic behavior in the interfacial conductance as a
575 > function of surface coverage. The maximum conductance (largest $G$)
576 > happens when the surfaces are about 75\% covered with butanethiol
577 > caps.  The reason for this behavior is not entirely clear.  One
578 > explanation is that incomplete butanethiol coverage allows small gaps
579 > between butanethiols to form. These gaps can be filled by transient
580 > solvent molecules.  These solvent molecules couple very strongly with
581 > the hot capping agent molecules near the surface, and can then carry
582 > away (diffusively) the excess thermal energy from the surface.
583  
584 < The toluene solvent does not exhibit the same behavior as hexane in that $G$ remains at approximately the same magnitude when the capping coverage increases from 25\% to 75\%.  Toluene, as a rigid planar molecule, cannot occupy the relatively small gaps between the capping agents as easily as the chain-like {\it n}-hexane.   The effect of solvent coupling to the capping agent is therefore weaker in toluene except at the very lowest coverage levels.  This effect counters the coverage-dependent conduction of heat away from the metal surface, leading to a much flatter $G$ vs. coverage trend than is observed in {\it n}-hexane.
584 > There appears to be a competition between the conduction of the
585 > thermal energy away from the surface by the capping agents (enhanced
586 > by greater coverage) and the coupling of the capping agents with the
587 > solvent (enhanced by interdigitation at lower coverages).  This
588 > competition would lead to the non-monotonic coverage behavior observed
589 > here.
590  
591 + Results for rigid body toluene solvent, as well as the UA hexane, are
592 + within the ranges expected from prior experimental
593 + work.\cite{Wilson:2002uq,cahill:793,PhysRevB.80.195406} This suggests
594 + that explicit hydrogen atoms might not be required for modeling
595 + thermal transport in these systems.  C-H vibrational modes do not see
596 + significant excited state population at low temperatures, and are not
597 + likely to carry lower frequency excitations from the solid layer into
598 + the bulk liquid.
599 +
600 + The toluene solvent does not exhibit the same behavior as hexane in
601 + that $G$ remains at approximately the same magnitude when the capping
602 + coverage increases from 25\% to 75\%.  Toluene, as a rigid planar
603 + molecule, cannot occupy the relatively small gaps between the capping
604 + agents as easily as the chain-like {\it n}-hexane.  The effect of
605 + solvent coupling to the capping agent is therefore weaker in toluene
606 + except at the very lowest coverage levels.  This effect counters the
607 + coverage-dependent conduction of heat away from the metal surface,
608 + leading to a much flatter $G$ vs. coverage trend than is observed in
609 + {\it n}-hexane.
610 +
611   \subsection{Effects due to Solvent \& Solvent Models}
612 < In addition to UA solvent and capping agent models, AA models have also been included in our simulations.  In most of this work, the same (UA or AA) model for solvent and capping agent was used, but it is also possible to utilize different models for different components.  We have also included isotopic substitutions (Hydrogen to Deuterium) to decrease the explicit vibrational overlap between solvent and capping agent. Table \ref{modelTest} summarizes the results of these studies.
612 > In addition to UA solvent and capping agent models, AA models have
613 > also been included in our simulations.  In most of this work, the same
614 > (UA or AA) model for solvent and capping agent was used, but it is
615 > also possible to utilize different models for different components.
616 > We have also included isotopic substitutions (Hydrogen to Deuterium)
617 > to decrease the explicit vibrational overlap between solvent and
618 > capping agent. Table \ref{modelTest} summarizes the results of these
619 > studies.
620  
621   \begin{table*}
622    \begin{minipage}{\linewidth}
# Line 529 | Line 625 | In addition to UA solvent and capping agent models, AA
625        \caption{Computed interfacial thermal conductance ($G$ and
626          $G^\prime$) values for interfaces using various models for
627          solvent and capping agent (or without capping agent) at
628 <        $\langle T\rangle\sim$200K. (D stands for deuterated solvent
629 <        or capping agent molecules; ``Avg.'' denotes results that are
630 <        averages of simulations under different applied thermal flux values $(J_z)$. Error
535 <        estimates are indicated in parentheses.)}
628 >        $\langle T\rangle\sim$200K.  Here ``D'' stands for deuterated
629 >        solvent or capping agent molecules. Error estimates are
630 >        indicated in parentheses.}
631        
632        \begin{tabular}{llccc}
633          \hline\hline
634 <        Butanethiol model & Solvent & $J_z$ & $G$ & $G^\prime$ \\
635 <        (or bare surface) & model & (GW/m$^2$) &
634 >        Butanethiol model & Solvent & $G$ & $G^\prime$ \\
635 >        (or bare surface) & model &
636          \multicolumn{2}{c}{(MW/m$^2$/K)} \\
637          \hline
638 <        UA    & UA hexane    & Avg. & 131(9)    & 87(10)    \\
639 <              & UA hexane(D) & 1.95 & 153(5)    & 136(13)   \\
640 <              & AA hexane    & Avg. & 131(6)    & 122(10)   \\
641 <              & UA toluene   & 1.96 & 187(16)   & 151(11)   \\
642 <              & AA toluene   & 1.89 & 200(36)   & 149(53)   \\
638 >        UA    & UA hexane    & 131(9)    & 87(10)    \\
639 >              & UA hexane(D) & 153(5)    & 136(13)   \\
640 >              & AA hexane    & 131(6)    & 122(10)   \\
641 >              & UA toluene   & 187(16)   & 151(11)   \\
642 >              & AA toluene   & 200(36)   & 149(53)   \\
643          \hline
644 <        AA    & UA hexane    & 1.94 & 116(9)    & 129(8)    \\
645 <              & AA hexane    & Avg. & 442(14)   & 356(31)   \\
646 <              & AA hexane(D) & 1.93 & 222(12)   & 234(54)   \\
647 <              & UA toluene   & 1.98 & 125(25)   & 97(60)    \\
648 <              & AA toluene   & 3.79 & 487(56)   & 290(42)   \\
644 >        AA    & UA hexane    & 116(9)    & 129(8)    \\
645 >              & AA hexane    & 442(14)   & 356(31)   \\
646 >              & AA hexane(D) & 222(12)   & 234(54)   \\
647 >              & UA toluene   & 125(25)   & 97(60)    \\
648 >              & AA toluene   & 487(56)   & 290(42)   \\
649          \hline
650 <        AA(D) & UA hexane    & 1.94 & 158(25)   & 172(4)    \\
651 <              & AA hexane    & 1.92 & 243(29)   & 191(11)   \\
652 <              & AA toluene   & 1.93 & 364(36)   & 322(67)   \\
650 >        AA(D) & UA hexane    & 158(25)   & 172(4)    \\
651 >              & AA hexane    & 243(29)   & 191(11)   \\
652 >              & AA toluene   & 364(36)   & 322(67)   \\
653          \hline
654 <        bare  & UA hexane    & Avg. & 46.5(3.2) & 49.4(4.5) \\
655 <              & UA hexane(D) & 0.98 & 43.9(4.6) & 43.0(2.0) \\
656 <              & AA hexane    & 0.96 & 31.0(1.4) & 29.4(1.3) \\
657 <              & UA toluene   & 1.99 & 70.1(1.3) & 65.8(0.5) \\
654 >        bare  & UA hexane    & 46.5(3.2) & 49.4(4.5) \\
655 >              & UA hexane(D) & 43.9(4.6) & 43.0(2.0) \\
656 >              & AA hexane    & 31.0(1.4) & 29.4(1.3) \\
657 >              & UA toluene   & 70.1(1.3) & 65.8(0.5) \\
658          \hline\hline
659        \end{tabular}
660        \label{modelTest}
# Line 567 | Line 662 | To facilitate direct comparison between force fields,
662    \end{minipage}
663   \end{table*}
664  
665 < To facilitate direct comparison between force fields, systems with the same capping agent and solvent were prepared with the same length scales for the simulation cells.  
665 > To facilitate direct comparison between force fields, systems with the
666 > same capping agent and solvent were prepared with the same length
667 > scales for the simulation cells.
668  
669 < On bare metal / solvent surfaces, different force field models for hexane yield similar results for both $G$ and $G^\prime$, and these two definitions agree with each other very well. This is primarily an indicator of weak interactions between the metal and the solvent, and is a typical case for acoustic impedance mismatch between these two phases.
669 > On bare metal / solvent surfaces, different force field models for
670 > hexane yield similar results for both $G$ and $G^\prime$, and these
671 > two definitions agree with each other very well. This is primarily an
672 > indicator of weak interactions between the metal and the solvent.
673  
674 < For the fully-covered surfaces, the choice of force field for the capping agent and solvent has a large impact on the calulated values of $G$ and $G^\prime$.  The AA thiol to AA solvent conductivities are much larger than their UA to UA counterparts, and these values exceed the experimental estimates by a large measure.  The AA force field allows significant energy to go into C-H (or C-D) stretching modes, and since these modes are high frequency, this non-quantum behavior is likely responsible for the overestimate of the conductivity.
674 > For the fully-covered surfaces, the choice of force field for the
675 > capping agent and solvent has a large impact on the calculated values
676 > of $G$ and $G^\prime$.  The AA thiol to AA solvent conductivities are
677 > much larger than their UA to UA counterparts, and these values exceed
678 > the experimental estimates by a large measure.  The AA force field
679 > allows significant energy to go into C-H (or C-D) stretching modes,
680 > and since these modes are high frequency, this non-quantum behavior is
681 > likely responsible for the overestimate of the conductivity.  Compared
682 > to the AA model, the UA model yields more reasonable conductivity
683 > values with much higher computational efficiency.
684  
576 The similarity in the vibrational modes available to solvent and capping agent can be reduced by deuterating one of the two components.  Once either the hexanes or the butanethiols are deuterated, one can see a significantly lower $G$ and $G^\prime$ (Figure \ref{aahxntln}).  Compared to the AA model, the UA model yields more reasonable conductivity values with much higher computational efficiency.
577
578 \begin{figure}
579 \includegraphics[width=\linewidth]{aahxntln}
580 \caption{Spectra obtained for all-atom (AA) Au / butanethiol / solvent
581  systems. When butanethiol is deuterated (lower left), its
582  vibrational overlap with hexane decreases significantly.  Since aromatic molecules and the butanethiol are vibrationally dissimilar, the change is not as dramatic when toluene is the solvent (right).}
583 \label{aahxntln}
584 \end{figure}
585
586 For the Au / butanethiol / toluene interfaces, having the AA butanethiol deuterated did not yield a significant change in the measured conductance. Compared to the C-H vibrational overlap between hexane and butanethiol, both of which have alkyl chains, the overlap between toluene and butanethiol is not as significant and thus does not contribute as much to the heat exchange process.  The presence of extra degrees of freedom in the AA force field for toluene yields higher heat exchange rates between the two phases and results in a much higher conductivity than in the UA force field.
587
685   \subsubsection{Are electronic excitations in the metal important?}
686 < Because they lack electronic excitations, the QSC and related embedded atom method (EAM) models for gold are known to predict unreasonably low values for bulk conductivity ($\lambda$).\cite{kuang:164101,ISI:000207079300006,Clancy:1992} If the conductance between the phases ($G$) is governed primarily by phonon excitation (and not electronic degrees of freedom), one would expect a classical model to capture most of the interfacial thermal conductance.  Our results for $G$ and $G^\prime$ indicate that this is indeed the case, and suggest that the modeling of interfacial thermal transport depends primarily on the description of the interactions between the various components at the interface.  When the metal is chemically capped, the primary barrier to thermal conductivity appears to be the interface between the capping agent and the surrounding solvent, so the excitations in the metal have little impact on the value of $G$.
686 > Because they lack electronic excitations, the QSC and related embedded
687 > atom method (EAM) models for gold are known to predict unreasonably
688 > low values for bulk conductivity
689 > ($\lambda$).\cite{kuang:164101,ISI:000207079300006,Clancy:1992} If the
690 > conductance between the phases ($G$) is governed primarily by phonon
691 > excitation (and not electronic degrees of freedom), one would expect a
692 > classical model to capture most of the interfacial thermal
693 > conductance.  Our results for $G$ and $G^\prime$ indicate that this is
694 > indeed the case, and suggest that the modeling of interfacial thermal
695 > transport depends primarily on the description of the interactions
696 > between the various components at the interface.  When the metal is
697 > chemically capped, the primary barrier to thermal conductivity appears
698 > to be the interface between the capping agent and the surrounding
699 > solvent, so the excitations in the metal have little impact on the
700 > value of $G$.
701  
702   \subsection{Effects due to methodology and simulation parameters}
703  
704 < START HERE
704 > We have varied the parameters of the simulations in order to
705 > investigate how these factors would affect the computation of $G$.  Of
706 > particular interest are: 1) the length scale for the applied thermal
707 > gradient (modified by increasing the amount of solvent in the system),
708 > 2) the sign and magnitude of the applied thermal flux, 3) the average
709 > temperature of the simulation (which alters the solvent density during
710 > equilibration), and 4) the definition of the interfacial conductance
711 > (Eqs. (\ref{discreteG}) or (\ref{derivativeG})) used in the
712 > calculation.
713  
714 < We have varied our protocol or other parameters of the simulations in order to investigate how these factors would affect the computation of $G$.
714 > Systems of different lengths were prepared by altering the number of
715 > solvent molecules and extending the length of the box along the $z$
716 > axis to accomodate the extra solvent.  Equilibration at the same
717 > temperature and pressure conditions led to nearly identical surface
718 > areas ($L_x$ and $L_y$) available to the metal and capping agent,
719 > while the extra solvent served mainly to lengthen the axis that was
720 > used to apply the thermal flux.  For a given value of the applied
721 > flux, the different $z$ length scale has only a weak effect on the
722 > computed conductivities (Table \ref{AuThiolHexaneUA}).
723  
597 We allowed $L_x$ and $L_y$ to change during equilibrating the liquid phase. Due to the stiffness of the crystalline Au structure, $L_x$ and $L_y$ would not change noticeably after equilibration. Although $L_z$ could fluctuate $\sim$1\% after a system is fully equilibrated in the NPT ensemble, this fluctuation, as well as those of $L_x$ and $L_y$ (which is significantly smaller), would not be magnified on the calculated $G$'s, as shown in Table \ref{AuThiolHexaneUA}. This insensivity to $L_i$ fluctuations allows reliable measurement of $G$'s without the necessity of extremely cautious equilibration process.
598
599 As stated in our computational details, the spacing filled with
600 solvent molecules can be chosen within a range. This allows some
601 change of solvent molecule numbers for the same Au-butanethiol
602 surfaces. We did this study on our Au-butanethiol/hexane
603 simulations. Nevertheless, the results obtained from systems of
604 different $N_{hexane}$ did not indicate that the measurement of $G$ is
605 susceptible to this parameter. For computational efficiency concern,
606 smaller system size would be preferable, given that the liquid phase
607 structure is not affected.
608
724   \subsubsection{Effects of applied flux}
725 < Our NIVS algorithm allows change of unphysical thermal flux both in
726 < direction and in quantity. This feature extends our investigation of
727 < interfacial thermal conductance. However, the magnitude of this
728 < thermal flux is not arbitary if one aims to obtain a stable and
729 < reliable thermal gradient. A temperature profile would be
730 < substantially affected by noise when $|J_z|$ has a much too low
731 < magnitude; while an excessively large $|J_z|$ that overwhelms the
732 < conductance capacity of the interface would prevent a thermal gradient
733 < to reach a stablized steady state. NIVS has the advantage of allowing
734 < $J$ to vary in a wide range such that the optimal flux range for $G$
735 < measurement can generally be simulated by the algorithm. Within the
736 < optimal range, we were able to study how $G$ would change according to
622 < the thermal flux across the interface. For our simulations, we denote
623 < $J_z$ to be positive when the physical thermal flux is from the liquid
624 < to metal, and negative vice versa. The $G$'s measured under different
625 < $J_z$ is listed in Table \ref{AuThiolHexaneUA} and
626 < \ref{AuThiolToluene}. These results do not suggest that $G$ is
627 < dependent on $J_z$ within this flux range. The linear response of flux
628 < to thermal gradient simplifies our investigations in that we can rely
629 < on $G$ measurement with only a couple $J_z$'s and do not need to test
630 < a large series of fluxes.
725 > The NIVS algorithm allows changes in both the sign and magnitude of
726 > the applied flux.  It is possible to reverse the direction of heat
727 > flow simply by changing the sign of the flux, and thermal gradients
728 > which would be difficult to obtain experimentally ($5$ K/\AA) can be
729 > easily simulated.  However, the magnitude of the applied flux is not
730 > arbitrary if one aims to obtain a stable and reliable thermal gradient.
731 > A temperature gradient can be lost in the noise if $|J_z|$ is too
732 > small, and excessive $|J_z|$ values can cause phase transitions if the
733 > extremes of the simulation cell become widely separated in
734 > temperature.  Also, if $|J_z|$ is too large for the bulk conductivity
735 > of the materials, the thermal gradient will never reach a stable
736 > state.  
737  
738 < \begin{table*}
739 <  \begin{minipage}{\linewidth}
740 <    \begin{center}
741 <      \caption{Computed interfacial thermal conductivity ($G$ and
742 <        $G^\prime$) values for the 100\% covered Au-butanethiol/hexane
743 <        interfaces with UA model and different hexane molecule numbers
744 <        at different temperatures using a range of energy
745 <        fluxes. Error estimates indicated in parenthesis.}
746 <      
747 <      \begin{tabular}{ccccccc}
748 <        \hline\hline
749 <        $\langle T\rangle$ & $N_{hexane}$ & Fixed & $\rho_{hexane}$ &
644 <        $J_z$ & $G$ & $G^\prime$ \\
645 <        (K) & & $L_x$ \& $L_y$? & (g/cm$^3$) & (GW/m$^2$) &
646 <        \multicolumn{2}{c}{(MW/m$^2$/K)} \\
647 <        \hline
648 <        200 & 266 & No  & 0.672 & -0.96 & 102(3)    & 80.0(0.8) \\
649 <            & 200 & Yes & 0.694 &  1.92 & 129(11)   & 87.3(0.3) \\
650 <            &     & Yes & 0.672 &  1.93 & 131(16)   & 78(13)    \\
651 <            &     & No  & 0.688 &  0.96 & 125(16)   & 90.2(15)  \\
652 <            &     &     &       &  1.91 & 139(10)   & 101(10)   \\
653 <            &     &     &       &  2.83 & 141(6)    & 89.9(9.8) \\
654 <            & 166 & Yes & 0.679 &  0.97 & 115(19)   & 69(18)    \\
655 <            &     &     &       &  1.94 & 125(9)    & 87.1(0.2) \\
656 <            &     & No  & 0.681 &  0.97 & 141(30)   & 78(22)    \\
657 <            &     &     &       &  1.92 & 138(4)    & 98.9(9.5) \\
658 <        \hline
659 <        250 & 200 & No  & 0.560 &  0.96 & 75(10)    & 61.8(7.3) \\
660 <            &     &     &       & -0.95 & 49.4(0.3) & 45.7(2.1) \\
661 <            & 166 & Yes & 0.570 &  0.98 & 79.0(3.5) & 62.9(3.0) \\
662 <            &     & No  & 0.569 &  0.97 & 80.3(0.6) & 67(11)    \\
663 <            &     &     &       &  1.44 & 76.2(5.0) & 64.8(3.8) \\
664 <            &     &     &       & -0.95 & 56.4(2.5) & 54.4(1.1) \\
665 <            &     &     &       & -1.85 & 47.8(1.1) & 53.5(1.5) \\
666 <        \hline\hline
667 <      \end{tabular}
668 <      \label{AuThiolHexaneUA}
669 <    \end{center}
670 <  \end{minipage}
671 < \end{table*}
738 > Within a reasonable range of $J_z$ values, we were able to study how
739 > $G$ changes as a function of this flux.  In what follows, we use
740 > positive $J_z$ values to denote the case where energy is being
741 > transferred by the method from the metal phase and into the liquid.
742 > The resulting gradient therefore has a higher temperature in the
743 > liquid phase.  Negative flux values reverse this transfer, and result
744 > in higher temperature metal phases.  The conductance measured under
745 > different applied $J_z$ values is listed in Tables 1 and 2 in the
746 > supporting information. These results do not indicate that $G$ depends
747 > strongly on $J_z$ within this flux range. The linear response of flux
748 > to thermal gradient simplifies our investigations in that we can rely
749 > on $G$ measurement with only a small number $J_z$ values.
750  
751 + The sign of $J_z$ is a different matter, however, as this can alter
752 + the temperature on the two sides of the interface. The average
753 + temperature values reported are for the entire system, and not for the
754 + liquid phase, so at a given $\langle T \rangle$, the system with
755 + positive $J_z$ has a warmer liquid phase.  This means that if the
756 + liquid carries thermal energy via diffusive transport, {\it positive}
757 + $J_z$ values will result in increased molecular motion on the liquid
758 + side of the interface, and this will increase the measured
759 + conductivity.
760 +
761   \subsubsection{Effects due to average temperature}
762  
763 < Furthermore, we also attempted to increase system average temperatures
764 < to above 200K. These simulations are first equilibrated in the NPT
765 < ensemble under normal pressure. As stated above, the TraPPE-UA model
766 < for hexane tends to predict a lower boiling point. In our simulations,
767 < hexane had diffculty to remain in liquid phase when NPT equilibration
768 < temperature is higher than 250K. Additionally, the equilibrated liquid
769 < hexane density under 250K becomes lower than experimental value. This
770 < expanded liquid phase leads to lower contact between hexane and
771 < butanethiol as well.[MAY NEED SLAB DENSITY FIGURE]
772 < And this reduced contact would
685 < probably be accountable for a lower interfacial thermal conductance,
686 < as shown in Table \ref{AuThiolHexaneUA}.
763 > We also studied the effect of average system temperature on the
764 > interfacial conductance.  The simulations are first equilibrated in
765 > the NPT ensemble at 1 atm.  The TraPPE-UA model for hexane tends to
766 > predict a lower boiling point (and liquid state density) than
767 > experiments.  This lower-density liquid phase leads to reduced contact
768 > between the hexane and butanethiol, and this accounts for our
769 > observation of lower conductance at higher temperatures.  In raising
770 > the average temperature from 200K to 250K, the density drop of
771 > $\sim$20\% in the solvent phase leads to a $\sim$40\% drop in the
772 > conductance.
773  
774 < A similar study for TraPPE-UA toluene agrees with the above result as
775 < well. Having a higher boiling point, toluene tends to remain liquid in
776 < our simulations even equilibrated under 300K in NPT
777 < ensembles. Furthermore, the expansion of the toluene liquid phase is
778 < not as significant as that of the hexane. This prevents severe
693 < decrease of liquid-capping agent contact and the results (Table
694 < \ref{AuThiolToluene}) show only a slightly decreased interface
695 < conductance. Therefore, solvent-capping agent contact should play an
696 < important role in the thermal transport process across the interface
697 < in that higher degree of contact could yield increased conductance.
774 > Similar behavior is observed in the TraPPE-UA model for toluene,
775 > although this model has better agreement with the experimental
776 > densities of toluene.  The expansion of the toluene liquid phase is
777 > not as significant as that of the hexane (8.3\% over 100K), and this
778 > limits the effect to $\sim$20\% drop in thermal conductivity.
779  
780 < \begin{table*}
781 <  \begin{minipage}{\linewidth}
782 <    \begin{center}
783 <      \caption{Computed interfacial thermal conductivity ($G$ and
784 <        $G^\prime$) values for a 90\% coverage Au-butanethiol/toluene
704 <        interface at different temperatures using a range of energy
705 <        fluxes. Error estimates indicated in parenthesis.}
706 <      
707 <      \begin{tabular}{ccccc}
708 <        \hline\hline
709 <        $\langle T\rangle$ & $\rho_{toluene}$ & $J_z$ & $G$ & $G^\prime$ \\
710 <        (K) & (g/cm$^3$) & (GW/m$^2$) & \multicolumn{2}{c}{(MW/m$^2$/K)} \\
711 <        \hline
712 <        200 & 0.933 &  2.15 & 204(12) & 113(12) \\
713 <            &       & -1.86 & 180(3)  & 135(21) \\
714 <            &       & -3.93 & 176(5)  & 113(12) \\
715 <        \hline
716 <        300 & 0.855 & -1.91 & 143(5)  & 125(2)  \\
717 <            &       & -4.19 & 135(9)  & 113(12) \\
718 <        \hline\hline
719 <      \end{tabular}
720 <      \label{AuThiolToluene}
721 <    \end{center}
722 <  \end{minipage}
723 < \end{table*}
780 > Although we have not mapped out the behavior at a large number of
781 > temperatures, is clear that there will be a strong temperature
782 > dependence in the interfacial conductance when the physical properties
783 > of one side of the interface (notably the density) change rapidly as a
784 > function of temperature.
785  
786 < Besides lower interfacial thermal conductance, surfaces in relatively
787 < high temperatures are susceptible to reconstructions, when
788 < butanethiols have a full coverage on the Au(111) surface. These
789 < reconstructions include surface Au atoms migrated outward to the S
790 < atom layer, and butanethiol molecules embedded into the original
791 < surface Au layer. The driving force for this behavior is the strong
792 < Au-S interactions in our simulations. And these reconstructions lead
793 < to higher ratio of Au-S attraction and thus is energetically
794 < favorable. Furthermore, this phenomenon agrees with experimental
795 < results\cite{doi:10.1021/j100035a033,doi:10.1021/la026493y}. Vlugt
796 < {\it et al.} had kept their Au(111) slab rigid so that their
797 < simulations can reach 300K without surface reconstructions. Without
798 < this practice, simulating 100\% thiol covered interfaces under higher
799 < temperatures could hardly avoid surface reconstructions. However, our
800 < measurement is based on assuming homogeneity on $x$ and $y$ dimensions
801 < so that measurement of $T$ at particular $z$ would be an effective
802 < average of the particles of the same type. Since surface
742 < reconstructions could eliminate the original $x$ and $y$ dimensional
743 < homogeneity, measurement of $G$ is more difficult to conduct under
744 < higher temperatures. Therefore, most of our measurements are
745 < undertaken at $\langle T\rangle\sim$200K.
786 > Besides the lower interfacial thermal conductance, surfaces at
787 > relatively high temperatures are susceptible to reconstructions,
788 > particularly when butanethiols fully cover the Au(111) surface. These
789 > reconstructions include surface Au atoms which migrate outward to the
790 > S atom layer, and butanethiol molecules which embed into the surface
791 > Au layer. The driving force for this behavior is the strong Au-S
792 > interactions which are modeled here with a deep Lennard-Jones
793 > potential. This phenomenon agrees with reconstructions that have been
794 > experimentally
795 > observed.\cite{doi:10.1021/j100035a033,doi:10.1021/la026493y}.  Vlugt
796 > {\it et al.} kept their Au(111) slab rigid so that their simulations
797 > could reach 300K without surface
798 > reconstructions.\cite{vlugt:cpc2007154} Since surface reconstructions
799 > blur the interface, the measurement of $G$ becomes more difficult to
800 > conduct at higher temperatures.  For this reason, most of our
801 > measurements are undertaken at $\langle T\rangle\sim$200K where
802 > reconstruction is minimized.
803  
804   However, when the surface is not completely covered by butanethiols,
805 < the simulated system is more resistent to the reconstruction
806 < above. Our Au-butanethiol/toluene system had the Au(111) surfaces 90\%
807 < covered by butanethiols, but did not see this above phenomena even at
808 < $\langle T\rangle\sim$300K. The empty three-fold sites not occupied by
809 < capping agents could help prevent surface reconstruction in that they
810 < provide other means of capping agent relaxation. It is observed that
811 < butanethiols can migrate to their neighbor empty sites during a
812 < simulation. Therefore, we were able to obtain $G$'s for these
756 < interfaces even at a relatively high temperature without being
757 < affected by surface reconstructions.
805 > the simulated system appears to be more resistent to the
806 > reconstruction. Our Au / butanethiol / toluene system had the Au(111)
807 > surfaces 90\% covered by butanethiols, but did not see this above
808 > phenomena even at $\langle T\rangle\sim$300K.  That said, we did
809 > observe butanethiols migrating to neighboring three-fold sites during
810 > a simulation.  Since the interface persisted in these simulations, we
811 > were able to obtain $G$'s for these interfaces even at a relatively
812 > high temperature without being affected by surface reconstructions.
813  
759
814   \section{Discussion}
815  
816 < \subsection{Capping agent acts as a vibrational coupler between solid
817 <  and solvent phases}
816 > The primary result of this work is that the capping agent acts as an
817 > efficient thermal coupler between solid and solvent phases.  One of
818 > the ways the capping agent can carry out this role is to down-shift
819 > between the phonon vibrations in the solid (which carry the heat from
820 > the gold) and the molecular vibrations in the liquid (which carry some
821 > of the heat in the solvent).
822 >
823   To investigate the mechanism of interfacial thermal conductance, the
824   vibrational power spectrum was computed. Power spectra were taken for
825   individual components in different simulations. To obtain these
826 < spectra, simulations were run after equilibration, in the NVE
827 < ensemble, and without a thermal gradient. Snapshots of configurations
828 < were collected at a frequency that is higher than that of the fastest
829 < vibrations occuring in the simulations. With these configurations, the
830 < velocity auto-correlation functions can be computed:
826 > spectra, simulations were run after equilibration in the
827 > microcanonical (NVE) ensemble and without a thermal
828 > gradient. Snapshots of configurations were collected at a frequency
829 > that is higher than that of the fastest vibrations occurring in the
830 > simulations. With these configurations, the velocity auto-correlation
831 > functions can be computed:
832   \begin{equation}
833   C_A (t) = \langle\vec{v}_A (t)\cdot\vec{v}_A (0)\rangle
834   \label{vCorr}
# Line 781 | Line 841 | symmetrized velocity autocorrelation function,
841   \label{fourier}
842   \end{equation}
843  
844 <
785 < \subsubsection{The role of specific vibrations}
844 > \subsection{The role of specific vibrations}
845   The vibrational spectra for gold slabs in different environments are
846   shown as in Figure \ref{specAu}. Regardless of the presence of
847 < solvent, the gold surfaces covered by butanethiol molecules, compared
848 < to bare gold surfaces, exhibit an additional peak observed at the
849 < frequency of $\sim$170cm$^{-1}$, which is attributed to the S-Au
850 < bonding vibration. This vibration enables efficient thermal transport
851 < from surface Au layer to the capping agents. Therefore, in our
852 < simulations, the Au/S interfaces do not appear major heat barriers
853 < compared to the butanethiol / solvent interfaces.
847 > solvent, the gold surfaces which are covered by butanethiol molecules
848 > exhibit an additional peak observed at a frequency of
849 > $\sim$165cm$^{-1}$.  We attribute this peak to the S-Au bonding
850 > vibration. This vibration enables efficient thermal coupling of the
851 > surface Au layer to the capping agents. Therefore, in our simulations,
852 > the Au / S interfaces do not appear to be the primary barrier to
853 > thermal transport when compared with the butanethiol / solvent
854 > interfaces.  This supports the results of Luo {\it et
855 >  al.}\cite{Luo20101}, who reported $G$ for Au-SAM junctions roughly
856 > twice as large as what we have computed for the thiol-liquid
857 > interfaces.
858  
796 \subsubsection{Overlap of power spectrum}
797 Simultaneously, the vibrational overlap between butanethiol and
798 organic solvents suggests higher thermal exchange efficiency between
799 these two components. Even exessively high heat transport was observed
800 when All-Atom models were used and C-H vibrations were treated
801 classically. Compared to metal and organic liquid phase, the heat
802 transfer efficiency between butanethiol and organic solvents is closer
803 to that within bulk liquid phase.
804
805 Furthermore, our observation validated previous
806 results\cite{hase:2010} that the intramolecular heat transport of
807 alkylthiols is highly effecient. As a combinational effects of these
808 phenomena, butanethiol acts as a channel to expedite thermal transport
809 process. The acoustic impedance mismatch between the metal and the
810 liquid phase can be effectively reduced with the presence of suitable
811 capping agents.
812
859   \begin{figure}
860   \includegraphics[width=\linewidth]{vibration}
861 < \caption{Vibrational spectra obtained for gold in different
862 <  environments.}
861 > \caption{The vibrational power spectrum for thiol-capped gold has an
862 >  additional vibrational peak at $\sim $165cm$^{-1}$.  Bare gold
863 >  surfaces (both with and without a solvent over-layer) are missing
864 >  this peak.   A similar peak at  $\sim $165cm$^{-1}$ also appears in
865 >  the vibrational power spectrum for the butanethiol capping agents.}
866   \label{specAu}
867   \end{figure}
868  
869 < \subsubsection{Isotopic substitution and vibrational overlap}
869 > Also in this figure, we show the vibrational power spectrum for the
870 > bound butanethiol molecules, which also exhibits the same
871 > $\sim$165cm$^{-1}$ peak.
872 >
873 > \subsection{Overlap of power spectra}
874   A comparison of the results obtained from the two different organic
875   solvents can also provide useful information of the interfacial
876 < thermal transport process. The deuterated hexane (UA) results do not
877 < appear to be substantially different from those of normal hexane (UA),
878 < given that butanethiol (UA) is non-deuterated for both solvents. The
879 < UA models, even though they have eliminated C-H vibrational overlap,
880 < still have significant overlap in the infrared spectra.  Because
881 < differences in the infrared range do not seem to produce an observable
882 < difference for the results of $G$ (Figure \ref{uahxnua}).
876 > thermal transport process.  In particular, the vibrational overlap
877 > between the butanethiol and the organic solvents suggests a highly
878 > efficient thermal exchange between these components.  Very high
879 > thermal conductivity was observed when AA models were used and C-H
880 > vibrations were treated classically. The presence of extra degrees of
881 > freedom in the AA force field yields higher heat exchange rates
882 > between the two phases and results in a much higher conductivity than
883 > in the UA force field. The all-atom classical models include high
884 > frequency modes which should be unpopulated at our relatively low
885 > temperatures.  This artifact is likely the cause of the high thermal
886 > conductance in all-atom MD simulations.
887  
888 + The similarity in the vibrational modes available to solvent and
889 + capping agent can be reduced by deuterating one of the two components
890 + (Fig. \ref{aahxntln}).  Once either the hexanes or the butanethiols
891 + are deuterated, one can observe a significantly lower $G$ and
892 + $G^\prime$ values (Table \ref{modelTest}).
893 +
894   \begin{figure}
895 + \includegraphics[width=\linewidth]{aahxntln}
896 + \caption{Spectra obtained for all-atom (AA) Au / butanethiol / solvent
897 +  systems. When butanethiol is deuterated (lower left), its
898 +  vibrational overlap with hexane decreases significantly.  Since
899 +  aromatic molecules and the butanethiol are vibrationally dissimilar,
900 +  the change is not as dramatic when toluene is the solvent (right).}
901 + \label{aahxntln}
902 + \end{figure}
903 +
904 + For the Au / butanethiol / toluene interfaces, having the AA
905 + butanethiol deuterated did not yield a significant change in the
906 + measured conductance. Compared to the C-H vibrational overlap between
907 + hexane and butanethiol, both of which have alkyl chains, the overlap
908 + between toluene and butanethiol is not as significant and thus does
909 + not contribute as much to the heat exchange process.
910 +
911 + Previous observations of Zhang {\it et al.}\cite{hase:2010} indicate
912 + that the {\it intra}molecular heat transport due to alkylthiols is
913 + highly efficient.  Combining our observations with those of Zhang {\it
914 +  et al.}, it appears that butanethiol acts as a channel to expedite
915 + heat flow from the gold surface and into the alkyl chain.  The
916 + vibrational coupling between the metal and the liquid phase can
917 + therefore be enhanced with the presence of suitable capping agents.
918 +
919 + Deuterated models in the UA force field did not decouple the thermal
920 + transport as well as in the AA force field.  The UA models, even
921 + though they have eliminated the high frequency C-H vibrational
922 + overlap, still have significant overlap in the lower-frequency
923 + portions of the infrared spectra (Figure \ref{uahxnua}).  Deuterating
924 + the UA models did not decouple the low frequency region enough to
925 + produce an observable difference for the results of $G$ (Table
926 + \ref{modelTest}).
927 +
928 + \begin{figure}
929   \includegraphics[width=\linewidth]{uahxnua}
930 < \caption{Vibrational spectra obtained for normal (upper) and
931 <  deuterated (lower) hexane in Au-butanethiol/hexane
932 <  systems. Butanethiol spectra are shown as reference. Both hexane and
933 <  butanethiol were using United-Atom models.}
930 > \caption{Vibrational power spectra for UA models for the butanethiol
931 >  and hexane solvent (upper panel) show the high degree of overlap
932 >  between these two molecules, particularly at lower frequencies.
933 >  Deuterating a UA model for the solvent (lower panel) does not
934 >  decouple the two spectra to the same degree as in the AA force
935 >  field (see Fig \ref{aahxntln}).}
936   \label{uahxnua}
937   \end{figure}
938  
939   \section{Conclusions}
940 < The NIVS algorithm we developed has been applied to simulations of
941 < Au-butanethiol surfaces with organic solvents. This algorithm allows
942 < effective unphysical thermal flux transferred between the metal and
943 < the liquid phase. With the flux applied, we were able to measure the
944 < corresponding thermal gradient and to obtain interfacial thermal
945 < conductivities. Under steady states, single trajectory simulation
946 < would be enough for accurate measurement. This would be advantageous
947 < compared to transient state simulations, which need multiple
849 < trajectories to produce reliable average results.
940 > The NIVS algorithm has been applied to simulations of
941 > butanethiol-capped Au(111) surfaces in the presence of organic
942 > solvents. This algorithm allows the application of unphysical thermal
943 > flux to transfer heat between the metal and the liquid phase. With the
944 > flux applied, we were able to measure the corresponding thermal
945 > gradients and to obtain interfacial thermal conductivities. Under
946 > steady states, 2-3 ns trajectory simulations are sufficient for
947 > computation of this quantity.
948  
949 < Our simulations have seen significant conductance enhancement with the
950 < presence of capping agent, compared to the bare gold / liquid
951 < interfaces. The acoustic impedance mismatch between the metal and the
952 < liquid phase is effectively eliminated by proper capping
953 < agent. Furthermore, the coverage precentage of the capping agent plays
954 < an important role in the interfacial thermal transport
955 < process. Moderately lower coverages allow higher contact between
956 < capping agent and solvent, and thus could further enhance the heat
957 < transfer process.
949 > Our simulations have seen significant conductance enhancement in the
950 > presence of capping agent, compared with the bare gold / liquid
951 > interfaces. The vibrational coupling between the metal and the liquid
952 > phase is enhanced by a chemically-bonded capping agent. Furthermore,
953 > the coverage percentage of the capping agent plays an important role
954 > in the interfacial thermal transport process. Moderately low coverages
955 > allow higher contact between capping agent and solvent, and thus could
956 > further enhance the heat transfer process, giving a non-monotonic
957 > behavior of conductance with increasing coverage.
958  
959 < Our measurement results, particularly of the UA models, agree with
960 < available experimental data. This indicates that our force field
863 < parameters have a nice description of the interactions between the
864 < particles at the interfaces. AA models tend to overestimate the
959 > Our results, particularly using the UA models, agree well with
960 > available experimental data.  The AA models tend to overestimate the
961   interfacial thermal conductance in that the classically treated C-H
962 < vibration would be overly sampled. Compared to the AA models, the UA
963 < models have higher computational efficiency with satisfactory
964 < accuracy, and thus are preferable in interfacial thermal transport
965 < modelings. Of the two definitions for $G$, the discrete form
962 > vibrations become too easily populated. Compared to the AA models, the
963 > UA models have higher computational efficiency with satisfactory
964 > accuracy, and thus are preferable in modeling interfacial thermal
965 > transport.
966 >
967 > Of the two definitions for $G$, the discrete form
968   (Eq. \ref{discreteG}) was easier to use and gives out relatively
969   consistent results, while the derivative form (Eq. \ref{derivativeG})
970   is not as versatile. Although $G^\prime$ gives out comparable results
971   and follows similar trend with $G$ when measuring close to fully
972 < covered or bare surfaces, the spatial resolution of $T$ profile is
973 < limited for accurate computation of derivatives data.
972 > covered or bare surfaces, the spatial resolution of $T$ profile
973 > required for the use of a derivative form is limited by the number of
974 > bins and the sampling required to obtain thermal gradient information.
975  
976 < Vlugt {\it et al.} has investigated the surface thiol structures for
977 < nanocrystal gold and pointed out that they differs from those of the
978 < Au(111) surface\cite{landman:1998,vlugt:cpc2007154}. This difference
979 < might lead to change of interfacial thermal transport behavior as
980 < well. To investigate this problem, an effective means to introduce
981 < thermal flux and measure the corresponding thermal gradient is
982 < desirable for simulating structures with spherical symmetry.
976 > Vlugt {\it et al.} have investigated the surface thiol structures for
977 > nanocrystalline gold and pointed out that they differ from those of
978 > the Au(111) surface.\cite{landman:1998,vlugt:cpc2007154} This
979 > difference could also cause differences in the interfacial thermal
980 > transport behavior. To investigate this problem, one would need an
981 > effective method for applying thermal gradients in non-planar
982 > (i.e. spherical) geometries.
983  
984   \section{Acknowledgments}
985   Support for this project was provided by the National Science
986   Foundation under grant CHE-0848243. Computational time was provided by
987   the Center for Research Computing (CRC) at the University of Notre
988   Dame.
989 +
990 + \section{Supporting Information}
991 + This information is available free of charge via the Internet at
992 + http://pubs.acs.org.
993 +
994   \newpage
995  
996   \bibliography{interfacial}

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines