trunk/chuckDissertation/nanoglass.tex

%!TEX root = /Users/charles/Documents/chuckDissertation/dissertation.tex
\chapter{\label{chap:nanoglass}GLASS FORMATION IN METALLIC NANOPARTICLES}

\section{INTRODUCTION}

Excitation of the plasmon resonance in metallic nanoparticles has
attracted enormous interest in the past several years. This is partly
due to the location of the plasmon band in the near IR for particles
in a wide range of sizes and geometries. Living tissue is nearly
transparent in the near IR, and for this reason, there is an
unrealized potential for metallic nanoparticles to be used in both
diagnostic and therapeutic settings.\cite{West:2003fk,Hu:2006lr} One
of the side effects of absorption of laser radiation at these
frequencies is the rapid (sub-picosecond) heating of the electronic
degrees of freedom in the metal. This hot electron gas quickly
transfers heat to the phonon modes of the particle, resulting in a
rapid heating of the lattice of the metal particles.  Since metallic
nanoparticles have a large surface area to volume ratio, many of the
metal atoms are at surface locations and experience relatively weak
bonding. This is observable in a lowering of the melting temperatures
of these particles when compared with bulk metallic
samples.\cite{Buffat:1976yq,Dick:2002qy} One of the side effects of
the excitation of small metallic nanoparticles at the plasmon
resonance is the facile creation of liquid metal
droplets.\cite{Mafune01,HartlandG.V._jp0276092,Link:2000lr,Plech:2003yq,plech:195423,Plech:2007rt}

Much of the experimental work on this subject has been carried out in
the Hartland, El-Sayed and Plech
groups.\cite{HartlandG.V._jp0276092,Hodak:2000rb,Hartland:2003lr,Petrova:2007qy,Link:2000lr,plech:195423,Plech:2007rt}
These experiments mostly use the technique of time-resolved optical
pump-probe spectroscopy, where a pump laser pulse serves to excite
conduction band electrons in the nanoparticle and a following probe
laser pulse allows observation of the time evolution of the
electron-phonon coupling. Hu and Hartland have observed a direct
relation between the size of the nanoparticle and the observed cooling
rate using such pump-probe techniques.\cite{Hu:2004lr} Plech {\it et
al.} have use pulsed x-ray scattering as a probe to directly access
changes to atomic structure following pump
excitation.\cite{plech:195423} They further determined that heat
transfer in nanoparticles to the surrounding solvent is goverened by
interfacial dynamics and not the thermal transport properties of the
solvent.  This is in agreement with Cahill,\cite{Wilson:2002uq}
but opposite to the conclusions in Reference \citen{Hu:2004lr}.

Since these experiments are carried out in condensed phase
surroundings, the large surface area to volume ratio makes the heat
transfer to the surrounding solvent a relatively rapid process. In our
recent simulation study of the laser excitation of gold
nanoparticles,\cite{VardemanC.F._jp051575r} we observed that the
cooling rate for these particles (10$^{11}$-10$^{12}$ K/s) is in
excess of the cooling rate required for glass formation in bulk
metallic alloys.\cite{Greer:1995qy} Given this fact, it may be
possible to use laser excitation to melt, alloy and quench metallic
nanoparticles in order to form glassy nanobeads.

To study whether or not glass nanobead formation is feasible, we have
chosen the bimetallic alloy of Silver (60\%) and Copper (40\%) as a
model system because it is an experimentally known glass former, and
has been used previously as a theoretical model for glassy
dynamics.\cite{Vardeman-II:2001jn} The Hume-Rothery rules suggest that
alloys composed of Copper and Silver should be miscible in the solid
state, because their lattice constants are within 15\% of each
another.\cite{Kittel:1996fk} Experimentally, however Ag-Cu alloys are
a well-known exception to this rule and are only miscible in the
liquid state given equilibrium conditions.\cite{Massalski:1986rt}
Below the eutectic temperature of 779 $^\circ$C and composition
(60.1\% Ag, 39.9\% Cu), the solid alloys of Ag and Cu will phase
separate into Ag and Cu rich $\alpha$ and $\beta$ phases,
respectively.\cite{Banhart:1992sv,Ma:2005fk} This behavior is due to a
positive heat of mixing in both the solid and liquid phases. For the
one-to-one composition fcc solid solution, $\Delta H_{\rm mix}$ is on
the order of +6~kJ/mole.\cite{Ma:2005fk} Non-equilibrium solid
solutions may be formed by undercooling, and under these conditions, a
compositionally-disordered $\gamma$ fcc phase can be
formed.\cite{najafabadi:3144}

Metastable alloys composed of Ag-Cu were first reported by Duwez in
1960 and were created by using a ``splat quenching'' technique in
which a liquid droplet is propelled by a shock wave against a cooled
metallic target.\cite{duwez:1136} Because of the small positive
$\Delta H_{\rm mix}$, supersaturated crystalline solutions are
typically obtained rather than an amorphous phase. Higher $\Delta
H_{\rm mix}$ systems, such as Ag-Ni, are immiscible even in liquid
states, but they tend to form metastable alloys much more readily than
Ag-Cu. If present, the amorphous Ag-Cu phase is usually seen as the
minority phase in most experiments. Because of this unique
crystalline-amorphous behavior, the Ag-Cu system has been widely
studied. Methods for creating such bulk phase structures include splat
quenching, vapor deposition, ion beam mixing and mechanical
alloying. Both structural \cite{sheng:184203} and
dynamic\cite{Vardeman-II:2001jn} computational studies have also been
performed on this system.

Although bulk Ag-Cu alloys have been studied widely, this alloy has
been mostly overlooked in nanoscale materials. The literature on
alloyed metallic nanoparticles has dealt with the Ag-Au system, which
has the useful property of being miscible on both solid and liquid
phases. Nanoparticles of another miscible system, Au-Cu, have been
successfully constructed using techniques such as laser
ablation,\cite{Malyavantham:2004cu} and the synthetic reduction of
metal ions in solution.\cite{Kim:2003lv} Laser induced alloying has
been used as a technique for creating Au-Ag alloy particles from
core-shell particles.\cite{Hartland:2003lr} To date, attempts at
creating Ag-Cu nanoparticles have used ion implantation to embed
nanoparticles in a glass matrix.\cite{De:1996ta,Magruder:1994rg} These
attempts have been largely unsuccessful in producing mixed alloy
nanoparticles, and instead produce phase segregated or core-shell
structures.

