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 3753 by skuang, Wed Jul 27 03:27:28 2011 UTC vs.
Revision 3763 by skuang, Tue Sep 27 21:02:48 2011 UTC

# Line 44 | Line 44 | Notre Dame, Indiana 46556}
44   \begin{doublespace}
45  
46   \begin{abstract}
47 <
48 < With the Non-Isotropic Velocity Scaling algorithm (NIVS) we have
49 < 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 <
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 the acoustic
54 >  impedance mismatch between the metal and liquid phases is
55 >  effectively reduced by the capping agents, leading to a greatly
56 >  enhanced conductivity at the interface.  Specifically, the chemical
57 >  bond between the metal and the capping agent introduces a
58 >  vibrational overlap that is not present without the capping agent,
59 >  and the overlap between the vibrational spectra (metal to cap, cap
60 >  to solvent) provides a mechanism for rapid thermal transport across
61 >  the interface. Our calculations also suggest that this is a
62 >  non-monotonic function of the fractional coverage of the surface, as
63 >  moderate coverages allow convective heat transport of solvent
64 >  molecules that have been in close contact with the capping agent.
65   \end{abstract}
66  
67   \newpage
# Line 73 | Line 73 | Due to the importance of heat flow in nanotechnology,
73   %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
74  
75   \section{Introduction}
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
76 > Due to the importance of heat flow (and heat removal) in
77 > nanotechnology, interfacial thermal conductance has been studied
78 > extensively both experimentally and computationally.\cite{cahill:793}
79 > Nanoscale materials have a significant fraction of their atoms at
80 > interfaces, and the chemical details of these interfaces govern the
81 > thermal transport properties.  Furthermore, the interfaces are often
82   heterogeneous (e.g. solid - liquid), which provides a challenge to
83 < traditional methods developed for homogeneous systems.
83 > computational methods which have been developed for homogeneous or
84 > bulk systems.
85  
86 < Experimentally, various interfaces have been investigated for their
87 < thermal conductance. Wang {\it et al.} studied heat transport through
88 < long-chain hydrocarbon monolayers on gold substrate at individual
89 < molecular level,\cite{Wang10082007} Schmidt {\it et al.} studied the
90 < role of CTAB on thermal transport between gold nanorods and
91 < solvent,\cite{doi:10.1021/jp8051888} and Juv\'e {\it et al.} studied
92 < the cooling dynamics, which is controlled by thermal interface
93 < resistence of glass-embedded metal
86 > Experimentally, the thermal properties of a number of interfaces have
87 > been investigated.  Cahill and coworkers studied nanoscale thermal
88 > transport from metal nanoparticle/fluid interfaces, to epitaxial
89 > TiN/single crystal oxides interfaces, and hydrophilic and hydrophobic
90 > interfaces between water and solids with different self-assembled
91 > monolayers.\cite{Wilson:2002uq,PhysRevB.67.054302,doi:10.1021/jp048375k,PhysRevLett.96.186101}
92 > Wang {\it et al.} studied heat transport through long-chain
93 > hydrocarbon monolayers on gold substrate at individual molecular
94 > level,\cite{Wang10082007} Schmidt {\it et al.} studied the role of
95 > cetyltrimethylammonium bromide (CTAB) on the thermal transport between
96 > gold nanorods and solvent,\cite{doi:10.1021/jp8051888} and Juv\'e {\it
97 >  et al.} studied the cooling dynamics, which is controlled by thermal
98 > interface resistance of glass-embedded metal
99   nanoparticles.\cite{PhysRevB.80.195406} Although interfaces are
100   normally considered barriers for heat transport, Alper {\it et al.}
101   suggested that specific ligands (capping agents) could completely
# Line 106 | Line 112 | difficult to measure with Equilibrium MD or forward NE
112   measurements for heat conductance of interfaces between the capping
113   monolayer on Au and a solvent phase have yet to be studied with their
114   approach. The comparatively low thermal flux through interfaces is
115 < difficult to measure with Equilibrium MD or forward NEMD simulation
115 > difficult to measure with Equilibrium
116 > MD\cite{doi:10.1080/0026897031000068578} or forward NEMD simulation
117   methods. Therefore, the Reverse NEMD (RNEMD)
118 < methods\cite{MullerPlathe:1997xw,kuang:164101} would have the
119 < advantage of applying this difficult to measure flux (while measuring
120 < the resulting gradient), given that the simulation methods being able
121 < to effectively apply an unphysical flux in non-homogeneous systems.
118 > methods\cite{MullerPlathe:1997xw,kuang:164101} would be advantageous
119 > in that they {\it apply} the difficult to measure quantity (flux),
120 > while {\it measuring} the easily-computed quantity (the thermal
121 > gradient).  This is particularly true for inhomogeneous interfaces
122 > where it would not be clear how to apply a gradient {\it a priori}.
123   Garde and coworkers\cite{garde:nl2005,garde:PhysRevLett2009} applied
124   this approach to various liquid interfaces and studied how thermal
125 < conductance (or resistance) is dependent on chemistry details of
126 < interfaces, e.g. hydrophobic and hydrophilic aqueous interfaces.
125 > conductance (or resistance) is dependent on chemical details of a
126 > number of hydrophobic and hydrophilic aqueous interfaces. {\bf And
127 >  Luo {\it et al.} studied the thermal conductance of Au-SAM-Au
128 >  junctions using the same approach, with comparison to a constant
129 >  temperature difference method\cite{Luo20101}. While this latter
130 >  approach establishes more thermal distributions compared to the
131 >  former RNEMD methods, it does not guarantee momentum or kinetic
132 >  energy conservations.}
133  
134   Recently, we have developed a Non-Isotropic Velocity Scaling (NIVS)
135   algorithm for RNEMD simulations\cite{kuang:164101}. This algorithm
# Line 129 | Line 143 | properties. Different models were used for both the ca
143   The work presented here deals with the Au(111) surface covered to
144   varying degrees by butanethiol, a capping agent with short carbon
145   chain, and solvated with organic solvents of different molecular
146 < properties. Different models were used for both the capping agent and
147 < the solvent force field parameters. Using the NIVS algorithm, the
148 < thermal transport across these interfaces was studied and the
149 < underlying mechanism for the phenomena was investigated.
146 > properties. {\bf To our knowledge, few previous MD inverstigations
147 >  have been found to address to these systems yet.} Different models
148 > were used for both the capping agent and the solvent force field
149 > parameters. Using the NIVS algorithm, the thermal transport across
150 > these interfaces was studied and the underlying mechanism for the
151 > phenomena was investigated.
152  
153   \section{Methodology}
154 < \subsection{Imposd-Flux Methods in MD Simulations}
154 > \subsection{Imposed-Flux Methods in MD Simulations}
155   Steady state MD simulations have an advantage in that not many
156   trajectories are needed to study the relationship between thermal flux
157   and thermal gradients. For systems with low interfacial conductance,
# Line 159 | Line 175 | can be applied between regions of particles of arbitar
175   kinetic energy fluxes without obvious perturbation to the velocity
176   distributions of the simulated systems. Furthermore, this approach has
177   the advantage in heterogeneous interfaces in that kinetic energy flux
178 < can be applied between regions of particles of arbitary identity, and
178 > can be applied between regions of particles of arbitrary identity, and
179   the flux will not be restricted by difference in particle mass.
180  
181   The NIVS algorithm scales the velocity vectors in two separate regions
182 < of a simulation system with respective diagonal scaling matricies. To
183 < determine these scaling factors in the matricies, a set of equations
182 > of a simulation system with respective diagonal scaling matrices. To
183 > determine these scaling factors in the matrices, a set of equations
184   including linear momentum conservation and kinetic energy conservation
185   constraints and target energy flux satisfaction is solved. With the
186   scaling operation applied to the system in a set frequency, bulk
# Line 187 | Line 203 | temperature of the two separated phases.
203   where ${E_{total}}$ is the total imposed non-physical kinetic energy
204   transfer during the simulation and ${\langle T_\mathrm{hot}\rangle}$
205   and ${\langle T_\mathrm{cold}\rangle}$ are the average observed
206 < temperature of the two separated phases.
206 > temperature of the two separated phases.  For an applied flux $J_z$
207 > operating over a simulation time $t$ on a periodically-replicated slab
208 > of dimensions $L_x \times L_y$, $E_{total} = J_z *(t)*(2 L_x L_y)$.
209  
210   When the interfacial conductance is {\it not} small, there are two
211   ways to define $G$. One common way is to assume the temperature is
212   discrete on the two sides of the interface. $G$ can be calculated
213   using the applied thermal flux $J$ and the maximum temperature
214   difference measured along the thermal gradient max($\Delta T$), which
215 < occurs at the Gibbs deviding surface (Figure \ref{demoPic}):
215 > occurs at the Gibbs dividing surface (Figure \ref{demoPic}). This is
216 > known as the Kapitza conductance, which is the inverse of the Kapitza
217 > resistance.
218   \begin{equation}
219    G=\frac{J}{\Delta T}
220   \label{discreteG}
# Line 205 | Line 225 | occurs at the Gibbs deviding surface (Figure \ref{demo
225   \caption{Interfacial conductance can be calculated by applying an
226    (unphysical) kinetic energy flux between two slabs, one located
227    within the metal and another on the edge of the periodic box.  The
228 <  system responds by forming a thermal response or a gradient.  In
229 <  bulk liquids, this gradient typically has a single slope, but in
230 <  interfacial systems, there are distinct thermal conductivity
231 <  domains.  The interfacial conductance, $G$ is found by measuring the
232 <  temperature gap at the Gibbs dividing surface, or by using second
233 <  derivatives of the thermal profile.}
228 >  system responds by forming a thermal gradient.  In bulk liquids,
229 >  this gradient typically has a single slope, but in interfacial
230 >  systems, there are distinct thermal conductivity domains.  The
231 >  interfacial conductance, $G$ is found by measuring the temperature
232 >  gap at the Gibbs dividing surface, or by using second derivatives of
233 >  the thermal profile.}
234   \label{demoPic}
235   \end{figure}
236  
# Line 249 | Line 269 | profile.
269  
270   \begin{figure}
271   \includegraphics[width=\linewidth]{gradT}
272 < \caption{A sample of Au-butanethiol/hexane interfacial system and the
273 <  temperature profile after a kinetic energy flux is imposed to
274 <  it. The 1st and 2nd derivatives of the temperature profile can be
275 <  obtained with finite difference approximation (lower panel).}
272 > \caption{A sample of Au (111) / butanethiol / hexane interfacial
273 >  system with the temperature profile after a kinetic energy flux has
274 >  been imposed.  Note that the largest temperature jump in the thermal
275 >  profile (corresponding to the lowest interfacial conductance) is at
276 >  the interface between the butanethiol molecules (blue) and the
277 >  solvent (grey).  First and second derivatives of the temperature
278 >  profile are obtained using a finite difference approximation (lower
279 >  panel).}
280   \label{gradT}
281   \end{figure}
282  
# Line 299 | Line 323 | between periodic images of the gold interfaces is $45
323   solvent molecules would change the normal behavior of the liquid
324   phase. Therefore, our $N_{solvent}$ values were chosen to ensure that
325   these extreme cases did not happen to our simulations. The spacing
326 < between periodic images of the gold interfaces is $45 \sim 75$\AA.
326 > between periodic images of the gold interfaces is $45 \sim 75$\AA in
327 > our simulations.
328  
329   The initial configurations generated are further equilibrated with the
330   $x$ and $y$ dimensions fixed, only allowing the $z$-length scale to
# Line 317 | Line 342 | gradient had stablized, the temperature profile of the
342   $\sim$200K. Therefore, thermal flux usually came from the metal to the
343   liquid so that the liquid has a higher temperature and would not
344   freeze due to lowered temperatures. After this induced temperature
345 < gradient had stablized, the temperature profile of the simulation cell
346 < was recorded. To do this, the simulation cell is devided evenly into
345 > gradient had stabilized, the temperature profile of the simulation cell
346 > was recorded. To do this, the simulation cell is divided evenly into
347   $N$ slabs along the $z$-axis. The average temperatures of each slab
348   are recorded for 1$\sim$2 ns. When the slab width $d$ of each slab is
349   the same, the derivatives of $T$ with respect to slab number $n$ can
# Line 349 | Line 374 | particles of different species.
374   \caption{Structures of the capping agent and solvents utilized in
375    these simulations. The chemically-distinct sites (a-e) are expanded
376    in terms of constituent atoms for both United Atom (UA) and All Atom
377 <  (AA) force fields.  Most parameters are from
378 <  Refs. \protect\cite{TraPPE-UA.alkanes,TraPPE-UA.alkylbenzenes,TraPPE-UA.thiols} (UA) and \protect\cite{OPLSAA} (AA). Cross-interactions with the Au atoms are given in Table \ref{MnM}.}
377 >  (AA) force fields.  Most parameters are from References
378 >  \protect\cite{TraPPE-UA.alkanes,TraPPE-UA.alkylbenzenes,TraPPE-UA.thiols}
379 >  (UA) and \protect\cite{OPLSAA} (AA). Cross-interactions with the Au
380 >  atoms are given in Table \ref{MnM}.}
381   \label{demoMol}
382   \end{figure}
383  
# Line 372 | Line 399 | However, the TraPPE-UA model for alkanes is known to p
399  
400   By eliminating explicit hydrogen atoms, the TraPPE-UA models are
401   simple and computationally efficient, while maintaining good accuracy.
402 < However, the TraPPE-UA model for alkanes is known to predict a slighly
402 > However, the TraPPE-UA model for alkanes is known to predict a slightly
403   lower boiling point than experimental values. This is one of the
404   reasons we used a lower average temperature (200K) for our
405   simulations. If heat is transferred to the liquid phase during the
# Line 470 | Line 497 | our simulations.
497    \end{minipage}
498   \end{table*}
499  
473 \subsection{Vibrational Power Spectrum}
500  
501 < To investigate the mechanism of interfacial thermal conductance, the
502 < vibrational power spectrum was computed. Power spectra were taken for
503 < individual components in different simulations. To obtain these
504 < spectra, simulations were run after equilibration, in the NVE
505 < ensemble, and without a thermal gradient. Snapshots of configurations
506 < were collected at a frequency that is higher than that of the fastest
507 < vibrations occuring in the simulations. With these configurations, the
508 < velocity auto-correlation functions can be computed:
483 < \begin{equation}
484 < C_A (t) = \langle\vec{v}_A (t)\cdot\vec{v}_A (0)\rangle
485 < \label{vCorr}
486 < \end{equation}
487 < The power spectrum is constructed via a Fourier transform of the
488 < symmetrized velocity autocorrelation function,
489 < \begin{equation}
490 <  \hat{f}(\omega) =
491 <  \int_{-\infty}^{\infty} C_A (t) e^{-2\pi it\omega}\,dt
492 < \label{fourier}
493 < \end{equation}
494 <
495 < \section{Results and Discussions}
496 < In what follows, how the parameters and protocol of simulations would
497 < affect the measurement of $G$'s is first discussed. With a reliable
498 < protocol and set of parameters, the influence of capping agent
499 < coverage on thermal conductance is investigated. Besides, different
500 < force field models for both solvents and selected deuterated models
501 < were tested and compared. Finally, a summary of the role of capping
502 < agent in the interfacial thermal transport process is given.
503 <
504 < \subsection{How Simulation Parameters Affects $G$}
505 < We have varied our protocol or other parameters of the simulations in
506 < order to investigate how these factors would affect the measurement of
507 < $G$'s. It turned out that while some of these parameters would not
508 < affect the results substantially, some other changes to the
509 < simulations would have a significant impact on the measurement
510 < results.
501 > \section{Results}
502 > There are many factors contributing to the measured interfacial
503 > conductance; some of these factors are physically motivated
504 > (e.g. coverage of the surface by the capping agent coverage and
505 > solvent identity), while some are governed by parameters of the
506 > methodology (e.g. applied flux and the formulas used to obtain the
507 > conductance). In this section we discuss the major physical and
508 > calculational effects on the computed conductivity.
509  
510 < In some of our simulations, we allowed $L_x$ and $L_y$ to change
513 < during equilibrating the liquid phase. Due to the stiffness of the
514 < crystalline Au structure, $L_x$ and $L_y$ would not change noticeably
515 < after equilibration. Although $L_z$ could fluctuate $\sim$1\% after a
516 < system is fully equilibrated in the NPT ensemble, this fluctuation, as
517 < well as those of $L_x$ and $L_y$ (which is significantly smaller),
518 < would not be magnified on the calculated $G$'s, as shown in Table
519 < \ref{AuThiolHexaneUA}. This insensivity to $L_i$ fluctuations allows
520 < reliable measurement of $G$'s without the necessity of extremely
521 < cautious equilibration process.
510 > \subsection{Effects due to capping agent coverage}
511  
512 < As stated in our computational details, the spacing filled with
513 < solvent molecules can be chosen within a range. This allows some
514 < change of solvent molecule numbers for the same Au-butanethiol
515 < surfaces. We did this study on our Au-butanethiol/hexane
516 < simulations. Nevertheless, the results obtained from systems of
517 < different $N_{hexane}$ did not indicate that the measurement of $G$ is
529 < susceptible to this parameter. For computational efficiency concern,
530 < smaller system size would be preferable, given that the liquid phase
531 < structure is not affected.
532 <
533 < Our NIVS algorithm allows change of unphysical thermal flux both in
534 < direction and in quantity. This feature extends our investigation of
535 < interfacial thermal conductance. However, the magnitude of this
536 < thermal flux is not arbitary if one aims to obtain a stable and
537 < reliable thermal gradient. A temperature profile would be
538 < substantially affected by noise when $|J_z|$ has a much too low
539 < magnitude; while an excessively large $|J_z|$ that overwhelms the
540 < conductance capacity of the interface would prevent a thermal gradient
541 < to reach a stablized steady state. NIVS has the advantage of allowing
542 < $J$ to vary in a wide range such that the optimal flux range for $G$
543 < measurement can generally be simulated by the algorithm. Within the
544 < optimal range, we were able to study how $G$ would change according to
545 < the thermal flux across the interface. For our simulations, we denote
546 < $J_z$ to be positive when the physical thermal flux is from the liquid
547 < to metal, and negative vice versa. The $G$'s measured under different
548 < $J_z$ is listed in Table \ref{AuThiolHexaneUA} and
549 < \ref{AuThiolToluene}. These results do not suggest that $G$ is
550 < dependent on $J_z$ within this flux range. The linear response of flux
551 < to thermal gradient simplifies our investigations in that we can rely
552 < on $G$ measurement with only a couple $J_z$'s and do not need to test
553 < a large series of fluxes.
554 <
555 < \begin{table*}
556 <  \begin{minipage}{\linewidth}
557 <    \begin{center}
558 <      \caption{Computed interfacial thermal conductivity ($G$ and
559 <        $G^\prime$) values for the 100\% covered Au-butanethiol/hexane
560 <        interfaces with UA model and different hexane molecule numbers
561 <        at different temperatures using a range of energy
562 <        fluxes. Error estimates indicated in parenthesis.}
563 <      
564 <      \begin{tabular}{ccccccc}
565 <        \hline\hline
566 <        $\langle T\rangle$ & $N_{hexane}$ & Fixed & $\rho_{hexane}$ &
567 <        $J_z$ & $G$ & $G^\prime$ \\
568 <        (K) & & $L_x$ \& $L_y$? & (g/cm$^3$) & (GW/m$^2$) &
569 <        \multicolumn{2}{c}{(MW/m$^2$/K)} \\
570 <        \hline
571 <        200 & 266 & No  & 0.672 & -0.96 & 102(3)    & 80.0(0.8) \\
572 <            & 200 & Yes & 0.694 &  1.92 & 129(11)   & 87.3(0.3) \\
573 <            &     & Yes & 0.672 &  1.93 & 131(16)   & 78(13)    \\
574 <            &     & No  & 0.688 &  0.96 & 125(16)   & 90.2(15)  \\
575 <            &     &     &       &  1.91 & 139(10)   & 101(10)   \\
576 <            &     &     &       &  2.83 & 141(6)    & 89.9(9.8) \\
577 <            & 166 & Yes & 0.679 &  0.97 & 115(19)   & 69(18)    \\
578 <            &     &     &       &  1.94 & 125(9)    & 87.1(0.2) \\
579 <            &     & No  & 0.681 &  0.97 & 141(30)   & 78(22)    \\
580 <            &     &     &       &  1.92 & 138(4)    & 98.9(9.5) \\
581 <        \hline
582 <        250 & 200 & No  & 0.560 &  0.96 & 75(10)    & 61.8(7.3) \\
583 <            &     &     &       & -0.95 & 49.4(0.3) & 45.7(2.1) \\
584 <            & 166 & Yes & 0.570 &  0.98 & 79.0(3.5) & 62.9(3.0) \\
585 <            &     & No  & 0.569 &  0.97 & 80.3(0.6) & 67(11)    \\
586 <            &     &     &       &  1.44 & 76.2(5.0) & 64.8(3.8) \\
587 <            &     &     &       & -0.95 & 56.4(2.5) & 54.4(1.1) \\
588 <            &     &     &       & -1.85 & 47.8(1.1) & 53.5(1.5) \\
589 <        \hline\hline
590 <      \end{tabular}
591 <      \label{AuThiolHexaneUA}
592 <    \end{center}
593 <  \end{minipage}
594 < \end{table*}
512 > A series of different initial conditions with a range of surface
513 > coverages was prepared and solvated with various with both of the
514 > solvent molecules. These systems were then equilibrated and their
515 > interfacial thermal conductivity was measured with the NIVS
516 > algorithm. Figure \ref{coverage} demonstrates the trend of conductance
517 > with respect to surface coverage.
518  
596 Furthermore, we also attempted to increase system average temperatures
597 to above 200K. These simulations are first equilibrated in the NPT
598 ensemble under normal pressure. As stated above, the TraPPE-UA model
599 for hexane tends to predict a lower boiling point. In our simulations,
600 hexane had diffculty to remain in liquid phase when NPT equilibration
601 temperature is higher than 250K. Additionally, the equilibrated liquid
602 hexane density under 250K becomes lower than experimental value. This
603 expanded liquid phase leads to lower contact between hexane and
604 butanethiol as well.[MAY NEED SLAB DENSITY FIGURE]
605 And this reduced contact would
606 probably be accountable for a lower interfacial thermal conductance,
607 as shown in Table \ref{AuThiolHexaneUA}.
608
609 A similar study for TraPPE-UA toluene agrees with the above result as
610 well. Having a higher boiling point, toluene tends to remain liquid in
611 our simulations even equilibrated under 300K in NPT
612 ensembles. Furthermore, the expansion of the toluene liquid phase is
613 not as significant as that of the hexane. This prevents severe
614 decrease of liquid-capping agent contact and the results (Table
615 \ref{AuThiolToluene}) show only a slightly decreased interface
616 conductance. Therefore, solvent-capping agent contact should play an
617 important role in the thermal transport process across the interface
618 in that higher degree of contact could yield increased conductance.
619
620 \begin{table*}
621  \begin{minipage}{\linewidth}
622    \begin{center}
623      \caption{Computed interfacial thermal conductivity ($G$ and
624        $G^\prime$) values for a 90\% coverage Au-butanethiol/toluene
625        interface at different temperatures using a range of energy
626        fluxes. Error estimates indicated in parenthesis.}
627      
628      \begin{tabular}{ccccc}
629        \hline\hline
630        $\langle T\rangle$ & $\rho_{toluene}$ & $J_z$ & $G$ & $G^\prime$ \\
631        (K) & (g/cm$^3$) & (GW/m$^2$) & \multicolumn{2}{c}{(MW/m$^2$/K)} \\
632        \hline
633        200 & 0.933 &  2.15 & 204(12) & 113(12) \\
634            &       & -1.86 & 180(3)  & 135(21) \\
635            &       & -3.93 & 176(5)  & 113(12) \\
636        \hline
637        300 & 0.855 & -1.91 & 143(5)  & 125(2)  \\
638            &       & -4.19 & 135(9)  & 113(12) \\
639        \hline\hline
640      \end{tabular}
641      \label{AuThiolToluene}
642    \end{center}
643  \end{minipage}
644 \end{table*}
645
646 Besides lower interfacial thermal conductance, surfaces in relatively
647 high temperatures are susceptible to reconstructions, when
648 butanethiols have a full coverage on the Au(111) surface. These
649 reconstructions include surface Au atoms migrated outward to the S
650 atom layer, and butanethiol molecules embedded into the original
651 surface Au layer. The driving force for this behavior is the strong
652 Au-S interactions in our simulations. And these reconstructions lead
653 to higher ratio of Au-S attraction and thus is energetically
654 favorable. Furthermore, this phenomenon agrees with experimental
655 results\cite{doi:10.1021/j100035a033,doi:10.1021/la026493y}. Vlugt
656 {\it et al.} had kept their Au(111) slab rigid so that their
657 simulations can reach 300K without surface reconstructions. Without
658 this practice, simulating 100\% thiol covered interfaces under higher
659 temperatures could hardly avoid surface reconstructions. However, our
660 measurement is based on assuming homogeneity on $x$ and $y$ dimensions
661 so that measurement of $T$ at particular $z$ would be an effective
662 average of the particles of the same type. Since surface
663 reconstructions could eliminate the original $x$ and $y$ dimensional
664 homogeneity, measurement of $G$ is more difficult to conduct under
665 higher temperatures. Therefore, most of our measurements are
666 undertaken at $\langle T\rangle\sim$200K.
667
668 However, when the surface is not completely covered by butanethiols,
669 the simulated system is more resistent to the reconstruction
670 above. Our Au-butanethiol/toluene system had the Au(111) surfaces 90\%
671 covered by butanethiols, but did not see this above phenomena even at
672 $\langle T\rangle\sim$300K. The empty three-fold sites not occupied by
673 capping agents could help prevent surface reconstruction in that they
674 provide other means of capping agent relaxation. It is observed that
675 butanethiols can migrate to their neighbor empty sites during a
676 simulation. Therefore, we were able to obtain $G$'s for these
677 interfaces even at a relatively high temperature without being
678 affected by surface reconstructions.
679
680 \subsection{Influence of Capping Agent Coverage on $G$}
681 To investigate the influence of butanethiol coverage on interfacial
682 thermal conductance, a series of different coverage Au-butanethiol
683 surfaces is prepared and solvated with various organic
684 molecules. These systems are then equilibrated and their interfacial
685 thermal conductivity are measured with our NIVS algorithm. Figure
686 \ref{coverage} demonstrates the trend of conductance change with
687 respect to different coverages of butanethiol. To study the isotope
688 effect in interfacial thermal conductance, deuterated UA-hexane is
689 included as well.
690
519   \begin{figure}
520   \includegraphics[width=\linewidth]{coverage}
521 < \caption{Comparison of interfacial thermal conductivity ($G$) values
522 <  for the Au-butanethiol/solvent interface with various UA models and
523 <  different capping agent coverages at $\langle T\rangle\sim$200K
524 <  using certain energy flux respectively.}
521 > \caption{The interfacial thermal conductivity ($G$) has a
522 >  non-monotonic dependence on the degree of surface capping.  This
523 >  data is for the Au(111) / butanethiol / solvent interface with
524 >  various UA force fields at $\langle T\rangle \sim $200K.}
525   \label{coverage}
526   \end{figure}
527  
528 < It turned out that with partial covered butanethiol on the Au(111)
529 < surface, the derivative definition for $G^\prime$
530 < (Eq. \ref{derivativeG}) was difficult to apply, due to the difficulty
531 < in locating the maximum of change of $\lambda$. Instead, the discrete
532 < definition (Eq. \ref{discreteG}) is easier to apply, as the Gibbs
533 < deviding surface can still be well-defined. Therefore, $G$ (not
706 < $G^\prime$) was used for this section.
528 > In partially covered surfaces, the derivative definition for
529 > $G^\prime$ (Eq. \ref{derivativeG}) becomes difficult to apply, as the
530 > location of maximum change of $\lambda$ becomes washed out.  The
531 > discrete definition (Eq. \ref{discreteG}) is easier to apply, as the
532 > Gibbs dividing surface is still well-defined. Therefore, $G$ (not
533 > $G^\prime$) was used in this section.
534  
535   From Figure \ref{coverage}, one can see the significance of the
536 < presence of capping agents. Even when a fraction of the Au(111)
537 < surface sites are covered with butanethiols, the conductivity would
538 < see an enhancement by at least a factor of 3. This indicates the
539 < important role cappping agent is playing for thermal transport
540 < phenomena on metal / organic solvent surfaces.
536 > presence of capping agents. When even a small fraction of the Au(111)
537 > surface sites are covered with butanethiols, the conductivity exhibits
538 > an enhancement by at least a factor of 3.  Capping agents are clearly
539 > playing a major role in thermal transport at metal / organic solvent
540 > surfaces.
541  
542 < Interestingly, as one could observe from our results, the maximum
543 < conductance enhancement (largest $G$) happens while the surfaces are
544 < about 75\% covered with butanethiols. This again indicates that
545 < solvent-capping agent contact has an important role of the thermal
546 < transport process. Slightly lower butanethiol coverage allows small
547 < gaps between butanethiols to form. And these gaps could be filled with
548 < solvent molecules, which acts like ``heat conductors'' on the
549 < surface. The higher degree of interaction between these solvent
550 < molecules and capping agents increases the enhancement effect and thus
724 < produces a higher $G$ than densely packed butanethiol arrays. However,
725 < once this maximum conductance enhancement is reached, $G$ decreases
726 < when butanethiol coverage continues to decrease. Each capping agent
727 < molecule reaches its maximum capacity for thermal
728 < conductance. Therefore, even higher solvent-capping agent contact
729 < would not offset this effect. Eventually, when butanethiol coverage
730 < continues to decrease, solvent-capping agent contact actually
731 < decreases with the disappearing of butanethiol molecules. In this
732 < case, $G$ decrease could not be offset but instead accelerated. [MAY NEED
733 < SNAPSHOT SHOWING THE PHENOMENA / SLAB DENSITY ANALYSIS]
542 > We note a non-monotonic behavior in the interfacial conductance as a
543 > function of surface coverage. The maximum conductance (largest $G$)
544 > happens when the surfaces are about 75\% covered with butanethiol
545 > caps.  The reason for this behavior is not entirely clear.  One
546 > explanation is that incomplete butanethiol coverage allows small gaps
547 > between butanethiols to form. These gaps can be filled by transient
548 > solvent molecules.  These solvent molecules couple very strongly with
549 > the hot capping agent molecules near the surface, and can then carry
550 > away (diffusively) the excess thermal energy from the surface.
551  
552 < A comparison of the results obtained from differenet organic solvents
553 < can also provide useful information of the interfacial thermal
554 < transport process. The deuterated hexane (UA) results do not appear to
555 < be much different from those of normal hexane (UA), given that
556 < butanethiol (UA) is non-deuterated for both solvents. These UA model
557 < studies, even though eliminating C-H vibration samplings, still have
741 < C-C vibrational frequencies different from each other. However, these
742 < differences in the infrared range do not seem to produce an observable
743 < difference for the results of $G$ (Figure \ref{uahxnua}).
552 > There appears to be a competition between the conduction of the
553 > thermal energy away from the surface by the capping agents (enhanced
554 > by greater coverage) and the coupling of the capping agents with the
555 > solvent (enhanced by interdigitation at lower coverages).  This
556 > competition would lead to the non-monotonic coverage behavior observed
557 > here.
558  
559 < \begin{figure}
560 < \includegraphics[width=\linewidth]{uahxnua}
561 < \caption{Vibrational spectra obtained for normal (upper) and
562 <  deuterated (lower) hexane in Au-butanethiol/hexane
563 <  systems. Butanethiol spectra are shown as reference. Both hexane and
564 <  butanethiol were using United-Atom models.}
565 < \label{uahxnua}
566 < \end{figure}
559 > Results for rigid body toluene solvent, as well as the UA hexane, are
560 > within the ranges expected from prior experimental
561 > work.\cite{Wilson:2002uq,cahill:793,PhysRevB.80.195406} This suggests
562 > that explicit hydrogen atoms might not be required for modeling
563 > thermal transport in these systems.  C-H vibrational modes do not see
564 > significant excited state population at low temperatures, and are not
565 > likely to carry lower frequency excitations from the solid layer into
566 > the bulk liquid.
567  
568 < Furthermore, results for rigid body toluene solvent, as well as other
569 < UA-hexane solvents, are reasonable within the general experimental
570 < ranges\cite{Wilson:2002uq,cahill:793,PhysRevB.80.195406}. This
571 < suggests that explicit hydrogen might not be a required factor for
572 < modeling thermal transport phenomena of systems such as
573 < Au-thiol/organic solvent.
568 > The toluene solvent does not exhibit the same behavior as hexane in
569 > that $G$ remains at approximately the same magnitude when the capping
570 > coverage increases from 25\% to 75\%.  Toluene, as a rigid planar
571 > molecule, cannot occupy the relatively small gaps between the capping
572 > agents as easily as the chain-like {\it n}-hexane.  The effect of
573 > solvent coupling to the capping agent is therefore weaker in toluene
574 > except at the very lowest coverage levels.  This effect counters the
575 > coverage-dependent conduction of heat away from the metal surface,
576 > leading to a much flatter $G$ vs. coverage trend than is observed in
577 > {\it n}-hexane.
578  
579 < However, results for Au-butanethiol/toluene do not show an identical
580 < trend with those for Au-butanethiol/hexane in that $G$ remains at
581 < approximately the same magnitue when butanethiol coverage differs from
582 < 25\% to 75\%. This might be rooted in the molecule shape difference
583 < for planar toluene and chain-like {\it n}-hexane. Due to this
584 < difference, toluene molecules have more difficulty in occupying
585 < relatively small gaps among capping agents when their coverage is not
586 < too low. Therefore, the solvent-capping agent contact may keep
587 < increasing until the capping agent coverage reaches a relatively low
770 < level. This becomes an offset for decreasing butanethiol molecules on
771 < its effect to the process of interfacial thermal transport. Thus, one
772 < can see a plateau of $G$ vs. butanethiol coverage in our results.
773 <
774 < \subsection{Influence of Chosen Molecule Model on $G$}
775 < In addition to UA solvent/capping agent models, AA models are included
776 < in our simulations as well. Besides simulations of the same (UA or AA)
777 < model for solvent and capping agent, different models can be applied
778 < to different components. Furthermore, regardless of models chosen,
779 < either the solvent or the capping agent can be deuterated, similar to
780 < the previous section. Table \ref{modelTest} summarizes the results of
781 < these studies.
579 > \subsection{Effects due to Solvent \& Solvent Models}
580 > In addition to UA solvent and capping agent models, AA models have
581 > also been included in our simulations.  In most of this work, the same
582 > (UA or AA) model for solvent and capping agent was used, but it is
583 > also possible to utilize different models for different components.
584 > We have also included isotopic substitutions (Hydrogen to Deuterium)
585 > to decrease the explicit vibrational overlap between solvent and
586 > capping agent. Table \ref{modelTest} summarizes the results of these
587 > studies.
588  
589   \begin{table*}
590    \begin{minipage}{\linewidth}
591      \begin{center}
592        
593 <      \caption{Computed interfacial thermal conductivity ($G$ and
593 >      \caption{Computed interfacial thermal conductance ($G$ and
594          $G^\prime$) values for interfaces using various models for
595          solvent and capping agent (or without capping agent) at
596 <        $\langle T\rangle\sim$200K. (D stands for deuterated solvent
597 <        or capping agent molecules; ``Avg.'' denotes results that are
598 <        averages of simulations under different $J_z$'s. Error
599 <        estimates indicated in parenthesis.)}
596 >        $\langle T\rangle\sim$200K.  Here ``D'' stands for deuterated
597 >        solvent or capping agent molecules; ``Avg.'' denotes results
598 >        that are averages of simulations under different applied
599 >        thermal flux $(J_z)$ values. Error estimates are indicated in
600 >        parentheses.}
601        
602        \begin{tabular}{llccc}
603          \hline\hline
# Line 825 | Line 632 | To facilitate direct comparison, the same system with
632    \end{minipage}
633   \end{table*}
634  
635 < To facilitate direct comparison, the same system with differnt models
636 < for different components uses the same length scale for their
637 < simulation cells. Without the presence of capping agent, using
831 < different models for hexane yields similar results for both $G$ and
832 < $G^\prime$, and these two definitions agree with eath other very
833 < well. This indicates very weak interaction between the metal and the
834 < solvent, and is a typical case for acoustic impedance mismatch between
835 < these two phases.
635 > To facilitate direct comparison between force fields, systems with the
636 > same capping agent and solvent were prepared with the same length
637 > scales for the simulation cells.
638  
639 < As for Au(111) surfaces completely covered by butanethiols, the choice
640 < of models for capping agent and solvent could impact the measurement
641 < of $G$ and $G^\prime$ quite significantly. For Au-butanethiol/hexane
642 < interfaces, using AA model for both butanethiol and hexane yields
643 < substantially higher conductivity values than using UA model for at
644 < least one component of the solvent and capping agent, which exceeds
843 < the general range of experimental measurement results. This is
844 < probably due to the classically treated C-H vibrations in the AA
845 < model, which should not be appreciably populated at normal
846 < temperatures. In comparison, once either the hexanes or the
847 < butanethiols are deuterated, one can see a significantly lower $G$ and
848 < $G^\prime$. In either of these cases, the C-H(D) vibrational overlap
849 < between the solvent and the capping agent is removed (Figure
850 < \ref{aahxntln}). Conclusively, the improperly treated C-H vibration in
851 < the AA model produced over-predicted results accordingly. Compared to
852 < the AA model, the UA model yields more reasonable results with higher
853 < computational efficiency.
639 > On bare metal / solvent surfaces, different force field models for
640 > hexane yield similar results for both $G$ and $G^\prime$, and these
641 > two definitions agree with each other very well. This is primarily an
642 > indicator of weak interactions between the metal and the solvent, and
643 > is a typical case for acoustic impedance mismatch between these two
644 > phases.  
645  
646 < \begin{figure}
647 < \includegraphics[width=\linewidth]{aahxntln}
648 < \caption{Spectra obtained for All-Atom model Au-butanethil/solvent
649 <  systems. When butanethiol is deuterated (lower left), its
650 <  vibrational overlap with hexane would decrease significantly,
651 <  compared with normal butanethiol (upper left). However, this
652 <  dramatic change does not apply to toluene as much (right).}
653 < \label{aahxntln}
654 < \end{figure}
646 > For the fully-covered surfaces, the choice of force field for the
647 > capping agent and solvent has a large impact on the calculated values
648 > of $G$ and $G^\prime$.  The AA thiol to AA solvent conductivities are
649 > much larger than their UA to UA counterparts, and these values exceed
650 > the experimental estimates by a large measure.  The AA force field
651 > allows significant energy to go into C-H (or C-D) stretching modes,
652 > and since these modes are high frequency, this non-quantum behavior is
653 > likely responsible for the overestimate of the conductivity.  Compared
654 > to the AA model, the UA model yields more reasonable conductivity
655 > values with much higher computational efficiency.
656  
657 < However, for Au-butanethiol/toluene interfaces, having the AA
658 < butanethiol deuterated did not yield a significant change in the
659 < measurement results. Compared to the C-H vibrational overlap between
660 < hexane and butanethiol, both of which have alkyl chains, that overlap
661 < between toluene and butanethiol is not so significant and thus does
662 < not have as much contribution to the heat exchange
663 < process. Conversely, extra degrees of freedom such as the C-H
664 < vibrations could yield higher heat exchange rate between these two
665 < phases and result in a much higher conductivity.
657 > \subsubsection{Are electronic excitations in the metal important?}
658 > Because they lack electronic excitations, the QSC and related embedded
659 > atom method (EAM) models for gold are known to predict unreasonably
660 > low values for bulk conductivity
661 > ($\lambda$).\cite{kuang:164101,ISI:000207079300006,Clancy:1992} If the
662 > conductance between the phases ($G$) is governed primarily by phonon
663 > excitation (and not electronic degrees of freedom), one would expect a
664 > classical model to capture most of the interfacial thermal
665 > conductance.  Our results for $G$ and $G^\prime$ indicate that this is
666 > indeed the case, and suggest that the modeling of interfacial thermal
667 > transport depends primarily on the description of the interactions
668 > between the various components at the interface.  When the metal is
669 > chemically capped, the primary barrier to thermal conductivity appears
670 > to be the interface between the capping agent and the surrounding
671 > solvent, so the excitations in the metal have little impact on the
672 > value of $G$.
673  
674 < Although the QSC model for Au is known to predict an overly low value
876 < for bulk metal gold conductivity\cite{kuang:164101}, our computational
877 < results for $G$ and $G^\prime$ do not seem to be affected by this
878 < drawback of the model for metal. Instead, our results suggest that the
879 < modeling of interfacial thermal transport behavior relies mainly on
880 < the accuracy of the interaction descriptions between components
881 < occupying the interfaces.
674 > \subsection{Effects due to methodology and simulation parameters}
675  
676 < \subsection{Role of Capping Agent in Interfacial Thermal Conductance}
677 < The vibrational spectra for gold slabs in different environments are
678 < shown as in Figure \ref{specAu}. Regardless of the presence of
679 < solvent, the gold surfaces covered by butanethiol molecules, compared
680 < to bare gold surfaces, exhibit an additional peak observed at the
681 < frequency of $\sim$170cm$^{-1}$, which is attributed to the S-Au
682 < bonding vibration. This vibration enables efficient thermal transport
683 < from surface Au layer to the capping agents. Therefore, in our
684 < simulations, the Au/S interfaces do not appear major heat barriers
892 < compared to the butanethiol / solvent interfaces.
676 > We have varied the parameters of the simulations in order to
677 > investigate how these factors would affect the computation of $G$.  Of
678 > particular interest are: 1) the length scale for the applied thermal
679 > gradient (modified by increasing the amount of solvent in the system),
680 > 2) the sign and magnitude of the applied thermal flux, 3) the average
681 > temperature of the simulation (which alters the solvent density during
682 > equilibration), and 4) the definition of the interfacial conductance
683 > (Eqs. (\ref{discreteG}) or (\ref{derivativeG})) used in the
684 > calculation.
685  
686 < Simultaneously, the vibrational overlap between butanethiol and
687 < organic solvents suggests higher thermal exchange efficiency between
688 < these two components. Even exessively high heat transport was observed
689 < when All-Atom models were used and C-H vibrations were treated
690 < classically. Compared to metal and organic liquid phase, the heat
691 < transfer efficiency between butanethiol and organic solvents is closer
692 < to that within bulk liquid phase.
686 > Systems of different lengths were prepared by altering the number of
687 > solvent molecules and extending the length of the box along the $z$
688 > axis to accomodate the extra solvent.  Equilibration at the same
689 > temperature and pressure conditions led to nearly identical surface
690 > areas ($L_x$ and $L_y$) available to the metal and capping agent,
691 > while the extra solvent served mainly to lengthen the axis that was
692 > used to apply the thermal flux.  For a given value of the applied
693 > flux, the different $z$ length scale has only a weak effect on the
694 > computed conductivities (Table \ref{AuThiolHexaneUA}).
695  
696 < Furthermore, our observation validated previous
697 < results\cite{hase:2010} that the intramolecular heat transport of
698 < alkylthiols is highly effecient. As a combinational effects of these
699 < phenomena, butanethiol acts as a channel to expedite thermal transport
700 < process. The acoustic impedance mismatch between the metal and the
701 < liquid phase can be effectively reduced with the presence of suitable
702 < capping agents.
696 > \subsubsection{Effects of applied flux}
697 > The NIVS algorithm allows changes in both the sign and magnitude of
698 > the applied flux.  It is possible to reverse the direction of heat
699 > flow simply by changing the sign of the flux, and thermal gradients
700 > which would be difficult to obtain experimentally ($5$ K/\AA) can be
701 > easily simulated.  However, the magnitude of the applied flux is not
702 > arbitrary if one aims to obtain a stable and reliable thermal gradient.
703 > A temperature gradient can be lost in the noise if $|J_z|$ is too
704 > small, and excessive $|J_z|$ values can cause phase transitions if the
705 > extremes of the simulation cell become widely separated in
706 > temperature.  Also, if $|J_z|$ is too large for the bulk conductivity
707 > of the materials, the thermal gradient will never reach a stable
708 > state.  
709  
710 < \begin{figure}
711 < \includegraphics[width=\linewidth]{vibration}
712 < \caption{Vibrational spectra obtained for gold in different
713 <  environments.}
714 < \label{specAu}
715 < \end{figure}
710 > Within a reasonable range of $J_z$ values, we were able to study how
711 > $G$ changes as a function of this flux.  In what follows, we use
712 > positive $J_z$ values to denote the case where energy is being
713 > transferred by the method from the metal phase and into the liquid.
714 > The resulting gradient therefore has a higher temperature in the
715 > liquid phase.  Negative flux values reverse this transfer, and result
716 > in higher temperature metal phases.  The conductance measured under
717 > different applied $J_z$ values is listed in Tables
718 > \ref{AuThiolHexaneUA} and \ref{AuThiolToluene}. These results do not
719 > indicate that $G$ depends strongly on $J_z$ within this flux
720 > range. The linear response of flux to thermal gradient simplifies our
721 > investigations in that we can rely on $G$ measurement with only a
722 > small number $J_z$ values.  
723  
724 < [MAY ADD COMPARISON OF AU SLAB WIDTHS BUT NOT MUCH TO TALK ABOUT...]
724 > \begin{table*}
725 >  \begin{minipage}{\linewidth}
726 >    \begin{center}
727 >      \caption{In the hexane-solvated interfaces, the system size has
728 >        little effect on the calculated values for interfacial
729 >        conductance ($G$ and $G^\prime$), but the direction of heat
730 >        flow (i.e. the sign of $J_z$) can alter the average
731 >        temperature of the liquid phase and this can alter the
732 >        computed conductivity.}
733 >      
734 >      \begin{tabular}{ccccccc}
735 >        \hline\hline
736 >        $\langle T\rangle$ & $N_{hexane}$  & $\rho_{hexane}$ &
737 >        $J_z$ & $G$ & $G^\prime$ \\
738 >        (K) &  & (g/cm$^3$) & (GW/m$^2$) &
739 >        \multicolumn{2}{c}{(MW/m$^2$/K)} \\
740 >        \hline
741 >        200 & 266 &  0.672 & -0.96 & 102(3)    & 80.0(0.8) \\
742 >            & 200 &  0.688 &  0.96 & 125(16)   & 90.2(15)  \\
743 >            &     &        &  1.91 & 139(10)   & 101(10)   \\
744 >            &     &        &  2.83 & 141(6)    & 89.9(9.8) \\
745 >            & 166 &  0.681 &  0.97 & 141(30)   & 78(22)    \\
746 >            &     &        &  1.92 & 138(4)    & 98.9(9.5) \\
747 >        \hline
748 >        250 & 200 &  0.560 &  0.96 & 75(10)    & 61.8(7.3) \\
749 >            &     &        & -0.95 & 49.4(0.3) & 45.7(2.1) \\
750 >            & 166 &  0.569 &  0.97 & 80.3(0.6) & 67(11)    \\
751 >            &     &        &  1.44 & 76.2(5.0) & 64.8(3.8) \\
752 >            &     &        & -0.95 & 56.4(2.5) & 54.4(1.1) \\
753 >            &     &        & -1.85 & 47.8(1.1) & 53.5(1.5) \\
754 >        \hline\hline
755 >      \end{tabular}
756 >      \label{AuThiolHexaneUA}
757 >    \end{center}
758 >  \end{minipage}
759 > \end{table*}
760  
761 < \section{Conclusions}
762 < The NIVS algorithm we developed has been applied to simulations of
763 < Au-butanethiol surfaces with organic solvents. This algorithm allows
764 < effective unphysical thermal flux transferred between the metal and
765 < the liquid phase. With the flux applied, we were able to measure the
766 < corresponding thermal gradient and to obtain interfacial thermal
767 < conductivities. Under steady states, single trajectory simulation
768 < would be enough for accurate measurement. This would be advantageous
769 < compared to transient state simulations, which need multiple
928 < trajectories to produce reliable average results.
761 > The sign of $J_z$ is a different matter, however, as this can alter
762 > the temperature on the two sides of the interface. The average
763 > temperature values reported are for the entire system, and not for the
764 > liquid phase, so at a given $\langle T \rangle$, the system with
765 > positive $J_z$ has a warmer liquid phase.  This means that if the
766 > liquid carries thermal energy via convective transport, {\it positive}
767 > $J_z$ values will result in increased molecular motion on the liquid
768 > side of the interface, and this will increase the measured
769 > conductivity.
770  
771 < Our simulations have seen significant conductance enhancement with the
931 < presence of capping agent, compared to the bare gold / liquid
932 < interfaces. The acoustic impedance mismatch between the metal and the
933 < liquid phase is effectively eliminated by proper capping
934 < agent. Furthermore, the coverage precentage of the capping agent plays
935 < an important role in the interfacial thermal transport
936 < process. Moderately lower coverages allow higher contact between
937 < capping agent and solvent, and thus could further enhance the heat
938 < transfer process.
771 > \subsubsection{Effects due to average temperature}
772  
773 < Our measurement results, particularly of the UA models, agree with
774 < available experimental data. This indicates that our force field
775 < parameters have a nice description of the interactions between the
776 < particles at the interfaces. AA models tend to overestimate the
773 > We also studied the effect of average system temperature on the
774 > interfacial conductance.  The simulations are first equilibrated in
775 > the NPT ensemble at 1 atm.  The TraPPE-UA model for hexane tends to
776 > predict a lower boiling point (and liquid state density) than
777 > experiments.  This lower-density liquid phase leads to reduced contact
778 > between the hexane and butanethiol, and this accounts for our
779 > observation of lower conductance at higher temperatures as shown in
780 > Table \ref{AuThiolHexaneUA}.  In raising the average temperature from
781 > 200K to 250K, the density drop of $\sim$20\% in the solvent phase
782 > leads to a $\sim$40\% drop in the conductance.
783 >
784 > Similar behavior is observed in the TraPPE-UA model for toluene,
785 > although this model has better agreement with the experimental
786 > densities of toluene.  The expansion of the toluene liquid phase is
787 > not as significant as that of the hexane (8.3\% over 100K), and this
788 > limits the effect to $\sim$20\% drop in thermal conductivity  (Table
789 > \ref{AuThiolToluene}).
790 >
791 > Although we have not mapped out the behavior at a large number of
792 > temperatures, is clear that there will be a strong temperature
793 > dependence in the interfacial conductance when the physical properties
794 > of one side of the interface (notably the density) change rapidly as a
795 > function of temperature.
796 >
797 > \begin{table*}
798 >  \begin{minipage}{\linewidth}
799 >    \begin{center}
800 >      \caption{When toluene is the solvent, the interfacial thermal
801 >        conductivity is less sensitive to temperature, but again, the
802 >        direction of the heat flow can alter the solvent temperature
803 >        and can change the computed conductance values.}
804 >      
805 >      \begin{tabular}{ccccc}
806 >        \hline\hline
807 >        $\langle T\rangle$ & $\rho_{toluene}$ & $J_z$ & $G$ & $G^\prime$ \\
808 >        (K) & (g/cm$^3$) & (GW/m$^2$) & \multicolumn{2}{c}{(MW/m$^2$/K)} \\
809 >        \hline
810 >        200 & 0.933 &  2.15 & 204(12) & 113(12) \\
811 >            &       & -1.86 & 180(3)  & 135(21) \\
812 >            &       & -3.93 & 176(5)  & 113(12) \\
813 >        \hline
814 >        300 & 0.855 & -1.91 & 143(5)  & 125(2)  \\
815 >            &       & -4.19 & 135(9)  & 113(12) \\
816 >        \hline\hline
817 >      \end{tabular}
818 >      \label{AuThiolToluene}
819 >    \end{center}
820 >  \end{minipage}
821 > \end{table*}
822 >
823 > Besides the lower interfacial thermal conductance, surfaces at
824 > relatively high temperatures are susceptible to reconstructions,
825 > particularly when butanethiols fully cover the Au(111) surface. These
826 > reconstructions include surface Au atoms which migrate outward to the
827 > S atom layer, and butanethiol molecules which embed into the surface
828 > Au layer. The driving force for this behavior is the strong Au-S
829 > interactions which are modeled here with a deep Lennard-Jones
830 > potential. This phenomenon agrees with reconstructions that have been
831 > experimentally
832 > observed.\cite{doi:10.1021/j100035a033,doi:10.1021/la026493y}.  Vlugt
833 > {\it et al.} kept their Au(111) slab rigid so that their simulations
834 > could reach 300K without surface
835 > reconstructions.\cite{vlugt:cpc2007154} Since surface reconstructions
836 > blur the interface, the measurement of $G$ becomes more difficult to
837 > conduct at higher temperatures.  For this reason, most of our
838 > measurements are undertaken at $\langle T\rangle\sim$200K where
839 > reconstruction is minimized.
840 >
841 > However, when the surface is not completely covered by butanethiols,
842 > the simulated system appears to be more resistent to the
843 > reconstruction. Our Au / butanethiol / toluene system had the Au(111)
844 > surfaces 90\% covered by butanethiols, but did not see this above
845 > phenomena even at $\langle T\rangle\sim$300K.  That said, we did
846 > observe butanethiols migrating to neighboring three-fold sites during
847 > a simulation.  Since the interface persisted in these simulations, we
848 > were able to obtain $G$'s for these interfaces even at a relatively
849 > high temperature without being affected by surface reconstructions.
850 >
851 > \section{Discussion}
852 >
853 > The primary result of this work is that the capping agent acts as an
854 > efficient thermal coupler between solid and solvent phases.  One of
855 > the ways the capping agent can carry out this role is to down-shift
856 > between the phonon vibrations in the solid (which carry the heat from
857 > the gold) and the molecular vibrations in the liquid (which carry some
858 > of the heat in the solvent).
859 >
860 > To investigate the mechanism of interfacial thermal conductance, the
861 > vibrational power spectrum was computed. Power spectra were taken for
862 > individual components in different simulations. To obtain these
863 > spectra, simulations were run after equilibration in the
864 > microcanonical (NVE) ensemble and without a thermal
865 > gradient. Snapshots of configurations were collected at a frequency
866 > that is higher than that of the fastest vibrations occurring in the
867 > simulations. With these configurations, the velocity auto-correlation
868 > functions can be computed:
869 > \begin{equation}
870 > C_A (t) = \langle\vec{v}_A (t)\cdot\vec{v}_A (0)\rangle
871 > \label{vCorr}
872 > \end{equation}
873 > The power spectrum is constructed via a Fourier transform of the
874 > symmetrized velocity autocorrelation function,
875 > \begin{equation}
876 >  \hat{f}(\omega) =
877 >  \int_{-\infty}^{\infty} C_A (t) e^{-2\pi it\omega}\,dt
878 > \label{fourier}
879 > \end{equation}
880 >
881 > \subsection{The role of specific vibrations}
882 > The vibrational spectra for gold slabs in different environments are
883 > shown as in Figure \ref{specAu}. Regardless of the presence of
884 > solvent, the gold surfaces which are covered by butanethiol molecules
885 > exhibit an additional peak observed at a frequency of
886 > $\sim$165cm$^{-1}$.  We attribute this peak to the S-Au bonding
887 > vibration. This vibration enables efficient thermal coupling of the
888 > surface Au layer to the capping agents. Therefore, in our simulations,
889 > the Au / S interfaces do not appear to be the primary barrier to
890 > thermal transport when compared with the butanethiol / solvent
891 > interfaces. {\bf This confirms the results from Luo {\it et
892 >    al.}\cite{Luo20101}, which reported $G$ for Au-SAM junctions
893 >  generally twice larger than what we have computed for the
894 >  thiol-liquid interfaces.}
895 >
896 > \begin{figure}
897 > \includegraphics[width=\linewidth]{vibration}
898 > \caption{The vibrational power spectrum for thiol-capped gold has an
899 >  additional vibrational peak at $\sim $165cm$^{-1}$.  Bare gold
900 >  surfaces (both with and without a solvent over-layer) are missing
901 >  this peak.   A similar peak at  $\sim $165cm$^{-1}$ also appears in
902 >  the vibrational power spectrum for the butanethiol capping agents.}
903 > \label{specAu}
904 > \end{figure}
905 >
906 > Also in this figure, we show the vibrational power spectrum for the
907 > bound butanethiol molecules, which also exhibits the same
908 > $\sim$165cm$^{-1}$ peak.
909 >
910 > \subsection{Overlap of power spectra}
911 > A comparison of the results obtained from the two different organic
912 > solvents can also provide useful information of the interfacial
913 > thermal transport process.  In particular, the vibrational overlap
914 > between the butanethiol and the organic solvents suggests a highly
915 > efficient thermal exchange between these components.  Very high
916 > thermal conductivity was observed when AA models were used and C-H
917 > vibrations were treated classically.  The presence of extra degrees of
918 > freedom in the AA force field yields higher heat exchange rates
919 > between the two phases and results in a much higher conductivity than
920 > in the UA force field. {\bf Due to the classical models used, this
921 >  even includes those high frequency modes which should be unpopulated
922 >  at our relatively low temperatures. This artifact causes high
923 >  frequency vibrations accountable for thermal transport in classical
924 >  MD simulations.}
925 >
926 > The similarity in the vibrational modes available to solvent and
927 > capping agent can be reduced by deuterating one of the two components
928 > (Fig. \ref{aahxntln}).  Once either the hexanes or the butanethiols
929 > are deuterated, one can observe a significantly lower $G$ and
930 > $G^\prime$ values (Table \ref{modelTest}).
931 >
932 > \begin{figure}
933 > \includegraphics[width=\linewidth]{aahxntln}
934 > \caption{Spectra obtained for all-atom (AA) Au / butanethiol / solvent
935 >  systems. When butanethiol is deuterated (lower left), its
936 >  vibrational overlap with hexane decreases significantly.  Since
937 >  aromatic molecules and the butanethiol are vibrationally dissimilar,
938 >  the change is not as dramatic when toluene is the solvent (right).}
939 > \label{aahxntln}
940 > \end{figure}
941 >
942 > For the Au / butanethiol / toluene interfaces, having the AA
943 > butanethiol deuterated did not yield a significant change in the
944 > measured conductance. Compared to the C-H vibrational overlap between
945 > hexane and butanethiol, both of which have alkyl chains, the overlap
946 > between toluene and butanethiol is not as significant and thus does
947 > not contribute as much to the heat exchange process.
948 >
949 > Previous observations of Zhang {\it et al.}\cite{hase:2010} indicate
950 > that the {\it intra}molecular heat transport due to alkylthiols is
951 > highly efficient.  Combining our observations with those of Zhang {\it
952 >  et al.}, it appears that butanethiol acts as a channel to expedite
953 > heat flow from the gold surface and into the alkyl chain.  The
954 > acoustic impedance mismatch between the metal and the liquid phase can
955 > therefore be effectively reduced with the presence of suitable capping
956 > agents.
957 >
958 > Deuterated models in the UA force field did not decouple the thermal
959 > transport as well as in the AA force field.  The UA models, even
960 > though they have eliminated the high frequency C-H vibrational
961 > overlap, still have significant overlap in the lower-frequency
962 > portions of the infrared spectra (Figure \ref{uahxnua}).  Deuterating
963 > the UA models did not decouple the low frequency region enough to
964 > produce an observable difference for the results of $G$ (Table
965 > \ref{modelTest}).
966 >
967 > \begin{figure}
968 > \includegraphics[width=\linewidth]{uahxnua}
969 > \caption{Vibrational power spectra for UA models for the butanethiol
970 >  and hexane solvent (upper panel) show the high degree of overlap
971 >  between these two molecules, particularly at lower frequencies.
972 >  Deuterating a UA model for the solvent (lower panel) does not
973 >  decouple the two spectra to the same degree as in the AA force
974 >  field (see Fig \ref{aahxntln}).}
975 > \label{uahxnua}
976 > \end{figure}
977 >
978 > \section{Conclusions}
979 > The NIVS algorithm has been applied to simulations of
980 > butanethiol-capped Au(111) surfaces in the presence of organic
981 > solvents. This algorithm allows the application of unphysical thermal
982 > flux to transfer heat between the metal and the liquid phase. With the
983 > flux applied, we were able to measure the corresponding thermal
984 > gradients and to obtain interfacial thermal conductivities. Under
985 > steady states, 2-3 ns trajectory simulations are sufficient for
986 > computation of this quantity.
987 >
988 > Our simulations have seen significant conductance enhancement in the
989 > presence of capping agent, compared with the bare gold / liquid
990 > interfaces. The acoustic impedance mismatch between the metal and the
991 > liquid phase is effectively eliminated by a chemically-bonded capping
992 > agent. Furthermore, the coverage percentage of the capping agent plays
993 > an important role in the interfacial thermal transport
994 > process. Moderately low coverages allow higher contact between capping
995 > agent and solvent, and thus could further enhance the heat transfer
996 > process, giving a non-monotonic behavior of conductance with
997 > increasing coverage.
998 >
999 > Our results, particularly using the UA models, agree well with
1000 > available experimental data.  The AA models tend to overestimate the
1001   interfacial thermal conductance in that the classically treated C-H
1002 < vibration would be overly sampled. Compared to the AA models, the UA
1003 < models have higher computational efficiency with satisfactory
1004 < accuracy, and thus are preferable in interfacial thermal transport
1005 < modelings. Of the two definitions for $G$, the discrete form
1002 > vibrations become too easily populated. Compared to the AA models, the
1003 > UA models have higher computational efficiency with satisfactory
1004 > accuracy, and thus are preferable in modeling interfacial thermal
1005 > transport.
1006 >
1007 > Of the two definitions for $G$, the discrete form
1008   (Eq. \ref{discreteG}) was easier to use and gives out relatively
1009   consistent results, while the derivative form (Eq. \ref{derivativeG})
1010   is not as versatile. Although $G^\prime$ gives out comparable results
1011   and follows similar trend with $G$ when measuring close to fully
1012 < covered or bare surfaces, the spatial resolution of $T$ profile is
1013 < limited for accurate computation of derivatives data.
1012 > covered or bare surfaces, the spatial resolution of $T$ profile
1013 > required for the use of a derivative form is limited by the number of
1014 > bins and the sampling required to obtain thermal gradient information.
1015  
1016 < Vlugt {\it et al.} has investigated the surface thiol structures for
1017 < nanocrystal gold and pointed out that they differs from those of the
1018 < Au(111) surface\cite{landman:1998,vlugt:cpc2007154}. This difference
1019 < might lead to change of interfacial thermal transport behavior as
1020 < well. To investigate this problem, an effective means to introduce
1021 < thermal flux and measure the corresponding thermal gradient is
1022 < desirable for simulating structures with spherical symmetry.
1016 > Vlugt {\it et al.} have investigated the surface thiol structures for
1017 > nanocrystalline gold and pointed out that they differ from those of
1018 > the Au(111) surface.\cite{landman:1998,vlugt:cpc2007154} This
1019 > difference could also cause differences in the interfacial thermal
1020 > transport behavior. To investigate this problem, one would need an
1021 > effective method for applying thermal gradients in non-planar
1022 > (i.e. spherical) geometries.
1023  
1024   \section{Acknowledgments}
1025   Support for this project was provided by the National Science
1026   Foundation under grant CHE-0848243. Computational time was provided by
1027   the Center for Research Computing (CRC) at the University of Notre
1028 < Dame. \newpage
1028 > Dame.
1029  
1030 + \section{Supporting Information}
1031 + This information is available free of charge via the Internet at
1032 + http://pubs.acs.org.
1033 +
1034 + \newpage
1035 +
1036   \bibliography{interfacial}
1037  
1038   \end{doublespace}

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines