| 77 |
|
|
| 78 |
|
Metallic particles have also been proposed for use in highly efficient |
| 79 |
|
thermal-transfer fluids, although careful experiments by Eapen {\it et |
| 80 |
< |
al.} have shown that metal-particle-based ``nanofluids'' have |
| 81 |
< |
thermal conductivities that match Maxwell |
| 82 |
< |
predictions.\cite{Eapen:2007th} The likely cause of previously |
| 83 |
< |
reported non-Maxwell |
| 80 |
> |
al.} have shown that metal-particle-based nanofluids have thermal |
| 81 |
> |
conductivities that match Maxwell predictions.\cite{Eapen:2007th} The |
| 82 |
> |
likely cause of previously reported non-Maxwell |
| 83 |
|
behavior\cite{Eastman:2001wb,Keblinski:2002bx,Lee:1999ct,Xue:2003ya,Xue:2004oa} |
| 84 |
|
is percolation networks of nanoparticles exchanging energy via the |
| 85 |
< |
solvent,\cite{Eapen:2007mw} so it is vital to get a detailed molecular |
| 86 |
< |
picture of particle-ligand and particle-solvent interactions in order |
| 87 |
< |
to understand the thermal behavior of complex fluids. To date, there |
| 88 |
< |
have been few reported values (either from theory or experiment) for |
| 89 |
< |
$G$ for ligand-protected nanoparticles embedded in liquids, and there |
| 90 |
< |
is a significant gap in knowledge about how chemically distinct |
| 91 |
< |
ligands or protecting groups will affect heat transport from the |
| 92 |
< |
particles. |
| 85 |
> |
solvent,\cite{Eapen:2007mw} so it is important to get a detailed |
| 86 |
> |
molecular picture of particle-ligand and ligand-solvent interactions |
| 87 |
> |
in order to understand the thermal behavior of complex fluids. To |
| 88 |
> |
date, there have been few reported values (either from theory or |
| 89 |
> |
experiment) for $G$ for ligand-protected nanoparticles embedded in |
| 90 |
> |
liquids, and there is a significant gap in knowledge about how |
| 91 |
> |
chemically distinct ligands or protecting groups will affect heat |
| 92 |
> |
transport from the particles. |
| 93 |
|
|
| 94 |
|
Experimentally, the thermal properties of various nanostructured |
| 95 |
|
interfaces have been investigated by a number of groups. Cahill and |
| 96 |
< |
coworkers studied nanoscale thermal transport from metal |
| 97 |
< |
nanoparticle/fluid interfaces, to epitaxial TiN/single crystal oxide |
| 98 |
< |
interfaces, and hydrophilic and hydrophobic interfaces between water |
| 99 |
< |
and solids with different self-assembled |
| 96 |
> |
coworkers studied thermal transport from metal nanoparticle/fluid |
| 97 |
> |
interfaces, epitaxial TiN/single crystal oxide interfaces, and |
| 98 |
> |
hydrophilic and hydrophobic interfaces between water and solids with |
| 99 |
> |
different self-assembled |
| 100 |
|
monolayers.\cite{cahill:793,Wilson:2002uq,PhysRevB.67.054302,doi:10.1021/jp048375k,PhysRevLett.96.186101} |
| 101 |
|
Wang {\it et al.} studied heat transport through long-chain |
| 102 |
|
hydrocarbon monolayers on gold substrate at the individual molecular |
| 107 |
|
interface resistance of glass-embedded metal |
| 108 |
|
nanoparticles.\cite{PhysRevB.80.195406} Although interfaces are |
| 109 |
|
normally considered barriers for heat transport, Alper {\it et al.} |
| 110 |
< |
suggested that specific ligands (capping agents) could completely |
| 110 |
> |
have suggested that specific ligands (capping agents) could completely |
| 111 |
|
eliminate this barrier |
| 112 |
|
($G\rightarrow\infty$).\cite{doi:10.1021/la904855s} |
| 113 |
|
|
| 115 |
|
non-equilibrium molecular dynamics (RNEMD) to calculate the |
| 116 |
|
interfacial thermal conductance at a metal / organic solvent interface |
| 117 |
|
that had been chemically protected by butanethiolate groups. Our |
| 118 |
< |
calculations suggest an explanation for the very large thermal |
| 118 |
> |
calculations suggested an explanation for the very large thermal |
| 119 |
|
conductivity at alkanethiol-capped metal surfaces when compared with |
| 120 |
|
bare metal/solvent interfaces. Specifically, the chemical bond |
| 121 |
|
between the metal and the ligand introduces a vibrational overlap that |
| 123 |
|
the vibrational spectra (metal to ligand, ligand to solvent) provides |
| 124 |
|
a mechanism for rapid thermal transport across the interface. |
| 125 |
|
|
| 126 |
< |
One notable result of our previous work was the observation of |
| 127 |
< |
non-monotonic dependence of the thermal conductance on the coverage of |
| 128 |
< |
a metal surface by a chemical protecting group. Our explanation for |
| 129 |
< |
this behavior was that gaps in surface coverage allowed solvent to |
| 130 |
< |
penetrate close to the capping molecules that had been heated by the |
| 131 |
< |
metal surface, to absorb thermal energy from these molecules, and then |
| 132 |
< |
diffuse away. The effect of surface coverage was relatively difficult |
| 133 |
< |
to study as the individual protecting groups have lateral mobility on |
| 135 |
< |
the surface and can aggregate to form domains on the timescale of the |
| 136 |
< |
simulation. |
| 126 |
> |
A notable result of the previous simulations was the non-monotonic |
| 127 |
> |
dependence of $G$ on the fractional coverage of the metal surface by |
| 128 |
> |
the chemical protecting group. Gaps in surface coverage allow the |
| 129 |
> |
solvent molecules come into direct contact with ligands that had been |
| 130 |
> |
heated by the metal surface, absorb thermal energy from the ligands, |
| 131 |
> |
and then diffuse away. Quantifying the role of surface coverage is |
| 132 |
> |
difficult as the ligands have lateral mobility on the surface and can |
| 133 |
> |
aggregate to form domains on the timescale of the simulation. |
| 134 |
|
|
| 135 |
< |
To prevent lateral mobility of the surface ligands, the current work |
| 136 |
< |
involves mixed-chain monolayers in which the length mismatch between |
| 137 |
< |
long and short chains is sufficient to accomodate solvent |
| 138 |
< |
molecules. These complete (but mixed-chain) surfaces are significantly |
| 139 |
< |
less prone to surface domain formation, and allow us to further |
| 140 |
< |
investigate the mechanism of heat transport to the solvent. A thermal |
| 141 |
< |
flux is created using velocity shearing and scaling reverse |
| 142 |
< |
non-equilibrium molecular dynamics (VSS-RNEMD), and the resulting |
| 143 |
< |
temperature profiles are analyzed to yield information about the |
| 144 |
< |
interfacial thermal conductance. |
| 135 |
> |
To isolate this effect without worrying about lateral mobility of the |
| 136 |
> |
surface ligands, the current work involves mixed-chain-length |
| 137 |
> |
monolayers in which the length mismatch between long and short chains |
| 138 |
> |
is sufficient to accomodate solvent molecules. These completely |
| 139 |
> |
covered (but mixed-chain) surfaces are significantly less prone to |
| 140 |
> |
surface domain formation, and allow us to further investigate the |
| 141 |
> |
mechanism of heat transport to the solvent. A thermal flux is created |
| 142 |
> |
using velocity shearing and scaling reverse non-equilibrium molecular |
| 143 |
> |
dynamics (VSS-RNEMD), and the resulting temperature profiles are |
| 144 |
> |
analyzed to yield information about the interfacial thermal |
| 145 |
> |
conductance. |
| 146 |
|
|
| 147 |
|
|
| 148 |
|
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
| 169 |
|
{\it Reverse} Non-Equilibrium Molecular Dynamics (RNEMD) methods adopt |
| 170 |
|
a different approach in that an unphysical {\it flux} is imposed |
| 171 |
|
between different regions or ``slabs'' of the simulation |
| 172 |
< |
box.\cite{MullerPlathe:1997xw,ISI:000080382700030,Kuang:2010uq} The |
| 172 |
> |
box.\cite{MullerPlathe:1997xw,ISI:000080382700030,Kuang2010} The |
| 173 |
|
system responds by developing a temperature or momentum {\it gradient} |
| 174 |
|
between the two regions. Since the amount of the applied flux is known |
| 175 |
|
exactly, and the measurement of a gradient is generally less |
| 205 |
|
|
| 206 |
|
The most useful alternative RNEMD approach developed so far utilizes a |
| 207 |
|
series of simultaneous velocity shearing and scaling exchanges between |
| 208 |
< |
the two slabs.\cite{2012MolPh.110..691K} This method provides a set of |
| 208 |
> |
the two slabs.\cite{Kuang2012} This method provides a set of |
| 209 |
|
conservation constraints while simultaneously creating a desired flux |
| 210 |
|
between the two slabs. Satisfying the constraint equations ensures |
| 211 |
|
that the new configurations are sampled from the same NVE ensemble. |
| 284 |
|
kept to a minimum. This ability to generate simultaneous thermal and |
| 285 |
|
shear fluxes has been previously utilized to map out the shear |
| 286 |
|
viscosity of SPC/E water over a wide range of temperatures (90~K) with |
| 287 |
< |
a single 1 ns simulation.\cite{2012MolPh.110..691K} |
| 287 |
> |
a single 1 ns simulation.\cite{Kuang2012} |
| 288 |
|
|
| 289 |
|
\begin{figure} |
| 290 |
|
\includegraphics[width=\linewidth]{figures/rnemd} |
| 346 |
|
|
| 347 |
|
In the particular case we are studying here, there are two interfaces |
| 348 |
|
involved in the transfer of heat from the gold slab to the solvent: |
| 349 |
< |
the gold/thiolate interface and the thiolate/solvent interface. We |
| 349 |
> |
the metal/thiolate interface and the thiolate/solvent interface. We |
| 350 |
|
could treat the temperature on each side of an interface as discrete, |
| 351 |
|
making the interfacial conductance the inverse of the Kaptiza |
| 352 |
|
resistance, or $G = \frac{1}{R_k}$. To model the total conductance |
| 381 |
|
Our simulations include a number of chemically distinct components. |
| 382 |
|
Figure \ref{fig:structures} demonstrates the sites defined for both |
| 383 |
|
the {\it n}-hexane and alkanethiolate ligands present in our |
| 384 |
< |
simulations. Force field parameters are needed for interactions both |
| 387 |
< |
between the same type of particles and between particles of different |
| 388 |
< |
species. |
| 384 |
> |
simulations. |
| 385 |
|
|
| 386 |
|
\begin{figure} |
| 387 |
|
\includegraphics[width=\linewidth]{figures/structures} |
| 408 |
|
sites are located at the carbon centers for alkyl groups. Bonding |
| 409 |
|
interactions, including bond stretches and bends and torsions, were |
| 410 |
|
used for intra-molecular sites closer than 3 bonds. For non-bonded |
| 411 |
< |
interactions, Lennard-Jones potentials are used. We have previously |
| 411 |
> |
interactions, Lennard-Jones potentials were used. We have previously |
| 412 |
|
utilized both united atom (UA) and all-atom (AA) force fields for |
| 413 |
< |
thermal conductivity in previous work,\cite{} and since the united |
| 413 |
> |
thermal conductivity,\cite{kuang:AuThl,Kuang2012} and since the united |
| 414 |
|
atom force fields cannot populate the high-frequency modes that are |
| 415 |
|
present in AA force fields, they appear to work better for modeling |
| 416 |
|
thermal conductivity. The TraPPE-UA model for alkanes is known to |
| 417 |
|
predict a slightly lower boiling point than experimental values. This |
| 418 |
|
is one of the reasons we used a lower average temperature (200K) for |
| 419 |
< |
our simulations. |
| 419 |
> |
our simulations. |
| 420 |
|
|
| 421 |
|
The TraPPE-UA force field includes parameters for thiol |
| 422 |
|
molecules\cite{TraPPE-UA.thiols} which were used for the |
| 441 |
|
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
| 442 |
|
\subsection{Simulation Protocol} |
| 443 |
|
|
| 444 |
< |
We have implemented the VSS-RNEMD algorithm in OpenMD, our |
| 445 |
< |
group molecular dynamics code. A 1188 atom gold slab was |
| 446 |
< |
equilibrated at 1 atm and 200 K. The periodic box was then expanded |
| 447 |
< |
in the $z$-direction to expose two Au(111) faces on either side of the 11-atom thick slab. |
| 444 |
> |
We have implemented the VSS-RNEMD algorithm in OpenMD, our group |
| 445 |
> |
molecular dynamics code.\cite{openmd} An 1188 atom gold slab was |
| 446 |
> |
equilibrated at 1 atm and 200 K. The periodic box was then expanded |
| 447 |
> |
in the $z$-direction to expose two Au(111) faces on either side of the |
| 448 |
> |
11-atom thick slab. |
| 449 |
|
|
| 450 |
< |
A full monolayer of thiolates (1/3 the number of surface gold atoms) was placed on three-fold hollow sites on the gold interfaces. To efficiently test the effect of thiolate binding sites on the thermal conductance, all systems had one gold interface with thiolates placed only on fcc hollow sites and the other interface with thiolates only on hcp hollow sites. To test the effect of thiolate chain length on interfacial thermal conductance, full coverages of five chain lengths were tested: butanethiolate, hexanethiolate, octanethiolate, decanethiolate, and dodecanethiolate. To test the effect of mixed chain lengths, full coverages of butanethiolate/decanethiolate and butanethiolate/dodecanethiolate mixtures were created in short/long chain percentages of 25/75, 50/50, 62.5/37.5, 75/25, and 87.5/12.5. The short and long chains were placed on the surface hollow sites in a random configuration. |
| 450 |
> |
A full monolayer of thiolates (1/3 the number of surface gold atoms) |
| 451 |
> |
was placed on three-fold hollow sites on the gold interfaces. The |
| 452 |
> |
effect of thiolate binding sites on the thermal conductance was tested |
| 453 |
> |
by placing thiolates at both fcc and hcp hollow sites. No appreciable |
| 454 |
> |
difference in the temperature profile was noted due to the location of |
| 455 |
> |
thiolate binding. |
| 456 |
|
|
| 457 |
< |
The simulation box $z$-dimension was set to roughly double the length of the gold/thiolate block. Hexane solvent molecules were placed in the vacant portion of the box using the packmol algorithm. Hexane, a straight chain flexible alkane, is very structurally similar to the thiolate alkane tails; previous work has shown that UA models of hexane and butanethiolate have a high degree of vibrational overlap.\cite{Kuang2011} This overlap should provide a mechanism for thermal energy transfer from the thiolates to the solvent. |
| 457 |
> |
To test the role of thiolate chain length on interfacial thermal |
| 458 |
> |
conductance, full coverages of each of five chain lengths were tested: |
| 459 |
> |
butanethiolate (C$_4$), hexanethiolate (C$_6$), octanethiolate |
| 460 |
> |
(C$_8$), decanethiolate (C$_{10}$), and dodecanethiolate |
| 461 |
> |
(C$_{12}$). To test the effect of mixed chain lengths, full coverages |
| 462 |
> |
of C$_4$ / C$_{10}$ and C$_4$ / C$_{12}$ mixtures were created in |
| 463 |
> |
short/long chain percentages of 25/75, 50/50, 62.5/37.5, 75/25, and |
| 464 |
> |
87.5/12.5. The short and long chains were placed on the surface hollow |
| 465 |
> |
sites in a random configuration. |
| 466 |
|
|
| 467 |
< |
The system was equilibrated to 220 K in the NVT ensemble, allowing the thiolates and solvent to warm gradually. Pressure correction to 1 atm was done in an NPT ensemble that allowed expansion or contraction only in the z direction, so as not to disrupt the crystalline structure of the gold. The diagonal elements of the pressure tensor were monitored during the pressure correction step. If the xx and/or yy elements had a mean above zero throughout the simulation -- indicating residual pressure in the plane of the gold slab -- an additional short NPT equilibration step was performed allowing all box dimensions to change. Once the pressure was stable at 1 atm, a final NVT simulation was performed. All systems were equilibrated in the microcanonical (NVE) ensemble before proceeding with the VSS-RNEMD step. |
| 467 |
> |
The simulation box $z$-dimension was set to roughly double the length |
| 468 |
> |
of the gold/thiolate block. Hexane solvent molecules were placed in |
| 469 |
> |
the vacant portion of the box using the packmol algorithm. |
| 470 |
|
|
| 471 |
< |
A kinetic energy flux was applied using VSS-RNEMD in the NVE ensemble. The total simulation time was 3 nanoseconds, with velocity scaling occurring every 10 femtoseconds. The hot slab was centered in the gold and the cold slab was placed in the center of the solvent region. The average temperature was 220 K, with a temperature difference between the hot and cold slabs of approximately 30 K. The average temperature and kinetic energy flux were carefully selected with two considerations in mind: 1) if the cold bin gets too cold (below ~180 K) the solvent may freeze or undergo a glassy transition, and 2) due to the deep sulfur-gold potential well, sulfur atoms routinely embed into the gold slab, particularly at temperatures above 250 K. Simulation conditions were chosen to avoid both of these undesirable situations. A reversed flux direction resulted in frozen long chain thiolates and solvent too near its boiling point. |
| 471 |
> |
The system was equilibrated to 220 K in the canonical (NVT) ensemble, |
| 472 |
> |
allowing the thiolates and solvent to warm gradually. Pressure |
| 473 |
> |
correction to 1 atm was done using an isobaric-isothermal (NPT) |
| 474 |
> |
integrator that allowed expansion or contraction only in the $z$ |
| 475 |
> |
direction, maintaining the crystalline structure of the gold as close |
| 476 |
> |
to the bulk result as possible. The diagonal elements of the pressure |
| 477 |
> |
tensor were monitored during the pressure equilibration stage. If the |
| 478 |
> |
$xx$ and/or $yy$ elements had a mean above zero throughout the |
| 479 |
> |
simulation -- indicating residual surface tension in the plane of the |
| 480 |
> |
gold slab -- an additional short NPT equilibration step was performed |
| 481 |
> |
allowing all box dimensions to change. Once the pressure was stable |
| 482 |
> |
at 1 atm, a final equilibration stage was performed at constant |
| 483 |
> |
temperature. All systems were equilibrated in the microcanonical (NVE) |
| 484 |
> |
ensemble before proceeding with the VSS-RNEMD and data collection |
| 485 |
> |
stages. |
| 486 |
|
|
| 487 |
< |
A temperature profile of the system was created by dividing the box into $\sim$ 3 \AA \, bins along the z axis and recording the average temperature of each bin. |
| 487 |
> |
A kinetic energy flux was applied using VSS-RNEMD. The total |
| 488 |
> |
simulation time was 3 nanoseconds, with velocity scaling occurring |
| 489 |
> |
every 10 femtoseconds. The hot slab was centered in the gold and the |
| 490 |
> |
cold slab was placed in the center of the solvent region. The entire |
| 491 |
> |
system has a (time-averaged) temperature of 220 K, with a temperature |
| 492 |
> |
difference between the hot and cold slabs of approximately 30 K. The |
| 493 |
> |
average temperature and kinetic energy flux were selected to prevent |
| 494 |
> |
solvent freezing (or glass formation) and to not allow the thiolates |
| 495 |
> |
to bury in the gold slab. The Au-S interaction has a deep potential |
| 496 |
> |
energy well, which allows sulfur atoms to embed into the gold slab at |
| 497 |
> |
temperatures above 250 K. Simulation conditions were chosen to avoid |
| 498 |
> |
both of these undesirable situations. |
| 499 |
|
|
| 500 |
+ |
Temperature profiles of the system were created by dividing the box |
| 501 |
+ |
into $\sim$ 3 \AA \, bins along the $z$ axis and recording the average |
| 502 |
+ |
temperature of each bin. |
| 503 |
|
|
| 504 |
+ |
|
| 505 |
|
|
| 506 |
|
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
| 507 |
|
% **RESULTS** |
| 508 |
|
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
| 509 |
|
\section{Results} |
| 510 |
|
|
| 511 |
+ |
The solvent, hexane, a straight chain flexible alkane, is structurally |
| 512 |
+ |
similar to the thiolate alkane tails, and previous work has shown that |
| 513 |
+ |
UA models of hexane and butanethiolate have a high degree of |
| 514 |
+ |
vibrational overlap.\cite{kuang:AuThl} This overlap provides a |
| 515 |
+ |
mechanism for thermal energy transfer from the thiolates to the |
| 516 |
+ |
solvent. |
| 517 |
|
|
| 518 |
|
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
| 519 |
|
% CHAIN LENGTH |