One of the more successful attempts at creating intermixed Ag-Cu
nanoparticles used alternate pulsed laser ablation and deposition in
an amorphous Al$_2$O$_3$ matrix.\cite{gonzalo:5163} Surface plasmon
resonance (SPR) of bimetallic core-shell structures typically show two
distinct resonance peaks where mixed particles show a single shifted
and broadened resonance.\cite{Hodak:2000rb} The SPR for pure silver
occurs at 400 nm and for copper at 570 nm.\cite{HengleinA._jp992950g}
On Al$_2$O$_3$ films, these resonances move to 424 nm and 572 nm for
the pure metals. For bimetallic nanoparticles with 40\% Ag an
absorption peak is seen between 400-550 nm. With increasing Ag
content, the SPR shifts towards the blue, with the peaks nearly
coincident at a composition of 57\% Ag. Gonzalo {\it et al.} cited the
existence of a single broad resonance peak as evidence of an alloyed
particle rather than a phase segregated system.  However, spectroscopy
may not be able to tell the difference between alloyed particles and
mixtures of segregated particles.  High-resolution electron microscopy
has so far been unable to determine whether the mixed nanoparticles
were an amorphous phase or a supersaturated crystalline phase.

Characterization of glassy behavior by molecular dynamics simulations
is typically done using dynamic measurements such as the mean squared
displacement, $\langle r^2(t) \rangle$. Liquids exhibit a mean squared
displacement that is linear in time (at long times). Glassy materials
deviate significantly from this linear behavior at intermediate times,
entering a sub-linear regime with a return to linear behavior in the
infinite time limit.\cite{Kob:1999fk} However, diffusion in
nanoparticles differs significantly from the bulk in that atoms are
confined to a roughly spherical volume and cannot explore any region
larger than the particle radius ($R$). In these confined geometries,
$\langle r^2(t) \rangle$ approaches a limiting value of
$3R^2/40$.\cite{ShibataT._ja026764r} This limits the utility of
dynamical measures of glass formation when studying nanoparticles.

However, glassy materials exhibit strong icosahedral ordering among
nearest-neghbors (in contrast with crystalline and liquid-like
configurations). Local icosahedral structures are the
three-dimensional equivalent of covering a two-dimensional plane with
5-sided tiles; they cannot be used to tile space in a periodic
fashion, and are therefore an indicator of non-periodic packing in
amorphous solids. Steinhart {\it et al.} defined an orientational bond
order parameter that is sensitive to icosahedral
ordering.\cite{Steinhardt:1983mo} This bond order parameter can
therefore be used to characterize glass formation in liquid and solid
solutions.\cite{wolde:9932}

Theoretical molecular dynamics studies have been performed on the
formation of amorphous single component nanoclusters of either
gold,\cite{Chen:2004ec,Cleveland:1997jb,Cleveland:1997gu} or
nickel,\cite{Gafner:2004bg,Qi:2001nn} by rapid cooling($\thicksim
10^{12}-10^{13}$ K/s) from a liquid state. All of these studies found
icosahedral ordering in the resulting structures produced by this
rapid cooling which can be evidence of the formation of an amorphous
structure.\cite{Strandburg:1992qy} The nearest neighbor information
was obtained from pair correlation functions, common neighbor analysis
and bond order parameters.\cite{Steinhardt:1983mo} It should be noted
that these studies used single component systems with cooling rates
that are only obtainable in computer simulations and particle sizes
less than 20\AA. Single component systems are known to form amorphous
states in small clusters,\cite{Breaux:rz} but do not generally form
amorphous structures in bulk materials.

Since the nanoscale Ag-Cu alloy has been largely unexplored, many
interesting questions remain about the formation and properties of
such a system. Does the large surface area to volume ratio aid Ag-Cu
nanoparticles in rapid cooling and formation of an amorphous state?
Nanoparticles have been shown to have a size dependent melting
transition ($T_m$),\cite{Buffat:1976yq,Dick:2002qy} so we might expect
a similar trend to follow for the glass transition temperature
($T_g$). By analogy, bulk metallic glasses exhibit a correlation
between $T_m$ and $T_g$ although such dependence is difficult to
establish because of the dependence of $T_g$ on cooling rate and the
process by which the glass is formed.\cite{Wang:2003fk} It has also
been demonstrated that there is a finite size effect depressing $T_g$
in polymer glasses in confined geometries.\cite{Alcoutlabi:2005kx}

In the sections below, we describe our modeling of the laser
excitation and subsequent cooling of the particles {\it in silico} to
mimic real experimental conditions. The simulation parameters have
been tuned to the degree possible to match experimental conditions,
and we discusss both the icosahedral ordering in the system, as well
as the clustering of icosahedral centers that we observed.

\section{COMPUTATIONAL METHODOLOGY} 
\label{nanoglass:sec:details} 

\subsection{INITIAL GEOMETRIES AND HEATING}


Cu-core / Ag-shell and random alloy structures were constructed on an
underlying FCC lattice (4.09 {\AA}) at the bulk eutectic composition
$\mathrm{Ag}_6\mathrm{Cu}_4$.  Both initial geometries were considered
although experimental results suggest that the random structure is the
most likely structure to be found following
synthesis.\cite{Jiang:2005lr,gonzalo:5163} Three different sizes of
nanoparticles corresponding to a 20 \AA radius (2382 atoms), 30 {\AA}
radius (6603 atoms) and 40 {\AA} radius (15683 atoms) were
constructed.  These initial structures were relaxed to their
equilibrium structures at 20 K for 20 ps and again at 300 K for 100 ps
sampling from a Maxwell-Boltzmann distribution at each
temperature. All simulations were conducted using the {\sc oopse}
molecular dynamics package.\cite{Meineke:2004uq}

To mimic the effects of the heating due to laser irradiation, the
particles were allowed to melt by sampling velocities from the Maxwell
Boltzmann distribution at a temperature of 900 K.  The particles were
run under microcanonical simulation conditions for 1 ns of simualtion
time prior to studying the effects of heat transfer to the solvent.
In all cases, center of mass translational and rotational motion of
the particles were set to zero before any data collection was
undertaken.  Structural features (pair distribution functions) were
used to verify that the particles were indeed liquid droplets before
cooling simulations took place.

