Table Of ContentNumerical study of anharmonic vibrational decay in amorphous and paracrystalline
silicon
3
0 Jaroslav Fabian
0 Institute for Theoretical Physics, Karl-Franzens University, Universit¨atsplatz 5, 8010 Graz, Austria
2
n Joseph L. Feldman and C. Stephen Hellberg
a Center for Computational Materials Science, Naval Research Laboratory, Washington DC 20375-5345
J
4 S. M. Nakhmanson
1 Department of Physics, North Carolina State University, Raleigh, NC 27695-8202
] The anharmonicdecay rates ofatomic vibrations in amorphoussilicon (a-Si)and paracrystalline
n
silicon (p-Si), containing small crystalline grains embedded in a disordered matrix, are calculated
n
- using realistic structural models. The models are 1000-atom four-coordinated networks relaxed
s to a local minimum of the Stillinger-Weber interatomic potential. The vibrational decay rates
i
d are calculated numerically by perturbation theory, taking into account cubic anharmonicity as the
. perturbation. The vibrational lifetimes for a-Si are found to be on picosecond time scales, in
t
a agreement with the previous perturbative and classical molecular dynamics calculations on a 216-
m atommodel. Thecalculateddecayratesforp-Siaresimilartothoseofa-Si. Nomodesinp-Sireside
entirelyonthecrystallinecluster,decoupledfromtheamorphousmatrix. Thelocalizedmodeswith
-
d thelargest(upto59%)weightontheclusterdecayprimarilytotwodiffusons. Thenumericalresults
n are discussed in relation to a recent suggestion by van der Voort et al. [Phys. Rev. B 62, 8072
o (2000)] that long vibrational relaxation inferred experimentally may be due to possible crystalline
c nanostructuresin some typesof a-Si.
[
PACSnumbers: 63.50.+x,65.60.+a
1
v
2
I. INTRODUCTION silicon contain nanoscale regions with correlated (if
4
not ordered) atoms, which, through enhanced size
2
1 The pioneering experiments by Dijkhuis and quantization and localization of vibrational frequencies,
0 coworkers1,2,3,4 explored transient dynamics of ex- inhibit anharmonic decay.
3 cited vibrational modes in a topologically disordered Theseexperimentalresultsareatoddswiththeknown
0
material–hydrogenated amorphous silicon. In these theories of anharmonic vibrational decay in disordered
/
t experiments nonequilibrium vibrational states were gen- materials.6,7,8 Intheirsocalled“fracton”model,Alexan-
a
m erated during relaxation and recombination of optically deretal.9assumedthatthemajorityofvibrationalstates
excited electrons, and monitored with a probe laser in disordered systems are localized. This seemed to ex-
-
d (anti-Stokes Raman spectroscopy) for transient behav- plain the above experimental findings since the anhar-
n ior. The experimental results are surprising: Scholten monic decay could be drastically reduced by the ex-
o et al.1,2 found that at low temperatures (2 K) and for tremely small likelihood of the overlapbetween three lo-
c vibrational frequencies greater than 10 meV (maximum calizedmodes.10,11Thatthesmallprobabilityoftheover-
:
v frequency in a-Si is about 70 meV) vibrations decay lap between three localized modes inhibits vibrational
i on time scales of tens of nanoseconds. Furthermore, decay was disputed by Fabian and Allen6 who put for-
X
the higher the vibrational frequency, the slower is the ward a probabilistic scaling argument that the interac-
r
a decay rate. In contrast, phonons in crystalline silicon tion between three localizedmodes would in fact be cru-
decay on time scales of tens of picoseconds5 and the cial for the anharmonic decay and cannot be neglected.
decay rates increase with increasing frequency. The Fabian later demonstrated12 the scaling argument on a
results of Scholten et al. were further confirmed by one dimensional anharmonic chain with random spring
van der Voort et al.,3 who suggested that the long constants, and similar conclusions were reached recently
lifetimes are due to the microstructure of amorphous by Leitner in a study of heat flow in a one dimen-
silicon. This suggestion was tested by van der Voort sionalglass13andvibrationalenergytransferinhelicesof
et al.4 by measuring the vibrational decay rates of a myoglobin.14 Thus the fracton model, even if true in its
mixed amorphous-nanocrystalline silicon, which was premise of localizationof the majority of the vibrational
an amorphous hydrogenated silicon with a sizable modes, does not explain the experiment. We note, how-
fraction of nanocrystallites (with the diameter of 1 5 ever,that even the premise of the model is questionable,
nm). Even the mixed sample displayed nanoseco−nd as it is in sharp contrast to what is found in finite-size
vibrational lifetimes, although the lifetimes appeared to realistic models of glasses, which normally exhibit local-
decrease with increasing frequency. A hypothesis was ization only in a small part of the spectrum.
put forward4 that the measured types of amorphous Numerical calculations of vibrational decay in glasses
2
have been performed both by evaluating a perturbation to nanocrystalline Si sample of Ref. 4 (24% crystalline
formula6 andby classicalmoleculardynamics.7,8 Pertur- fractionand4.5nmaveragegraindiameter,versus 10%
≈
bation theory was appliedto the problemof anharmonic and 1 nm in the model) nor can it account for various
decay in glasses by Fabian and Allen6 who computed othertopologicaldefectspresentinrealmaterial. Still, if
the decay rates for a 216 atom model of amorphous Si. van der Voort’s supposition were correct, we would ob-
The decay rates were found to be fractions ofmeV (that serve inhibited decay rates of the modes which would be
is, lifetimes are picoseconds), in general increasing with predominantly localized on the crystalline cluster in the
increasing frequency. The anharmonic lifetimes of local- model. However,we do not find any modes localized ex-
ized modes were similar to those of the extended modes, clusivelyonthecluster: oneofthemostcluster-localized
even in the case of a model alloy SixGe1−x, where lo- modes has only 59% weight on the cluster, and is there-
calized modes span more than a half of the spectrum fore well coupled to the disordered matrix. It is not sur-
andthe overlapbetweenlocalizedstatesbecomes impor- prisingthatsuchmodeshavedecayratessimilartoother
tant.6 BickhamandFeldman7 reportedvibrationaldecay localized states.
rates for selected modes of 216 and4096 atommodels of We remark that the names amorphous and paracrys-
amorphousSi, usingclassicalmoleculardynamics. Their tallineinreferencetoourmodelsareamatteroftermino-
results agree with the perturbative calculation, though logicalconveniencerather than anattempt in classifying
the computed decay rates are somewhat greater due to real materials. We refer to a-Si as describing a homoge-
the fact that molecular dynamics takes into account all neous,continuousrandomnetworkofsiliconatoms,while
the anharmonic interaction, while the perturbative cal- p-Simodelsaresuchnetworksfilledwithcrystallineclus-
culation in Ref. 6 only cubic anharmonicity. In the cal- ters. Realmaterials—whicharenormallytermedamor-
culation of Bickham and Feldman, a chosen vibrational phous silicon — are likely of the p-Si type, containing
modewasgivenagreaterthanaveragekineticenergyand nanoscale crystallites with a distribution of sizes.16
wasallowedto equilibratewhile keepingthe overalltem- In the following we first introduce the structuralmod-
perature constant. From the decay of the kinetic energy els of a-Si and p-Si and their harmonic vibrationalprop-
in time, the mode decay rate was obtained. While the erties,thendiscusstheperturbativecalculationofanhar-
advantage of molecular dynamics over perturbation the- monic decay rates and present the results for the 1000-
ory in calculating vibrational decay rates is that the full atom models of a-Si and p-Si. Finally, we discuss our
anharmonicinteractionisconsidered,thedisadvantageis results with respect to the experiment.
that the classical dynamics does not capture accurately
the low temperature decay rates (for example, the rates
computedbyclassicalmoleculardynamicsvanishatzero
II. MODELS
temperature,7whileinrealitytheyarefiniteduetoquan-
tum effects6).
Both models employed in our studies were created
The purpose of this paper is twofold: (i) to extend with similar techniques: the homogeneous model for a-
the previous numerical studies of perturbative anhar- Si was constructed using the WWW method17 and the
monic decay in homogeneous amorphous silicon (a-Si) paracrystallinemodelwithavariation15 oftheBarkema-
to a larger system, and (ii) to calculate vibrational de- Mousseau method.18 (For a recent review of modeling
cay rates for a model of amorphous silicon — paracrys- continuous random networks see Ref. 19.) The for-
talline silicon (p-Si) — that includes nanocrystallites. mer model was also studied20 for its harmonic proper-
Thelargersystemisa1000-atommodelofa-Si,prepared ties within the framework of the Stillinger-Weber (SW)
similarly to the previously used 216-atom model.6 The potential21 prior to the present work. Apart from the
calculated decay rates display smaller statistical fluctu- computational efficiency issues, the major difference be-
ations and agree, on average, with those of the smaller tween the methods of WWW and Barkema-Mousseauis
model. Also, the larger system has a lower minimum the starting configuration used for the model construc-
frequency, allowing to calculate decay rates at the fre- tion: crystalline silicon (c-Si) is used in the former ap-
quencies previously unattainable. Studying paracrys- proach and a random close packed configuration in the
talline silicon, a material where small crystalline grains latter. Since neither of the starting models of this work
are embedded in a disordered matrix, allows us to test employthe SWpotential,bothmodelshadtoberelaxed
the hypothesis of van der Voort4 regarding the struc- with respect to SW prior to the decay times calculation.
tural origin of the anomalous long vibrational lifetimes This was done in an earlier work for the homogeneous
in a mixed amorphous-nanocrystallineSi system. In our model and in the present work for the paracrystalline
calculations we have used a 1000-atom (86 out of which model through the use of a molecular dynamics quench,
belong to a single crystalline grain) model created by where both temperature and virial pressure are set to
Nakhmansonet. al.15 tosimulatemedium-rangeorderin zerointhisstageofmodelpreparation. Thedensitiesfor
amorphous silicon. We should point out, that, although the two SW relaxedmodels are thus found to be slightly
providing a more realistic subject for the verification of (three to four percent) less than that of the density of
van der Voort’s hypothesis than “regular” models for a- c-Si, with the density of the paracrystalline model being
Si,thissimplemodelisneitheranexactstructuralmatch slightly (two percent) higher than that of the homoge-
3
4
3
S
O
2
D
V
1
0
0 10 20 30 40 50 60 70 80 90
PHONON ENERGY (meV)
FIG.2: Vibrationaldensityofstatesofthe1000-atommodels
ofa-Si(shadedarea) andp-Si(linewithout shading). Modes
with ω>72 meV (indicated by thevertical line) are locons.
notdecayexponentially(asatrulylocalizedmodewould)
FIG. 1: The 1000-atom model of p-Si. Atoms belonging to
the crystalline grain are shown in grey; those in white form outsideofthedefectregion. Inmodelsofa-Siresonances
theamorphous matrix. oftenhavefrequenciesbelowthelowestfrequency(which
is not zero due to the finite size) of the corresponding
models of crystalline silicon. The low frequency spec-
neous model. Changes in the atomic positions resulting tral region of propagons and resonances is a subject of
from the SW relaxation were found to be quite small great interest in relation to the so called boson peak,30
and the view of the SW relaxed model of p-Si shown in which was also investigatedin models of a-Si.31,32 Diffu-
Fig. 1 is very close to the unrelaxed structure. In gen- sons make up the majority of the spectrum. They have
eral it is known that the SW potential produces relaxed frequencies above the Ioffe-Regel limit (in principle, as
structuresthathavetwotothreepercentfive-foldcoordi- suggested in Ref. 33, there can be two such limits, one
natedatoms,evenifthestartingstructureswereperfectly for transverse and one for longitudinal acoustic modes),
four-foldcoordinated(whichis thecaseforthe paracrys- so that they cannot be characterized by a wave vector,
talline model). At variance with the optical properties but nevertheless are extended modes and contribute to
of the models, this deviation from the perfect four-fold thermalconductionwithinKubotheoryasshownby the
coordination does not noticeably alter their vibrational numericalcalculations ofRefs. 20,34. Finally, loconsare
properties and will not be discussed here in more detail. localized modes in the sense of strong (Anderson) local-
ization. Inmoststudiedmodelsloconsformonlyasmall,
high-frequencypartofthe spectrum. Inourmodelsofa-
III. HARMONIC VIBRATIONS Si many of the most localized locons reside in regions
where there are five-fold coordinated atoms.
In the harmonic approximation vibrational eigenfre- Experimentally the character of the atomic vibrations
quencies ω(i)and eigenvectorsei are computedby diag- in glasses has been studied by inelastic x-ray scattering
onalizing the corresponding dynaamical matrix (through- in various glassy systems.35,36,37,38,39 The general con-
out the paper symbols j, k, and l will represent vibra- census seems to be that sound waves cease to propagate
tional modes, while a, b, and c atoms). The results of early in the spectrum (for example, at 9 meV in den-
numericalcalculationsfromvariousgroups22,23,24,25,26,27 sified silica38 and at 5 meV in amorphous selenium35).
indicate that vibrational eigenstates in glasses belong to ThisconfirmstheexistenceoftheIoffe-Regelcrossoverin
one of four groups6,28,29: propagons, resonant modes, glassesandthepropagonpicture. Theissueofthenature
diffusons, and locons. Propagons are essentially sound of the vibrational modes above the crossover frequency
waves scattered by structural disorder. Resonances (where the models predict diffusons) has not been re-
(sometimes called quasilocalizedmodes) are modes tem- solved experimentally yet. The recent experimental and
porarilytrappedintopologicaldefects(inSimodelsthese theoretical progress is reviewed in Ref. 29.
are the groups of undercoordinated atoms). Spatially, a In Fig. 2 we plot the calculated vibrational density of
resonant mode has a large vibrational amplitude at the states (VDOS) for the models of a-Si and p-Si. Both
defect, but although being rather weak, the mode does curves look very similar, which is in agreement with
4
0.6
0.06
0.25
0.5
s
m
0.04
o
0.2 at0.4
-
c
1/p0.15 0.02 ht at 0.3
g
0.1 0 wei
0.2
4 5 6 7 8 e
d
o
0.05 m
0.1
0
0
0 10 20 30 40 50 60 70 80 90
0 10 20 30 40 50 60 70 80 90
PHONON ENERGY (meV)
PHONON ENERGY (meV)
FIG. 3: Inverse participation ratio 1/p of the vibrational
FIG. 4: Weight of the modes at the crystalline cluster as a
states in the 1000-atom models of a-Si (shaded area) and p-
function of mode frequency in the 1000-atom model of p-Si.
Si (line without shading). The modes with the frequencies
Plotted is thesquare of the atom displacement summed over
above 72 meV (thevertical line) can be considered localized.
the atoms forming the cluster. The horizontal line shows a
Quasilocalization occurs at low frequencies (resonant modes)
weight of 0.086 (8.6%) indicating an unbiased displacement
and around 30 meV which corresponds to the band edges.
pattern. The histogram is VDOSfor thecluster (see text).
The inset is a detailed view of thelow-frequency region.
the VDOS calculation of Ref. 15 made with a modified as would be expected from their idiosyncratic character,
version of the SW potential.40 The calculated spectrum canbelocalized(fullyorpartially)on,through,oroffthe
agrees rather well with the experimental one,20 except cluster. None of the modes is localized fully on the clus-
that the calculation overestimates the highest frequen- ter. There are 4 locons with the weight on the cluster of
cies by about 15%. This is a known artifact of the SW 30% or greater, the maximum weight being that of 59%
potential. for a mode with frequency ω = 73.05 meV and partici-
Localization properties of the modes can be judged pation ratio p=13. The second most localized mode on
from the participation ratio p(j), which indicates how the cluster has the frequency of 72.69 meV, the weight
many atoms “participate” in vibrational eigenmodes j. of 55% and p = 12. The third and fourth modes are
Inverseparticipationratio1/pfora-Siandp-Si,asafunc- more delocalized, having frequency (weight,p) ω =70.67
tion ofmode frequency is shownin Fig. 3. The majority (31%,160) and ω = 71.12 (30%,117), respectively. All
of vibrations in both models is delocalized, with the lo- fourmodeslieinthemobilityedgeregion. Inadditionto
calizationtransitiontaking place ataround72 meV (the these, there are modes with frequencies around 30 meV
mobility edge). The modes around 30 meV and some which have enhanced affinity for the cluster (see Fig. 4).
modes below 10 meV appear to be localized too. The Theweightofthese modesatthe crystallinecluster does
latter are resonant modes. The extended modes below not exceed 30%, but six modes have the weight between
about 15 meV are propagons, while all the rest are dif- 20 and 30%.
fusons (with possibly some longitudinal propagonsleft33 Harmonic vibrations in a-Si explain well many ob-
at small frequencies). The localization character in both served thermodynamic properties41 of the material, as
a-Siandp-Simodelsissimilar. Thepresenceofthecrys- well as kinetics such as heat flow.20,34 Anharmonicity
tallineclusterdoesnotleadtoadditionallocalizedmodes doesnotdirectlyaffectheatflowindielectricglasses,but
elsewhere in the spectrum. isveryimportantinrelaxingthe perturbedvibrationsto
In order to understand what fractionof eachmode re- maintain local equilibrium (temperature gradient, to be
sidesonthecrystallineclusterinthep-Simodel,wecom- specific). More directly, anharmonicity affects thermal
putetheweighteachmodehasonthecluster(thatis,we expansionandsoundattenuation. The 1000-atommodel
sum ej 2 for each j over all atoms a from the cluster). of a-Si was employed to demonstrate the importance of
The r|esau|lt is shown in Fig. 4, together with a histogram thermal vibrations in both of these phenomena.42,43 It
ofVDOSofthe cluster calculatedbysolvingthe dynam- was found that anharmonicity is rather weak in a-Si, al-
ical equations for the cluster atoms with the surround- thoughsomewhatstrongerthaninc-Si, primarilydue to
ing atoms held fixed. An unbiased mode has a weight strong anharmonicity of resonant modes. Indeed, res-
of 0.086 (8.6%), corresponding to the percentage of the onant modes show giant Gru¨neisen parameters in the
atomsmakingupthecluster. Forallthemodesbelowthe model, strongly enhancing the effects of anharmonicity,
mobility edge the weight fluctuates around 0.086, show- although still within the limits of perturbation theory
ing no special affinity for the cluster. Localized modes, based on cubic anharmonicity.
5
de 1 de 0.7
o o
m m
V V 0.6
e e
m m
91 91 0.5
4 4
2. 2.
1 1
meV) of 0.1 ww==00..122, 0m.4e,V ..., 2.4 meV) of 00..43 T1300=001 0KK K (times 20)
E ( E (
T T0.2
A A
R R
Y Y 0.1
A A
C C
DE0.01 E 0
D
0 100 200 300 400
0 0.4 0.8 1.2 1.6 2 2.4
TEMPERATURE (K)
w (meV)
FIG.5: Calculateddecayrateofthemode(apropagon)with
FIG. 6: Calculated decay rate of the mode with frequency
frequency ω = 12.49 meV in a-Si as a function of temper-
ω=12.49 meV in a-Siasa function of w at 10, 100, and 300
ature for different widths w of the rectangle function θw(ω) K.
representingthedeltafunctioninEq. 1. Thecurvesareforw
equal0.12(dashedline)and0.2,0.4,...,2.4meV(solidlines),
theorder beingnot mirrored in themagnitudeof thecurves.
In Eq. 1 the term with the temperature factor 1 +
The greatest decay rate is for w=0.2 meV, while thelowest
n(k)+n(l) corresponds to the decay j k +l, while
for w=1.4 meV. The curve representing w=1 meV chosen
→
in thecalculation is in themiddle of the bunch. the term with n(k) n(l) representsthe “difference”de-
−
cay j +k l. Energy conservation is ensured by the
→
delta functions. At low (down to zero) temperatures the
IV. VIBRATIONAL LIFETIMES
first term in Eq. 1 dominates, giving rise to a constant
2Γ, while both terms are generally equally important at
Usingcubicanharmonicityasthesmallperturbationto large temperatures, where Γ T. In crystals V(j,k,l)
∼
the harmonicHamiltonian,anharmonicdecayrate2Γ(j) vanishes unless the modes’ momentum is conserved in
of mode j can be obtained from the formula44,45 the decay process. In glasses, where lattice momentum
itself is not a valid concept (except for propagons and
¯h2π V(j,k,l)2
2Γ(j) = X| | resonances), all the modes k and l from the spectrum
4ω(j) ω(k)ω(l) contribute to V(j,k,l) for a given j.
kl
1 Crucial for determining decay rates perturbatively
( [1+n(k)+n(l)]δ[ω(j) ω(k) ω(l)]
from a finite-size model is the delta-function regular-
× 2 − −
ization. We approximate δ(ω) θ (ω), where θ (ω)
+[n(k) n(l)]δ[ω(j)+ω(k) ω(l)]). (1) w w
≈
− − is a rectangle of width w and height 1/w centered at
Here ω(j) is the frequency of mode j, n(j) is the mode ω = 0. In our calculations with 1000 atoms we choose
occupation number given by n(j) = exp[h¯ω(j)/kBT] w =1 meV, which fits about 40 modes in the rectangle.
1 −1 with T denoting temperature, a{nd V(j,k,l) is th−e The choice ofw needs to be a compromise betweengood
}
matrix element of the cubic anharmonicity of the inter- statistics and computer power. The statistics is deter-
atomic potential V in the harmonic representation: minedbyboththenumberand“similarity”ofthemodes
in a rectangle. If w is too large, the rectangle function
∂3V ej ek el
V(j,k,l)=XX aα bβ cγ . (2) will sample modes with distinct characteristics, not rep-
∂uaα∂ubβ∂ucγ √ma√mb√mc resenting faithfully the modes of the chosen frequency.
abc αβγ
This problem is likely to be absent for diffusons, which
Greek symbols α, β, and γ stand for the cartesian co- do not differ much on small spectral scales due to the
ordinates of both the atomic displacements u from the absenceofdegeneracy(cf. Ref. 42), butmaybe relevant
equilibrium positions, and of the normalized vibrational for propagons (which are mixed with resonant modes)
eigenvectors e. The atomic masses are denoted as m. and locons (which are idiosyncratic). Fortunately, the
Anharmonic vibrational lifetimes are the inverse of the averaging,firstwithinthe rectangleandsecond,overthe
rates: whole spectrum (see Eq. 1) makes the decay rates quite
insensitive to the choice of w, for a reasonable interval
τ(j)=1/2Γ(j). (3)
of values. In the earlier calculation6 w was chosen to be
Inthispaperwepresentdecayratesinthe units ofmeV. 0.4 meV for a 216-atom model, fitting about 4 modes in
For conversion into lifetimes, a decay rate of 1 meV is the rectangle. As we will see from the comparisonof the
equivalent to a lifetime of about 0.7 ps. two calculations in the following section, this choice was
6
1 3
T=10 K 216 T=300 K 216
0.8 1000 1000
V) V)
me me2
E (0.6 E (
T T
A A
R R
Y 0.4 Y
A A
EC C1
D E
0.2 D
0
0
0 10 20 30 40 50 60 70 80 90
0 10 20 30 40 50 60 70 80 90
PHONON ENERGY (meV)
PHONON ENERGY (meV)
FIG.7: Calculateddecayratesversusfrequencyforthe1000-
FIG.8: Calculateddecayratesversusfrequencyforthe1000-
atom model at 10 K (thick line). For comparison the rates
atom model at 300 K. For comparison the rates of the 216-
of the 216-atom model from Ref. 6 (not smoothed) are also
atom model from Ref. 6 (not smoothed) are also shown
shown (shaded area).
(shaded area).
already good enough, although it may have contributed
decayratesasafunctionoffrequencyat10K(andatlow
somewhat to the statistical noise, especially at low tem-
temperatures,generally)followthejointdensityofstates
peraturesandsmallfrequencies. Toillustratetheeffectw
[ δ(ω(j) ω(k) ω(l))] which counts, for a chosen
has on the decay rates, we show in Fig. 5 the calculated Pkl − −
mode j, the number of combination decay possibilities
ratesofthe mode with ω =12.49meVina-Si,as afunc-
j k+l with the constraintof energy conservation. At
tion of temperature, for selected w, ranging from 0.12 →
larger temperatures one must add the number of differ-
meV (corresponding to about 4 modes per rectangle) to
encedecaychannelsj k ltoreproduce,qualitatively,
2.4 meV (80 modes/rectangle). Except for the smallest → −
the calculated 2Γ(ω). More detailed physics of the an-
w, the results are grouped together with the dispersion
harmonic decay in glassesand especially the statistics of
of less than 10% above 100 K. The greatestdispersion is
the decay matrix elements can be found in Ref. 6.
at the lowest temperatures, where it reaches 25%. (The
There are several features which make the calculated
low temperature properties of the model do not describe
decay rates for the 216-atom and 1000-atom models
well the real a-Si structure, because of the existence of
somewhat different. The first is the overall reduction
the minimum frequency of 4 meV in the model). Figure
in noise for the 1000-atom model (the data are not
6 shows the decay rate for the same mode as a function
smoothed as was done in Ref. 6). The reason is both
of w, for selected temperatures. The rates become rea-
the greater model size (spectral averaging) and greater
sonably insensitive to w above 0.4 meV. The dispersion
w (rectangle averaging). Note that the observed noise
duetothe sensitivityonw isafactorcontributingtothe
in the spectrally resolved 2Γ is about 10% or less, con-
uncertainty of the calculated values.
sistent with a dispersion of the decay rates with w, dis-
cussed in the previous section. Second, the calculated
ratesforthe1000-atommodelaresomewhatsmallerthan
V. AMORPHOUS SILICON those of the 216-atom model, that is, the latter model
appears to be slightly more anharmonic. This is at vari-
Wenowpresentthecalculateddecayratesforthe1000- ance with the calculation of thermal expansion42 where
atommodelofa-Si. Duetothecomputationalpowerlim- the 216-atom model seems less anharmonic. The latter
itations we have sampled the spectrum uniformly with differenceprobablycanbe explainedbythe anomalously
about 200 modes for which we computed 2Γ. The cal- largenegativemodeGruneisenparametersofthelowfre-
culated decay rates are presented as a function of the quency resonance modes of the 1000 atom model, as the
mode frequency for two different temperatures: 10 K in thermodynamic Gruneisen parameter depends on an av-
Fig. 7 and300K in Fig. 8. For comparisonthe previous erage mode Gruneisen parameter at high temperatures.
calculations6 on a 216-atom model of a-Si are included. We note that the structural models differ in other ways:
Overall, the decay rates for the two models agree. The the smaller model is more topologically constrained,42
ratesareontheorderofmeV(picosecondlifetimes). Per- has smaller energy/atom, and has higher density than
turbation theory is thus valid for all the sampled modes the 1000-atom model. Third, the calculated rates of
with the exception of few in the lowest part of the spec- the 1000-atommodel extend to a lower frequency region
trum at 300 K (see below). As was shown in Ref. 6 the as the minimum frequency of the model is smaller than
7
10 1
0.4
a-Si
0.8 p-Si
V) 1 V)
me me 0.2
TE ( E (0.6 T=10 K
A 0.1 T
R 5.40 meV A
CAY 2505..7254 AY R0.4 025 30 35
E 78.78 C
D 0.01 E
D
0.2
0.001
0
0 200 400 600 800 1000
TEMPERATURE (K) 0 10 20 30 40 50 60 70 80 90
PHONON ENERGY (meV)
FIG. 9: Calculated decay ratesof selected modes in a-Siver-
sus temperature. The lines are labeled according to modes’ FIG.10: Calculateddecayratesofthe1000-atommodelofp-
frequenciesinmeV.Thelowestfrequencymodeisapropagon, Siat10K.Forcomparison,theratesforthe1000-atommodel
the following two are diffusons (acoustic-like and optic–like) of a-Si are also shown (shaded area). The empty circles are
and the highest frequency modeis a locon. for three modes in the mobility region with more than 30%
weight on the crystalline cluster, while the inset shows the
decay rates (filled circles) of modes around 30 meV, which
havelarge affinity (weight up to 30 %) for the cluster.
thatofthe216-atommodel. Finally,somelow-frequency
modes (resonances) at 300 K exhibit giant decay rates,
comparable to the modes’ frequencies. These rates are
in fact invalid, since they are not consistentwith pertur- circles. The decay rates of these modes have the same
bation theory. However, they indicate what may be ex- magnitude as those of the other locons. Finally, in the
pected from a full anharmonic calculation (for example insetsofthetwofiguresweshowthedecayratesofmodes
by molecular dynamics). This important physics issue with frequencies around30 meV, the regionof especially
will be discussed elsewhere. high affinity for the cluster (see Fig. 4). Decay rates of
In Fig. 9 we plot the temperature dependence of the more than 80 modes in that spectral region are plotted.
decayratesofselectedmodes. We showthe temperature Although many of the modes have large weight (some of
dependenceforapropagon,anacoustic-likeandanoptic- them up to 30%)on the cluster, most are unbiased. The
like diffuson, and a locon. The low-frequency propagon fact that 2Γ of all of these modes are similar in magni-
has a divergentlifetime (decay rate vanishes) as temper- tude at different temperatures implies no special decay
ature decreases to zero, since there are no two modes behavior for the modes of strong affinity for the cluster.
intowhichitcoulddecay,duetotheenergyconservation
Figure 12 shows the temperature dependence of three
constraint and the existence of the minimum-frequency
modes with more than 30% weight on the crystalline
mode. All the other modes have constant decay rates
cluster. The modes have frequencies (weight,p) 70,67
at small temperatures. The constant goes smoothly to
meV (31%,160), 72.68 meV (55%,12), and 73.05 meV
a linear function at large temperatures, which is due to
(59%,13). In addition, the figure plots the decay rate of
the fact that the population density of thermal phonons
a“normal”loconwithω =77.76meV(0.02%,8),residing
increases linearly with temperature.
outside of the cluster. No unusual features are observed
inthetemperaturedependenceofthedecayratesofthese
modes. Finally, in Fig. 13 we plot the anharmonic ma-
VI. PARACRYSTALLINE SILICON trixelementsV(j,k,l)ofthecombinationdecayj k+l
→
forthemaximallylocalizedmodeonthecluster,withfre-
The resultsforthe 1000-atommodelofp-Siareshown quency 73.05 meV to visualize the decay channels. The
in Figs. 10 and 11, which plot 2Γ as a function of mode figure shows that the dominant channel is a decay into
frequency. Forcomparisonwealsopresentthedatafora- twodiffusons. Decayintoapropagonandadiffuson(the
Sidiscussedintheprevioussection. Theresultsarequan- points in Fig. 13 below 15 meV and above 58 meV) is
titatively similar for both models. There are no anoma- somewhat less important. The corresponding matrix el-
lousdecayratesappearinginthe spectrumofp-Siwhich ements are much smaller. This may be related to the
would be due to the crystalline cluster. In addition to fact that propagon’s weight on the cluster is systemati-
the sampling modes, we have computed the decay rates callylowerthan8.5%(see Fig. 4). The diffusons’weight
specifically for three modes in the mobility edge region on the cluster is much more more scattered, with a sig-
with the weight at the crystalline cluster greater than nificant number of diffusons having the weight of 8.5%
30%. They are presented in Figs. 10 and 11 by empty and more. Decay into another locon and a propagon is
8
5 2
2)
2 3/
a-Si -Si
m
)4 p-Si -3A 1
V J
me 1 T=300 K 20
ATE (3 -NT (10 0
Y R2 0 ME
A 25 30 35 E
C L
E E -1
D X
1 RI
T
A
M
-2
0
0 10 20 30 40 50 60 70 80 90
0 10 20 30 40 50 60 70 80 90 PHONON ENERGY (meV)
PHONON ENERGY (meV)
FIG.13: Calculated matrix element V(j,k,l) ofp-Sifor mode
FIG. 11: Calculated decay rates of the 1000-atom model of
jmaximallylocalized(59%)onthecrystallinecluster[ω(j)=
p-Si at 300 K. For comparison, the rates for the 1000-atom
73.05 meV] as a function of ω(k). Shown are data for k and
model of a-Si are also shown (shaded area). The empty cir-
l obeying energy conservation: |ω(j)−ω(k)−ω(l)| < w/2,
cles are for three modes with more than 30% weight on the
where w = 0.2 meV (taken to be smaller than w = 1 meV
crystalline cluster, and the inset plots the decay rates (filled
used in thecalculation, toget a manageable graphics size).
circles)ofthemodesaround30meVwithlargeweightonthe
cluster.
3 are on the picosecond time scales, generally increasing
with increasing frequency. The decay rates of locons are
70.67 meV (31%) idiosyncratic,but are by no means inhibited. Calculated
72.69 (55%)
) decay rates of p-Si are similar to those of a-Si, showing
V
me2 73.05 (59%) little sensitivity to structural properties. These findings
E ( 77.76 (0.02%) disagree with the interpretation of recent experiments
T which find decay rates on the order of nanoseconds and
A
R somewhat greater sensitivity to structural properties.
Y
CA1 Theexplanationthatweoffertoaccountforthese dis-
E crepancies is that the calculation and experiment refer
D
to two different things. First, as we have pointed out
earlier, simple (and at such scale usually over-relaxed)
models like a continuous random network type WWW
0
model or a similar model containing a crystalline grain
0 200 400 600 800 1000
used in this study can not faithfully reproduce a broad
TEMPERATURE (K)
range of various topological features — some or combi-
nations of which may be responsible for increased decay
FIG. 12: Calculated decay rates for selected locons in p-Si.
timesobservedintheexperiment—presentinarealma-
The curves are labeled according to frequency in meV. The
numbersin thebracketsshowthemodes’weight onthecrys- terial. Second, in our calculations only “perturbative”
talline cluster. decayrates,whereasmall(infinitesimal)populationofa
single mode goes out of equilibrium are computed. The
experiments measure the decay of vibrational states ex-
forbidden by energy conservation. cited overa large portionof the spectrum. Furthermore,
the laserexcitationproducesphononpopulationstoo far
off the equilibrium to be called small perturbations. In
VII. CONCLUSIONS Ref.3,forexample,theexcitedphononpopulationnlies
between 0.03 and 0.3, which, for a mode with frequency,
We have calculated anharmonic decay rates of 1000- say, 50 meV corresponds to an effective temperature of
atom models of a-Si and p-Si using perturbation theory 160 and 400 K, respectively. This is huge compared to 2
with cubic anharmonicity in the interatomic potential. K at which the samples are held. A numerical investiga-
The results for a-Si are in agreement with the previous tionofBickham8indeedshowsthatastrongperturbation
perturbative calculations on a smaller model, as well as of the vibrationalspectrum of a-Si can relax on a 100 ps
with a molecular dynamics calculation. The results reit- time scale, compared to 10 ps for a weak perturbation.
erate the previous findings that the vibrational lifetimes Inadditiontopurevibrationalrelaxation,itisalsolikely,
9
assuggestedbyBickhamandFeldman,7thatcorrespond- different models (of sizes from 216 to 4096 atoms) and
ingly large local deviations in the atomic displacements different techniques (perturbation theory and molecular
cause local structural rearrangements which may relax dynamics) it is rather unlikely that the interpretation of
tolocalmetastableminimawhileemittingphonons. Itis the experimental findings in terms of the (pure) vibra-
possible that these rearrangementsoccur on nanosecond tional lifetimes (that is, without considering structural
timescales. relaxationandpossiblyphotoexcited-electronrecombina-
We conclude that current simple models for a-Si in tion processes) will be validated by investigating larger
combinationwithpresentedabovemethodsofanalysisdo models.
not provide an answer to the question why vibrations in
real material under the specific experimental conditions
appear to decay oversuch long time scales. A better an-
Acknowledgments
swer may be obtained by studying larger and more real-
isticmodelswhencomputationalpowerrequiredtomake
such calculations tractable becomes available. However, We thank Phil Allen for useful discussionsandM. van
considering the consistent results obtained thus far with der Voort for suggesting this calculation.
1 A.Scholten,A.V.Akimov,and J. I.Dijkhuis,Phys.Rev. 24 H.R.SchoberandC.Oligschleger,Phys.Rev.B53,11469
B 47, 13910 (1993). (1996).
2 A. Scholten and J. I. Dijkhuis, Phys. Rev. B 53, 3837 25 S. N. Taraskin and S. R. Elliott, Phys. Rev. B 61, 12031
(1996). (2000).
3 M.vanderVoort,A.V.Akimov,andJ.I.Dijkhuis,Phys. 26 S. N. Taraskin and S. R. Elliott, Phys. Rev. B 61, 12017
Rev.B 62, 8072 (2000). (2000).
4 M. van der Voort, O. L. Muskens, A. V. Akimov, A. B. 27 M. Marinov and N.Zotov, Phys. Rev.B 55, 2938 (1999).
Pevtsov, and J. I. Dijkhuis, Phys. Rev. B 64, 045203 28 P.B.Allen,J.L.Feldman,J.Fabian,andF.Wooten,Phil.
(2001). Mag. B 79, 1715 (1999).
5 J. Menendez and M. Cardona, Phys. Rev. B 29, 2051 29 E. Courtens, M. Foret, B. Hehlen, and R. Vacher, Solid
(1984). State Commun. 117, 187 (2001).
6 J.FabianandP.B.Allen,Phys.Rev.Lett.77,3839(1996). 30 T. Nakayama, Rep.Prog. Phys. 65, 1195 (2002).
7 S.R.BickhamandJ.L.Feldman,Phys.Rev.B57,12234 31 J.L.Feldman,P.B.Allen,andS.Bickham,Phys.Rev.B
(1999). 59, 3551 (1999).
8 S. R.Bickham, Phys.Rev.B 59, 4894 (1999). 32 S. M. Nakhmanson and D. A. Drabold, Phys. Rev. B 61,
9 S.Alexander,C. Laermans, R.Orbach,and H.M. Rosen- 5376 (2000).
berg, Phys.Rev.B 28, 4615 (1983). 33 J. L. Feldman, J. Non-Cryst.Solids 307-310, 128 (2002).
10 A. Jagannathan, R. Orbach, and O. Entin-Wohlman, 34 P. B. Allen and J. L. Feldman, Phys. Rev. B 48, 12581
Phys. Rev.B 39, 13465 (1989). (1993).
11 R. Orbach and A. Jagannathan, J. Phys. Chem. 98, 7411 35 M.Foret,B.Hehlen,G.Taillades,E.Courtens,R.Vacher,
(1994). H. Casalta, and B. Dorner, Phys. Rev. Lett. 81, 2100
12 J. Fabian, Phys.Rev.B 55, 3328 (1997). (1998).
13 D. Leitner, Phys. Rev.B 64, 094201 (2001). 36 G. Ruocco, F. Sette, R. D. Leonardo, D. Fioretto,
14 D. Leitner, Phys. Rev.Lett. 87, 188102 (2001). M. Krisch, M. Lorenzen, C. Masciovecchio, G. Monaco,
15 S. M. Nakhmanson, P. M. Voyles, N. Mousseau, G. T. F. Pignon, and T. Scopigno, Phys. Rev. Lett. 83, 5583
Barkema, and D. A. Drabold, Phys. Rev. B 63, 235207 (1999).
(2001). 37 O. Pilla, A. Cunsolo, A. Fontana, C. Masciovecchio,
16 P. M. Voyles, J. E. Gerbi, M. M. J. Treacy, J. M. Gibson, G. Monaco, M. Montagna, G. Ruocco, T. Scopigno, and
and J. R.Abelson, Phys.Rev.Lett. 86, 5514 (2001). F. Sette, Phys.Rev.Lett. 85, 2136 (2000).
17 F.Wooten,K.Winer,andD.Weaire,Phys.Rev.Lett.54, 38 E. Rat, M. Foret, E. Courtens, R. Vacher, and M. Arai,
1392 (1985). Phys. Rev.Lett. 83, 1355 (1999).
18 G. T. Barkema and N. Mousseau, Phys. Rev. Lett. 77, 39 M. Foret, R. Vacher, E. Courtens, and G. Monaco, Phys.
4358 (1996). Rev.B 66, 024204 (2002).
19 N. Mousseau, G. T. Barkema, and S. M. Nakhmanson, 40 R. L. C. Vink, G. T. Barkema, W. F. van der Weg, and
Phil. Mag. B 82, 171 (2002). N. Mousseau, J. Non-Cryst.Solids 282, 248 (2001).
20 J. L. Feldman, M. D. Kluge, P. B. Allen, and F. Wooten, 41 J. L.Feldman, C. S.Hellberg, G. Viliani, W.Garber, and
Phys. Rev.B 48, 12589 (1993). F. M. Tangerman, Phil. Mag. B 82, 133 (2002).
21 F. H. Stillinger and T. A. Weber, Phys. Rev. B 31, 5262 42 J.FabianandP.B.Allen,Phys.Rev.Lett.79,1885(1997).
(1985). 43 J.FabianandP.B.Allen,Phys.Rev.Lett.82,1478(1999).
22 R. Biswas, A. M. Bouchard, W. A. Kamikatahara, G. S. 44 A. A. Maradudin and A. E. Fein, Phys. Rev. 128, 2589
Grest, and C. M. Soukolis, Phys. Rev. Letters 60, 2280 (1961).
(1988). 45 R.A. Cowley, Adv.Phys. 12, 421 (1963).
23 H. R. Schober and B. B. Laird, Phys. Rev. B 44, 6746
(1991).