\subsection{MODELING RANDOM ALLOY AND CORE SHELL PARTICLES IN SOLUTION PHASE ENVIRONMENTS}


To approximate the effects of rapid heat transfer to the solvent
following a heating at the plasmon resonance, we utilized a
methodology in which atoms contained in the outer $4$ {\AA} radius of
the nanoparticle evolved under Langevin Dynamics,
\begin{equation}
m \frac{\partial^2 \vec{x}}{\partial t^2} = F_\textrm{sys}(\vec{x}(t))
- 6 \pi a \eta \vec{v}(t)  + F_\textrm{ran}
\label{eq:langevin}
\end{equation}
with a solvent friction ($\eta$) approximating the contribution from
the solvent and capping agent.  Atoms located in the interior of the
nanoparticle evolved under Newtonian dynamics.  The set-up of our
simulations is nearly identical with the ``stochastic boundary
molecular dynamics'' ({\sc sbmd}) method that has seen wide use in the
protein simulation
community.\cite{BROOKS:1985kx,BROOKS:1983uq,BRUNGER:1984fj} A sketch
of this setup can be found in Fig. \ref{fig:langevinSketch}.  In
Eq. (\ref{eq:langevin}) the frictional forces of a spherical atom
of radius $a$ depend on the solvent viscosity.  The random forces are
usually taken as gaussian random variables with zero mean and a
variance tied to the solvent viscosity and temperature,
\begin{equation}
\langle F_\textrm{ran}(t) \cdot F_\textrm{ran} (t')
\rangle = 2 k_B T (6 \pi \eta a) \delta(t - t')
\label{eq:stochastic}
\end{equation}
Due to the presence of the capping agent and the lack of details about
the atomic-scale interactions between the metallic atoms and the
solvent, the effective viscosity is a essentially a free parameter
that must be tuned to give experimentally relevant simulations.
\begin{figure}[htbp]
\centering
\includegraphics[width=5in]{images/stochbound.pdf}
\caption{Methodology used to mimic the experimental cooling conditions
of a hot nanoparticle surrounded by a solvent.  Atoms in the core of
the particle evolved under Newtonian dynamics, while atoms that were
in the outer skin of the particle evolved under Langevin dynamics.
The radius of the spherical region operating under Newtonian dynamics,
$r_\textrm{Newton}$ was set to be 4 {\AA} smaller than the original
radius ($R$) of the liquid droplet.}
\label{fig:langevinSketch}
\end{figure}

The viscosity ($\eta$) can be tuned by comparing the cooling rate that
a set of nanoparticles experience with the known cooling rates for
similar particles obtained via the laser heating experiments.
Essentially, we tune the solvent viscosity until the thermal decay
profile matches a heat-transfer model using reasonable values for the
interfacial conductance and the thermal conductivity of the solvent.

Cooling rates for the experimentally-observed nanoparticles were
calculated from the heat transfer equations for a spherical particle
embedded in a ambient medium that allows for diffusive heat transport.
Following Plech {\it et al.},\cite{plech:195423} we use a heat
transfer model that consists of two coupled differential equations
in the Laplace domain,
\begin{eqnarray}
Mc_{P}\cdot(s\cdot T_{p}(s)-T_{0})+4\pi R^{2} G\cdot(T_{p}(s)-T_{f}(r=R,s)=0\\
\left(\frac{\partial}{\partial r} T_{f}(r,s)\right)_{r=R} +
\frac{G}{K}(T_{p}(s)-T_{f}(r,s) = 0 
\label{eq:heateqn} 
\end{eqnarray}
where $s$ is the time-conjugate variable in Laplace space. The
variables in these equations describe a nanoparticle of radius $R$,
mass $M$, and specific heat $c_{p}$ at an initial temperature
$T_0$. The surrounding solvent has a thermal profile $T_f(r,t)$,
thermal conductivity $K$, density $\rho$, and specific heat $c$. $G$
is the interfacial conductance between the nanoparticle and the
surrounding solvent, and contains information about heat transfer to
the capping agent as well as the direct metal-to-solvent heat loss.
The temperature of the nanoparticle as a function of time can then
obtained by the inverse Laplace transform,
\begin{equation}
T_{p}(t)=\frac{2 k R^2 g^2
T_0}{\pi}\int_{0}^{\infty}\frac{\exp(-\kappa u^2
t/R^2)u^2}{(u^2(1 + R g) - k R g)^2+(u^3 - k R g u)^2}\mathrm{d}u. 
\label{eq:laplacetransform}
\end{equation}
For simplicity, we have introduced the thermal diffusivity $\kappa =
K/(\rho c)$,  and defined $k=4\pi R^3 \rho c /(M c_p)$ and $g = G/K$ in
Eq. (\ref{eq:laplacetransform}).

Eq. (\ref{eq:laplacetransform}) was solved numerically for the Ag-Cu
system using mole-fraction weighted values for $c_p$ and $\rho_p$ of
0.295 $(\mathrm{J g^{-1} K^{-1}})$ and $9.826\times 10^6$ $(\mathrm{g
m^{-3}})$ respectively. Since most of the laser excitation experiments
have been done in aqueous solutions, parameters used for the fluid are
$K$ of $0.6$ $(\mathrm{Wm^{-1}K^{-1}})$, $\rho$ of $1.0\times10^6$
$(\mathrm{g m^{-3}})$ and $c$ of $4.184$ $(\mathrm{J g^{-1} K^{-1}})$.

Values for the interfacial conductance have been determined by a
number of groups for similar nanoparticles and range from a low
$87.5\times 10^{6} (\mathrm{Wm^{-2}K^{-1}})$ to $130\times 10^{6}
(\mathrm{Wm^{-2}K^{-1}})$.\cite{Wilson:2002uq} Wilson {\it
et al.}  worked with Au, Pt, and AuPd nanoparticles and obtained an
estimate for the interfacial conductance of $G=130
(\mathrm{Wm^{-2}K^{-1}})$.\cite{Wilson:2002uq} Similarly, Plech {\it
et al.}  reported a value for the interfacial conductance of $G=105\pm
15 (\mathrm{Wm^{-2}K^{-1}})$ for Au nanoparticles.\cite{plech:195423}

We conducted our simulations at both ends of the range of
experimentally-determined values for the interfacial conductance.
This allows us to observe both the slowest and fastest heat transfers
from the nanoparticle to the solvent that are consistent with
experimental observations.  For the slowest heat transfer, a value for
G of $87.5\times 10^{6}$ $(\mathrm{Wm^{-2}K^{-1}})$ was used and for
the fastest heat transfer, a value of $117\times 10^{6}$
$(\mathrm{Wm^{-2}K^{-1}})$ was used.  Based on calculations we have
done using raw data from the Hartland group's thermal half-time
experiments on Au nanospheres,\cite{HuM._jp020581+} the true G values
are probably in the faster regime: $117\times 10^{6}$
$(\mathrm{Wm^{-2}K^{-1}})$.

The rate of cooling for the nanoparticles in a molecular dynamics
simulation can then be tuned by changing the effective solvent
viscosity ($\eta$) until the nanoparticle cooling rate matches the
cooling rate described by the heat-transfer Eq.
(\ref{eq:heateqn}). The effective solvent viscosity (in Pa s) for a G
of $87.5\times 10^{6}$ $(\mathrm{Wm^{-2}K^{-1}})$ is $4.2 \times
10^{-6}$, $5.0 \times 10^{-6}$, and
$5.5 \times 10^{-6}$ for 20 {\AA}, 30 {\AA}, and 40 {\AA} particles, respectively. The
effective solvent viscosity (again in Pa s) for an interfacial
conductance of $117\times 10^{6}$ $(\mathrm{Wm^{-2}K^{-1}})$ is $5.7
\times 10^{-6}$, $7.2 \times 10^{-6}$, and $7.5 \times 10^{-6}$
for 20 {\AA}, 30 {\AA} and 40 {\AA} particles.  These viscosities are
essentially gas-phase values, a fact which is consistent with the
initial temperatures of the particles being well into the
super-critical region for the aqueous environment.  Gas bubble
generation has also been seen experimentally around gold nanoparticles
in water.\cite{kotaidis:184702} Instead of a single value for the
effective viscosity, a time-dependent parameter might be a better
mimic of the cooling vapor layer that surrounds the hot particles.
This may also be a contributing factor to the size-dependence of the
effective viscosities in our simulations.

Cooling traces for each particle size are presented in
Fig. \ref{fig:images_cooling_plot}. It should be noted that the
Langevin thermostat produces cooling curves that are consistent with
Newtonian (single-exponential) cooling, which cannot match the cooling
profiles from Eq. (\ref{eq:laplacetransform}) exactly. Fitting the
Langevin cooling profiles to a single-exponential produces
$\tau=25.576$ ps, $\tau=43.786$ ps, and $\tau=56.621$ ps for the 20,
30 and 40 {\AA} nanoparticles and a G of $87.5\times 10^{6}$
$(\mathrm{Wm^{-2}K^{-1}})$.  For comparison's sake, similar
single-exponential fits with an interfacial conductance of G of
$117\times 10^{6}$ $(\mathrm{Wm^{-2}K^{-1}})$ produced a $\tau=13.391$
ps, $\tau=30.426$ ps, $\tau=43.857$ ps for the 20, 30 and 40 {\AA}
nanoparticles.

\begin{figure}[htbp]
\centering
\includegraphics[width=5in]{images/cooling_plot.pdf}
\caption{Thermal cooling curves obtained from the inverse Laplace
transform heat model in Eq. (\ref{eq:laplacetransform}) (solid line) as
well as from molecular dynamics simulations (circles).  Effective
solvent viscosities of 4.2-7.5 $\times 10^{-6}$ Pa s (depending on the
radius of the particle) give the best fit to the experimental cooling
curves.  This viscosity suggests that the nanoparticles in these
experiments are surrounded by a vapor layer (which is a reasonable
assumptions given the initial temperatures of the particles).  }
\label{fig:images_cooling_plot}
\end{figure}

\subsection{POTENIALS FOR CLASSICAL SIMULATIONS OF BIMETALLIC NANOPARTICLES}

Several different potential models have been developed that reasonably
describe interactions in transition metals. In particular, the
Embedded Atom Model ({\sc eam})~\cite{PhysRevB.33.7983} and
Sutton-Chen ({\sc sc})~\cite{Chen90} potential have been used to study
a wide range of phenomena in both bulk materials and
nanoparticles.\cite{Vardeman-II:2001jn,ShibataT._ja026764r,Sankaranarayanan:2005lr,Chui:2003fk,Wang:2005qy,Medasani:2007uq} Both
potentials are based on a model of a metal which treats the nuclei and
core electrons as pseudo-atoms embedded in the electron density due to
the valence electrons on all of the other atoms in the system. The
{\sc sc} potential has a simple form that closely resembles that of
the ubiquitous Lennard Jones potential,
\begin{equation}
\label{eq:SCP1} 
U_{tot}=\sum _{i}\left[ \frac{1}{2}\sum _{j\neq i}D_{ij}V^{pair}_{ij}(r_{ij})-c_{i}D_{ii}\sqrt{\rho_{i}}\right] , 
\end{equation}
where $V^{pair}_{ij}$ and $\rho_{i}$ are given by 
\begin{equation}
\label{eq:SCP2} 
V^{pair}_{ij}(r)=\left( \frac{\alpha_{ij}}{r_{ij}}\right)^{n_{ij}}, \rho_{i}=\sum_{j\neq i}\left( \frac{\alpha_{ij}}{r_{ij}}\right) ^{m_{ij}}. 
\end{equation}
$V^{pair}_{ij}$ is a repulsive pairwise potential that accounts for
interactions between the pseudoatom cores. The $\sqrt{\rho_i}$ term in
Eq. (\ref{eq:SCP1}) is an attractive many-body potential that models
the interactions between the valence electrons and the cores of the
pseudo-atoms. $D_{ij}$, $D_{ii}$ set the appropriate overall energy
scale, $c_i$ scales the attractive portion of the potential relative
to the repulsive interaction and $\alpha_{ij}$ is a length parameter
that assures a dimensionless form for $\rho$. These parameters are
tuned to various experimental properties such as the density, cohesive
energy, and elastic moduli for FCC transition metals. The quantum
Sutton-Chen ({\sc q-sc}) formulation matches these properties while
including zero-point quantum corrections for different transition
metals.\cite{PhysRevB.59.3527} This particular parametarization has
been shown to reproduce the experimentally available heat of mixing
data for both FCC solid solutions of Ag-Cu and the high-temperature
liquid.\cite{sheng:184203} In contrast, the {\sc eam} potential does
not reproduce the experimentally observed heat of mixing for the
liquid alloy.\cite{MURRAY:1984lr} In this work, we have utilized the
{\sc q-sc} formulation for our potential energies and forces.
Combination rules for the alloy were taken to be the arithmetic
average of the atomic parameters with the exception of $c_i$ since its
values is only dependent on the identity of the atom where the density
is evaluated.  For the {\sc q-sc} potential, cutoff distances are
traditionally taken to be $2\alpha_{ij}$ and include up to the sixth
coordination shell in FCC metals.

%\subsection{Sampling single-temperature configurations from a cooling
%trajectory} 

To better understand the structural changes occurring in the
nanoparticles throughout the cooling trajectory, configurations were
sampled at regular intervals during the cooling trajectory. These
configurations were then allowed to evolve under NVE dynamics to
sample from the proper distribution in phase space. Fig.
\ref{fig:images_cooling_time_traces} illustrates this sampling.


\begin{figure}[htbp]
        \centering
                \includegraphics[height=3in]{images/cooling_time_traces.pdf}
        \caption{Illustrative cooling profile for the 40 {\AA}
nanoparticle evolving under stochastic boundary conditions
corresponding to $G=$$87.5\times 10^{6}$
$(\mathrm{Wm^{-2}K^{-1}})$. At temperatures along the cooling
trajectory, configurations were sampled and allowed to evolve in the
NVE ensemble. These subsequent trajectories were analyzed for
structural features associated with bulk glass formation.}
        \label{fig:images_cooling_time_traces}
\end{figure}


\begin{figure}[htbp]
\centering
\includegraphics[width=5in]{images/cross_section_array.jpg}
\caption{Cutaway views of 30 \AA\ Ag-Cu nanoparticle structures for
random alloy (top) and Cu (core) / Ag (shell) initial conditions
(bottom).  Shown from left to right are the crystalline, liquid
droplet, and final glassy bead configurations.}
\label{fig:cross_sections}
\end{figure}

\section{ANALYSIS}

Frank first proposed local icosahedral ordering of atoms as an
explanation for supercooled atomic (specifically metallic) liquids,
and further showed that a 13-atom icosahedral cluster has a 8.4\%
higher binding energy the either a face centered cubic ({\sc fcc}) or
hexagonal close-packed ({\sc hcp}) crystal structures.\cite{19521106}
Icosahedra also have six five-fold symmetry axes that cannot be
extended indefinitely in three dimensions; long-range translational
order is therefore incommensurate with local icosahedral ordering.
This does not preclude icosahedral clusters from possessing long-range
{\it orientational} order. The ``frustrated'' packing of these
icosahedral structures into dense clusters has been proposed as a
model for glass formation.\cite{19871127} The size of the icosahedral
clusters is thought to increase until frustration prevents any further
growth.\cite{HOARE:1976fk} Molecular dynamics simulations of a
two-component Lennard-Jones glass showed that clusters of face-sharing
icosahedra are distributed throughout the
material.\cite{PhysRevLett.60.2295} Simulations of freezing of single
component metalic nanoclusters have shown a tendency for icosohedral
structure formation particularly at the surfaces of these
clusters.\cite{Gafner:2004bg,PhysRevLett.89.275502,Chen:2004ec}
Experimentally, the splitting (or shoulder) on the second peak of the
X-ray structure factor in binary metallic glasses has been attributed
to the formation of tetrahedra that share faces of adjoining
icosahedra.\cite{Waal:1995lr}

Various structural probes have been used to characterize structural
order in molecular systems including: common neighbor analysis,
Voronoi tesselations, and orientational bond-order
parameters.\cite{HoneycuttJ.Dana_j100303a014,Iwamatsu:2007lr,hsu:4974,nose:1803}
The method that has been used most extensively for determining local
and extended orientational symmetry in condensed phases is the
bond-orientational analysis formulated by Steinhart {\it et
al.}\cite{Steinhardt:1983mo} In this model, a set of spherical
harmonics is associated with each of the near neighbors of a central
atom.  Neighbors (or ``bonds'') are defined as having a distance from
the central atom that is within the first peak in the radial
distribution function. The spherical harmonic between a central atom
$i$ and a neighboring atom $j$ is
\begin{equation}
 Y_{lm}\left(\vec{r}_{ij}\right)=Y_{lm}\left(\theta(\vec{r}_{ij}),\phi(\vec{r}_{ij})\right)
\label{eq:spharm}
\end{equation}
where, $\{ Y_{lm}(\theta,\phi)\}$ are spherical harmonics, and
$\theta(\vec{r})$ and $\phi(\vec{r})$ are the polar and azimuthal
angles made by the bond vector $\vec{r}$ with respect to a reference
coordinate system. We chose for simplicity the origin as defined by
the coordinates for our nanoparticle. (Only even-$l$ spherical
harmonics are considered since permutation of a pair of identical
particles should not affect the bond-order parameter.) The local
environment surrounding atom $i$ can be defined by
the average over all neighbors, $N_b(i)$, surrounding that atom,
\begin{equation}
\bar{q}_{lm}(i) = \frac{1}{N_{b}(i)}\sum_{j=1}^{N_b(i)} Y_{lm}(\vec{r}_{ij}).
\label{eq:local_avg_bo}
\end{equation}
We can further define a global average orientational-bond order over
all $\bar{q}_{lm}$ by calculating the average of $\bar{q}_{lm}(i)$
over all $N$ particles 
\begin{equation}
\bar{Q}_{lm} = \frac{\sum^{N}_{i=1}N_{b}(i)\bar{q}_{lm}(i)}{\sum^{N}_{i=1}N_{b}(i)}
\label{eq:sys_avg_bo}
\end{equation}
The $\bar{Q}_{lm}$ contained in Eq. (\ref{eq:sys_avg_bo}) is not
necessarily invariant under rotations of the arbitrary reference
coordinate system.  Second- and third-order rotationally invariant
combinations, $Q_l$ and $W_l$, can be taken by summing over $m$ values
of $\bar{Q}_{lm}$,
\begin{equation}
Q_{l} = \left( \frac{4\pi}{2 l+1} \sum^{l}_{m=-l} \left| \bar{Q}_{lm} \right|^2\right)^{1/2}
\label{eq:sec_ord_inv}
\end{equation}
and
\begin{equation} 
\hat{W}_{l} = \frac{W_{l}}{\left( \sum^{l}_{m=-l} \left| \bar{Q}_{lm} \right|^2\right)^{3/2}}
\label{eq:third_ord_inv}
\end{equation}
where
\begin{equation} 
W_{l} = \sum_{\substack{m_1,m_2,m_3 \\m_1 + m_2 + m_3 = 0}} \left( \begin{array}{ccc} l & l & l \\ m_1 & m_2 & m_3 \end{array}\right) \\ \times \bar{Q}_{lm_1}\bar{Q}_{lm_2}\bar{Q}_{lm_3}.
\label{eq:third_inv}
\end{equation}
The factor in parentheses in Eq. (\ref{eq:third_inv}) is the Wigner-3$j$
symbol, and the sum is over all valid ($|m| \leq l$) values of $m_1$,
$m_2$, and $m_3$ which sum to zero.

\begin{table}
\caption{Values of bond orientational order parameters for
simple structures cooresponding to $l=4$ and $l=6$ spherical harmonic
functions.\cite{wolde:9932} $Q_4$ and $\hat{W}_4$ have values for {\it
individual} icosahedral clusters, but these values are not invariant
under rotations of the reference coordinate systems.  Similar behavior
is observed in the bond-orientational order parameters for individual
liquid-like structures.}
\begin{center}
\begin{tabular}{ccccc}
\hline
\hline
   & $Q_4$ & $Q_6$ & $\hat{W}_4$ & $\hat{W}_6$\\

fcc & 0.191 & 0.575 & -0.159 & -0.013\\

hcp & 0.097 & 0.485 & 0.134 & -0.012\\

bcc & 0.036 & 0.511 & 0.159 & 0.013\\

sc & 0.764 & 0.354 & 0.159 & 0.013\\

Icosahedral & - & 0.663 & - & -0.170\\

(liquid) & - & - & - &  -\\
\hline
\end{tabular}
\end{center}
\label{table:bopval}
\end{table} 
 
For ideal face-centered cubic ({\sc fcc}), body-centered cubic ({\sc
bcc}) and simple cubic ({\sc sc}) as well as hexagonally close-packed
({\sc hcp}) structures, these rotationally invariant bond order
parameters have fixed values independent of the choice of coordinate
reference frames.  For ideal icosahedral structures, the $l=6$
invariants, $Q_6$ and $\hat{W}_6$ are also independent of the
coordinate system.  $Q_4$ and $\hat{W}_4$ will have non-vanishing
values for {\it individual} icosahedral clusters, but these values are
not invariant under rotations of the reference coordinate systems.
Similar behavior is observed in the bond-orientational order
parameters for individual liquid-like structures.  Additionally, both
$Q_6$ and $\hat{W}_6$ are thought to have extreme values for the
icosahedral clusters.\cite{Steinhardt:1983mo} This makes the $l=6$
bond-orientational order parameters particularly useful in identifying
the extent of local icosahedral ordering in condensed phases.  For
example, a local structure which exhibits $\hat{W}_6$ values near
-0.17 is easily identified as an icosahedral cluster and cannot be
mistaken for distorted cubic or liquid-like structures.

One may use these bond orientational order parameters as an averaged
property to obtain the extent of icosahedral ordering in a supercooled
liquid or cluster.  It is also possible to accumulate information
about the {\it distributions} of local bond orientational order
parameters, $p(\hat{W}_6)$ and $p(Q_6)$, which provide information
about individual atomic sites that are central to local icosahedral
structures.

The distributions of atomic $Q_6$ and $\hat{W}_6$ values are plotted
as a function of temperature for our nanoparticles in Fig.
\ref{fig:q6} and \ref{fig:w6}.   At high temperatures, the
distributions are unstructured and are broadly distributed across the
entire range of values.  As the particles are cooled, however, there
is a dramatic increase in the fraction of atomic sites which have
local icosahedral ordering around them.  (This corresponds to the
sharp peak appearing in Fig. \ref{fig:w6} at $\hat{W}_6=-0.17$ and
to the broad shoulder appearing in Fig. \ref{fig:q6} at $Q_6 =
0.663$.)

\begin{figure}[htbp]
\centering
\includegraphics[width=5in]{images/w6_stacked_plot.pdf}
\caption{Distributions of the bond orientational order parameter
($\hat{W}_6$) at different temperatures.  The upper, middle, and lower
panels are for 20, 30, and 40 \AA\ particles, respectively.  The
left-hand column used cooling rates commensurate with a low
interfacial conductance ($87.5 \times 10^{6}$
$\mathrm{Wm^{-2}K^{-1}}$), while the right-hand column used a more
physically reasonable value of $117 \times 10^{6}$
$\mathrm{Wm^{-2}K^{-1}}$.  The peak at $\hat{W}_6 \approx -0.17$ is
due to local icosahedral structures.  The different curves in each of
the panels indicate the distribution of $\hat{W}_6$ values for samples
taken at different times along the cooling trajectory.  The initial
and final temperatures (in K) are indicated on the plots adjacent to
their respective distributions.}
\label{fig:w6}
\end{figure}

\begin{figure}[htbp]
\centering
\includegraphics[width=5in]{images/q6_stacked_plot.pdf}
\caption{Distributions of the bond orientational order parameter
($Q_6$) at different temperatures.  The curves in the six panels in
this figure were computed at identical conditions to the same panels in
figure \ref{fig:w6}.}
\label{fig:q6}
\end{figure}

The probability distributions of local order can be used to generate
free energy surfaces using the local orientational ordering as a
reaction coordinate.  By making the simple statistical equivalence
between the free energy and the probabilities of occupying certain
states,
\begin{equation}
g(\hat{W}_6) = - k_B T \ln p(\hat{W}_6),
\end{equation}
we can obtain a sequence of free energy surfaces (as a function of
temperature) for the local ordering around central atoms within our
particles.  Free energy surfaces for the 40 \AA\ particle at a range
of temperatures are shown in Fig. \ref{fig:freeEnergy}.  Note that
at all temperatures, the liquid-like structures are global minima on
the free energy surface, while the local icosahedra appear as local
minima once the temperature has fallen below 528 K.  As the
temperature falls, it is possible for substructures to become trapped
in the local icosahedral well, and if the cooling is rapid enough,
this trapping leads to vitrification.  A similar analysis of the free
energy surface for orientational order in bulk glass formers can be
found in the work of van~Duijneveldt and
Frenkel.\cite{duijneveldt:4655}


\begin{figure}[htbp]
\centering
\includegraphics[width=5in]{images/freeEnergyVsW6.pdf}
\caption{Free energy as a function of the orientational order
parameter ($\hat{W}_6$) for 40 {\AA} bimetallic nanoparticles as they
are cooled from 902 K to 310 K.  As the particles cool below 528 K, a
local minimum in the free energy surface appears near the perfect
icosahedral ordering ($\hat{W}_6 = -0.17$).  At all temperatures,
liquid-like structures are a global minimum on the free energy
surface, but if the cooling rate is fast enough, substructures
may become trapped with local icosahedral order, leading to the
formation of a glass.}
\label{fig:freeEnergy}
\end{figure}

We have also calculated the fraction of atomic centers which have
strong local icosahedral order:
\begin{equation}
f_\textrm{icos} = \int_{-\infty}^{w_i} p(\hat{W}_6) d \hat{W}_6
\label{eq:ficos}
\end{equation}
where $w_i$ is a cutoff value in $\hat{W}_6$ for atomic centers that
are displaying icosahedral environments.  We have chosen a (somewhat
arbitrary) value of $w_i= -0.15$ for the purposes of this work.  A
plot of $f_\textrm{icos}(T)$ as a function of temperature of the
particles is given in Fig. \ref{fig:ficos}.  As the particles cool,
the fraction of local icosahedral ordering rises smoothly to a plateau
value.  The smaller particles (particularly the ones that were cooled
in a higher viscosity solvent) show a slightly larger tendency towards
icosahedral ordering.

\begin{figure}[htbp]
\centering
\includegraphics[width=5in]{images/fraction_icos.pdf}
\caption{Temperautre dependence of the fraction of atoms with local
icosahedral ordering, $f_\textrm{icos}(T)$ for 20, 30, and 40 \AA\
particles cooled at two different values of the interfacial
conductance.} 
\label{fig:ficos}
\end{figure}

Since we have atomic-level resolution of the local bond-orientational
ordering information, we can also look at the local ordering as a
function of the identities of the central atoms.  In figure
\ref{fig:AgVsCu} we display the distributions of $\hat{W}_6$ values
for both the silver and copper atoms, and we note a strong
predilection for the copper atoms to be central to icosahedra.  This
is probably due to local packing competition of the larger silver
atoms around the copper, which would tend to favor icosahedral
structures over the more densely packed cubic structures.

\begin{figure}[htbp]
\centering
\includegraphics[width=5in]{images/w6_stacked_bytype_plot.pdf}
\caption{Distributions of the bond orientational order parameter
($\hat{W}_6$) for the two different elements present in the
nanoparticles.  This distribution was taken from the fully-cooled 40
\AA\ nanoparticle.  Local icosahedral ordering around copper atoms is
much more prevalent than around silver atoms.}
\label{fig:AgVsCu}
\end{figure}

The locations of these icosahedral centers are not uniformly
distrubted throughout the particles.  In Fig. \ref{fig:icoscluster}
we show snapshots of the centers of the local icosahedra (i.e. any
atom which exhibits a local bond orientational order parameter
$\hat{W}_6 < -0.15$).  At high temperatures, the icosahedral centers
are transitory, existing only for a few fs before being reabsorbed
into the liquid droplet.  As the particle cools, these centers become
fixed at certain locations, and additional icosahedra develop
throughout the particle, clustering around the sites where the
structures originated.  There is a strong preference for icosahedral
ordering near the surface of the particles.  Identification of these
structures by the type of atom shows that the silver-centered
icosahedra are evident only at the surface of the particles.

\begin{figure}[htbp]
\centering
\begin{tabular}{c c c}
\includegraphics[width=2.1in]{images/Cu_Ag_random_30A_liq_icosonly.pdf}
\includegraphics[width=2.1in]{images/Cu_Ag_random_30A__0007_icosonly.pdf}
\includegraphics[width=2.1in]{images/Cu_Ag_random_30A_glass_icosonly.pdf}
\end{tabular}
\caption{Centers of local icosahedral order ($\hat{W}_6<0.15$) at 900 
K, 471 K and 315 K for the 30 \AA\ nanoparticle cooled with an
interfacial conductance $G = 87.5 \times 10^{6}$
$\mathrm{Wm^{-2}K^{-1}}$. Silver atoms (blue) exhibit icosahedral
order at the surface of the nanoparticle while copper icosahedral
centers (green) are distributed throughout the nanoparticle.  The
icosahedral centers appear to cluster together and these clusters
increase in size with decreasing temperature.}
\label{fig:icoscluster}
\end{figure}

In contrast with the silver ordering behavior, the copper atoms which
have local icosahedral ordering are distributed more evenly throughout
the nanoparticles.  Fig. \ref{fig:Surface} shows this tendency as a
function of distance from the center of the nanoparticle.  Silver,
since it has a lower surface free energy than copper, tends to coat
the skins of the mixed particles.\cite{Zhu:1997lr} This is true even
for bimetallic particles that have been prepared in the Ag (core) / Cu
(shell) configuration.  Upon forming a liquid droplet, approximately 1
monolayer of Ag atoms will rise to the surface of the particles.  This
can be seen visually in Fig. \ref{fig:cross_sections} as well as in
the density plots in the bottom panel of Fig. \ref{fig:Surface}.
This observation is consistent with previous experimental and
theoretical studies on bimetallic alloys composed of noble
metals.\cite{MainardiD.S._la0014306,HuangS.-P._jp0204206,Ramirez-Caballero:2006lr}
Bond order parameters for surface atoms are averaged only over the
neighboring atoms, so packing constraints that may prevent icosahedral
ordering around silver in the bulk are removed near the surface.  It
would certainly be interesting to see if the relative tendency of
silver and copper to form local icosahedral structures in a bulk glass
differs from our observations on nanoparticles.

\begin{figure}[htbp]
\centering
\includegraphics[width=5in]{images/dens_fracr_stacked_plot.pdf}
\caption{Appearance of icosahedral clusters around central silver atoms
is largely due to the presence of these silver atoms at or near the
surface of the nanoparticle. The upper panel shows the fraction of
icosahedral atoms ($f_\textrm{icos}(r)$ for each of the two metallic
atoms as a function of distance from the center of the nanoparticle
($r$).  The lower panel shows the radial density of the two
constituent metals (relative to the overall density of the
nanoparticle).  Icosahedral clustering around copper atoms are more
evenly distributed throughout the particle, while icosahedral
clustering around silver is largely confined to the silver atoms at
the surface.}
\label{fig:Surface}
\end{figure}

The methods used by Sheng, He, and Ma to estimate the glass transition
temperature, $T_g$, in bulk Ag-Cu alloys involve finding
discontinuities in the slope of the average atomic volume, $\langle V
\rangle / N$, or enthalpy when plotted against the temperature of the
alloy.  They obtained a bulk glass transition temperature, $T_g$ = 510
K for a quenching rate of $2.5 \times 10^{13}$ K/s.

For simulations of nanoparticles, there is no periodic box, and
therefore, no easy way to compute the volume exactly.  Instead, we
estimate the volume of our nanoparticles using Barber {\it et al.}'s
very fast quickhull algorithm to obtain the convex hull for the
collection of 3-d coordinates of all of atoms at each point in
time.~\cite{barber96quickhull,qhull} The convex hull is the smallest convex
polyhedron which includes all of the atoms, so the volume of this
polyhedron is an excellent estimate of the volume of the nanoparticle.
This method of estimating the volume will be problematic if the
nanoparticle breaks into pieces (i.e. if the bounding surface becomes
concave), but for the relatively short trajectories used in this
study, it provides an excellent measure of particle volume as a
function of time (and temperature).

Using the discontinuity in the slope of the average atomic volume
vs. temperature, we arrive at an estimate of $T_g$ that is
approximately 488 K.  We note that this temperature is somewhat below
the onset of icosahedral ordering exhibited in the bond orientational
order parameters. It appears that icosahedral ordering sets in while
the system is still somewhat fluid, and is locked in place once the
temperature falls below $T_g$.  We did not observe any dependence of
our estimates for $T_g$ on either the nanoparticle size or the value
of the interfacial conductance.  However, the cooling rates and size
ranges we utilized are all sampled from a relatively narrow range, and
it is possible that much larger particles would have substantially
different values for $T_g$.  Our estimates for the glass transition
temperatures for all three particle sizes and both interfacial
conductance values are shown in table \ref{table:Tg}.

\begin{table}
\caption{Estimates of the glass transition temperatures $T_g$ for
three different sizes of bimetallic Ag$_6$Cu$_4$ nanoparticles cooled
under two different values of the interfacial conductance, $G$.} 
\begin{center}
\begin{tabular}{ccccc}
\hline
\hline
Radius (\AA\ ) & Interfacial conductance & Effective cooling rate
(K/s $\times 10^{13}$) &  & $T_g$ (K) \\
20 & 87.5 & 2.4 & 477 \\
20 & 117  & 4.5 & 502 \\
30 & 87.5 & 1.3 & 491 \\
30 & 117  & 1.9 & 493 \\
40 & 87.5 & 1.0 & 476 \\
40 & 117  & 1.3 & 487 \\
\hline
\end{tabular}
\end{center}
\label{table:Tg}
\end{table}

\section{CONCLUSIONS}
\label{metglass:sec:conclusion}

Our heat-transfer calculations have utilized the best current
estimates of the interfacial heat transfer coefficient (G) from recent
experiments.  Using reasonable values for thermal conductivity in both
the metallic particle and the surrounding solvent, we have obtained
cooling rates for laser-heated nanoparticles that are in excess of
10$^{13}$ K / s.  To test whether or not this cooling rate can form
glassy nanoparticles, we have performed a mixed molecular dynamics
simulation in which the atoms in contact with the solvent or capping
agent are evolved under Langevin dynamics while the remaining atoms
are evolved under Newtonian dynamics.  The effective solvent viscosity
($\eta$) is a free parameter which we have tuned so that the particles
in the simulation follow the same cooling curve as their experimental
counterparts.  From the local icosahedral ordering around the atoms in
the nanoparticles (particularly Copper atoms), we deduce that it is
likely that glassy nanobeads are created via laser heating of
bimetallic nanoparticles, particularly when the initial temperature of
the particles approaches the melting temperature of the bulk metal
alloy.

Improvements to our calculations would require: 1) explicit treatment
of the capping agent and solvent, 2) another radial region to handle
the heat transfer to the solvent vapor layer that is likely to form
immediately surrounding the hot
particle,\cite{Hu:2004lr,kotaidis:184702} and 3) larger particles in
the size range most easily studied via laser heating experiments.

The local icosahedral ordering we observed in these bimetallic
particles is centered almost completely around the copper atoms, and
this is likely due to the size mismatch leading to a more efficient
packing of 5-membered rings of silver around a central copper atom.
This size mismatch should be reflected in bulk calculations, and work
is ongoing in our lab to confirm this observation in bulk
glass-formers.

The physical properties (bulk modulus, frequency of the breathing
mode, and density) of glassy nanobeads should be somewhat different
from their crystalline counterparts.  However, observation of these
differences may require single-particle resolution of the ultrafast
vibrational spectrum of one particle both before and after the
crystallite has been converted into a glassy bead.
Revision:	3483
Committed:	Tue Jan 13 14:39:50 2009 UTC (16 years, 10 months ago) by chuckv
Content type:	application/x-tex
File size:	48122 byte(s)
Log Message:	Chuck's dissertation for PhD Jan 2009