Ab initio calculations of third-order elastic constants and related properties for selected semiconductors
Michal Lopuszy? ski n
Interdisciplinary Centre for Mathematical and Computational Modelling, University of Warsaw, Pawi?skiego 5A, 02-106 Warsaw, Poland n
arXiv:cond-mat/0701410v3 [cond-mat.mtrl-sci] 4 Sep 2007
Jacek A. Majewski
Institute of Theoretical Physics, Faculty of Physics, University of Warsaw, Ho˙ a 69, 00-681 Warsaw, Poland z (Dated: February 6, 2008) We present theoretical studies for the third-order elastic constants Cijk in zinc-blende nitrides AlN, GaN, and InN. Our predictions for these compounds are based on detailed ab initio calculations of strain-energy and strain-stress relations in the framework of the density functional theory. To judge the computational accuracy, we compare the ab initio calculated results for Cijk with experimental data available for Si and GaAs. We also underline the relation of the third-order elastic constants to other quantities characterizing anharmonic behaviour of materials, such as pressure derivatives of the second-order elastic constants c′ and the mode Gr¨neisen constants for long-wavelength u ij acoustic modes γ(q, j). Paper accepted to Phys. Rev. B.
PACS numbers: 62.20.Dc,43.25.+y,62.50.+p
I.
INTRODUCTION
Third-order elastic constants Cijk are important quantities characterizing nonlinear elastic properties of materials and the interest in them dates back to the beginning of modern solid state physics.1,2,3,4,5 Third- and higherorder elastic constants are useful not only in describing mechanical phenomena when large stresses and strains are involved (e.g., in heterostructures of optoelectronic devices), but they can also serve as a basis for discussion of other anharmonic properties. The applications include phenomena such as thermal expansion, temperature dependence of elastic properties, phonon-phonon interactions etc.6 As far as theoretical studies are concerned, at the beginning the third-order elastic constants were calculated in the framework of the valence force Keating model.7 Later on, many other more sophisticated microscopic theories were employed to describe and predict nonlinear elastic properties of crystals on the basis of their atomic composition.6 Nowadays, precise ab initio calculations seem to be the most promising approach to handle this task. Such applications of density functional theory (DFT) on the local density approximation level (LDA) were already reported.8,9 Recently, one observes increased interest in nonlinear e?ects in elastic10,11,12 and piezoelectric properties.13,14 This is strongly connected to the fact that research focuses nowadays on the semiconductor nanostructures. In such systems these nonlinear e?ects are not only more pronounced than in bulk materials, but very often their reliable quantitative description is a prerequisite for correct theoretical explanation of the experimental data.12,14,15,16,17,18 In this paper, we perform ab initio
calculations of the unknown third-order elastic constants in cubic nitrides. The nitrides are technologically important group of materials for which the nonlinear effects are particularly signi?cant.10,12,19,20 Therefore, the knowledge of the third-order elastic moduli will de?nitely improve the modeling of nitride based nanostructures. In this work we also brie?y discuss the applications of Cijk to determination of other anharmonic properties, namely, pressure derivatives of second-order elastic moduli c′ and mode Gr¨ neisen constants γ(q, j). Since the u ij third-order e?ects are rather subtle, their computational determination can also serve as a precise test of accuracy for modern ab initio codes based on DFT approach. The paper is organized as follows. In Sec. II we give a general overview of the nonlinear elasticity theory. Sec. III contains a description of employed methodology. Also results for third-order elastic constants obtained from ab initio calculations are presented there. Our ?ndings for Si and GaAs are compared with previous numerical calculations and measurements, later on theoretical predictions for zinc-blende nitrides AlN, GaN, and InN are given. Secs. IV and V deal with the determination of quantities related to third-order elastic constants, namely, the pressure dependent elastic constants and mode Gr¨ neisen u constants, respectively. Finally, we conclude the paper in Sec. VI.
II.
OVERVIEW OF NONLINEAR ELASTICITY THEORY
Here we will recall some basic facts from nonlinear theory of elasticity.1,2,3,4,5,6 Let us consider point a which, after applying strain to a crystal, moves to the position
2 x. After introducing the Jacobian matrix J Jij = ?xi ?aj (1)
A. Methodology and computational details III. DETERMINATION OF THIRD-ORDER ELASTIC CONSTANTS
we may de?ne the Lagrangian strain 1 T (J J ? 1), (2) 2 which is a convenient measure of deformation for an elastic body. The energy per unit mass E(η) corresponding to the applied strain may be developed in power series with respect to η. This leads to the expression η= ρ0 E(η) = 1 1 cij ηi ηj + 2! i,j=1,6 3! Cijk ηi ηj ηk + . . . ,
i,j,k=1,6
In this work, we have determined third-order elastic constants for Si, GaAs, and zinc-blende nitrides (AlN, GaN, and InN) on the basis of quantum DFT calculations for deformed crystals. The results were obtained in two ways - employing strain-energy formula [Eq. (5)] and from strain-stress relation [Eqs. (6) and (7)]. The detailed procedure was as follows. We considered six sets of deformations parametrized by η ηA ηB ηC ηD ηE ηF = = = = = = (η, 0, 0, 0, 0, 0), (η, η, 0, 0, 0, 0), (η, 0, 0, η, 0, 0), (η, 0, 0, 0, η, 0), (η, η, η, 0, 0, 0), (0, 0, 0, η, η, η).
(3) where we applied Voigt convention (η11 → η1 , η22 → η2 , η33 → η3 , η23 → η4 /2, η13 → η5 /2, η12 → η6 /2) and introduced the density of unstrained crystal ρ0 . The cij and Cijk denote here second- and third-order elastic constants respectively.41 If we introduce J = (1 + ?) and assume that ? is symmetric (rotation free) linear strain tensor, the de?nition of η [Eq. (2)] yields
(8)
In every case, η was varied between ?0.08 and 0.08 with step 0.008. For every deformed con?guration, the po1 sitions of atoms were optimized and both energy and η = ? + ?2 . (4) 2 stress tensors were calculated on the basis of quantum Substituting the above result to the expansion in Eq. (3) DFT formalism. In this way, for each type of distortion, and leaving only terms up to second order with respect dependencies of energy E(η) and stress tensor t(η) on to components of ? recover the in?nitesimal theory of strain parameter η were obtained. The numerical results elasticity. have been in turn compared with the expressions from Naturally, the general expression for energy of strained the nonlinear theory of elasticity, which are summarized crystal, as given by Eq. (3), can be simpli?ed by emin Table I. This allows to extract the values of the secondploying symmetry considerations. For cubic crystals, this and third-order elastic constants, by performing suitable procedure yields the following formula: polynomial ?ts. ` 2 ? 1 ` 2 ? 1 2 2 2 2 The DFT calculations have been performed using the ρ0 E(η) = c11 η1 + η2 + η3 + c44 η4 + η5 + η6 + 2 2 ab initio total energy code VASP developed at the Insti+c12 (η1 η2 + η3 η2 + η1 η3 ) + tut f¨ r Materialphysik of Universit¨t Wien.21,22,23 The u a ` 3 ? 24 1 3 3 (5) projector augmented wave (PAW) approach has been C111 η1 + η2 + η3 + 6 used in its variant available in the VASP package.25 For ? ` 1 2 2 2 2 2 2 C112 η2 η1 + η3 η1 + η2 η1 + η3 η1 + η2 η3 + η2 η3 + the exchange-correlation functional generalized gradient 2 approximation (GGA) according to Perdew, Burke and ` ? 1 2 2 2 Ernzerhof (PBE)26,27 has been applied. For Ga and In, C123 η1 η2 η3 + C144 η1 η4 + η2 η5 + η3 η6 + 2 semicore 3d and 4d electrons have been explicitly in` ? 1 2 2 2 2 2 2 C155 η2 η4 + η3 η4 + η1 η5 + η3 η5 + η1 η6 + η2 η6 + cluded in the calculations. 2 Since the determination of subtle third-order e?ects C456 η4 η5 η6 + . . . requires high precision, we have performed careful convergence tests for parameters governing the accuracy of Another fundamental quantity in the theory of ?nite computations. On the basis of our tests we have chosen Si GaAs deformations is Lagrangian stress the following energy cuto?s Ecuto? = 600 eV, Ecuto? = AlN GaN InN 700 eV, and Ecuto? = Ecuto? = Ecuto? = 800 eV. For the ?E , (6) tij = ρ0 Brillouin zone integrals we have followed the Monkhorst?ηij Pack scheme28 , in Si and GaAs we have used 13 × 13 × 13 which can be expressed in terms of linear stress tensor σ mesh, whereas for AlN, GaN, and InN we have applied using the following formula 11×11×11 sampling. One example of performed tests for ?1 GaN is presented in Fig. 1. It illustrates the dependence t = det(J)J ?1 σ J T . (7) of two sample elastic moduli C111 and C144 on the energy cuto? and density of Monkhorst-Pack k-point mesh. For Again, Voigt convention (t11 → t1 , t22 → t2 , t33 → t3 , GaN the chosen parameters (Ecuto? = 800 eV and 11 × 11 × 11 t23 → t4 , t13 → t5 , t12 → t6 ) is used here.
3
TABLE I: Dependencies of energy and stress on deformation parameter η for considered types of deformation ηA , . . . , ηF , which have been used to determine second- and third-order elastic constants. Energy: . 1 ρ0 E(ηA ) = ` C111 η 3 + 1 c11 η 2 = fA (η) 6 2 ? . 3 1 ρ0 E(ηB ) = ` 3 C111 + C112 η + (c` + c12 )η 2 = fB (η) ? 3 11 ? . 1 1 1 1 ρ0 E(ηC ) = ` 6 C111 + 2 C144 ? η + ` 2 c11 + 2 c44 ? η 2 = fC (η) 3 2 . 1 1 1 1 ρ0 E(ηD ) = ` 6 C111 + 2 C155 η + ?2 c11 + 2 c44 η = ? D (η) f ` . 1 ρ0 E(ηE ) = 2 C111 + 3C112 + C123 η 3 + 3 c11 + 3c12 η 2 = fE (η) 2 2 . 3 ρ0 E(ηF ) = C456 η + 3 c44 η = fF (η) 2 . t1 (ηA ) = 1 C111 η 2 + c11 η = gA1 (η) 2 . t2 (ηA ) = 1 C112 η 2 + c12 η = gA2 (η) 2 . 2 t3 (ηB ) = (C123 + C112 ) η + 2c12 η = gB (η) . 2 t4 (ηC ) = C144 η + c44 η = gC (η) . t5 (ηD ) = C155 η 2 + c44 η = gD (η) . 2 t4 (ηF ) = C456 η + c44 η = gF (η)
Stress:
C111 [GPa]
?1180
?45
B.
Results and discussion
?1200
?50
Results are presented in Tables II and III. Table II contains our ?ndings for benchmark materials Si and GaAs, accompanied by available experimental data and previous theoretical ?ndings within LDA-DFT theory. Table III gives our prediction for the unknown values of Cijk for cubic nitrides. For completeness, we also provide there our prediction for second-order elastic moduli and compare them with previous calculations.12 For cij values, sometimes it was possible to determine one constant from a few ?ts [e.g., c44 from coe?cients in fC (η), fD (η), and fF (η)], obtaining slightly di?erent results [e.g., for GaN, c44 = 145, 151, 147 GPa from fC (η), fD (η), and fF (η) respectively]. In such cases the average of all obtained values was given in the tables. The sample plots of both energy and stress dependencies for GaN together with ?tted polynomials are depicted in Figs. 2 and 3. When analyzing the above results, one has to bear in mind that both measurements and calculations of the third-order elastic constants are di?cult. The reported experimental results for Cijk are determined with signi?cant uncertainties and quite often exhibit discrepancies between ?ndings of di?erent groups (see, e.g., GaAs in Table II). On the other hand, calculations of subtle third-order e?ects require reaching the limits of accuracy of modern quantum codes. When comparing experimental values with DFT results, it is also worth noticing that ab initio calculations are strictly valid for perfect crystalline structure and in
?1220 300 ?1204
400
500
600 700 800 Energy cutoff [eV]
900
1000
?55 1100
(b)
?1206
C111 C144
?44
C111 [GPa]
?1208
?1210
?45
?1212 7 9 11 k?points grid size 13
FIG. 1: Sample convergence tests for the third-order elastic constants in zinc-blende GaN. Panel (a) illustrates the dependence of C111 and C144 on energy cuto? (Monkhorst-Pack sampling 11 × 11 × 11 was applied for all points). Panel (b) shows the analogous dependence on the density of k-points mesh (energy cuto? 800 eV was used for all points). Note di?erent scales for C111 and C144 . Later on, all calculations for GaN have been performed with the energy cuto? 800 eV and the Monkhorst-Pack sampling 11 × 11 × 11.
C144 [GPa]
C144 [GPa]
k-point mesh) the di?erence between successive values of examined constants in our test is lower than 1 GPa. This di?erence is smaller than e.g. discrepancies observed between results obtained from strain-energy and strainstress approach which, in the opinion of the authors, indicates that the convergence with respect to parameters responsible for numerical accuracy is very reasonable.
?1140
(a)
?1160
C111 C144
?35
?40
4
0.030 fA(η) fB(η) fC(η) gA1(η) gA2(η) gB(η)
20
0.025 Lagrangian Stress [GPa] ?0.08 0.060 ?0.06 ?0.04 ?0.02 0 0.02 Lagrangian Strain η 0.04 0.06 0.08 10
0.020 Energy [eV/?3]
0
0.015
0.010
?10
0.005
?20
0.000
?30 ?0.08 ?0.06 ?0.04 ?0.02 0 0.02 Lagrangian Strain η 0.04 0.06 0.08
fD(η) fE(η) fF(η)
20
gC(η) gD(η) gF(η)
0.050 Lagrangian Stress [GPa] ?0.08 ?0.06 ?0.04 ?0.02 0 0.02 0.04 0.06 0.08 10
0.040 Energy [eV/?3]
0
0.030
0.020
?10
0.010
?20
0.000 Lagrangian Strain η
?30 ?0.08 ?0.06 ?0.04 ?0.02 0 0.02 0.04 0.06 0.08 Lagrangian Strain η
FIG. 2: Plots of strain-energy relations for GaN. Squares denote results of DFT computations, solid line represents a polynomial ?t. See main text for details and Table I for de?nitions of fA (η),fB (η),. . . ,fF (η).
FIG. 3: Plots of strain-stress relations for GaN. Circles denote results of DFT computations, solid line represents a polynomial ?t. See main text for details and Table I for de?nitions of gA1 (η),gA2 (η),. . . , gF (η)). Curves for gC (η) and gF (η) coincide because of very similar values of C144 and C456 in GaN.
the limit of 0K temperature. The experiments, however, are often performed in conditions far from this idealized case. Particularly, the importance of temperature factor can be veri?ed when comparing the results of measurements for Cijk of Si in temperatures T = 298K and T = 4K given in Table II (see Ref. 29 for detailed experimental study). One can observe that for this semiconductor the values of constants C144 and C123 even change their sign, when the material is cooled down. As far as calculations of third-order elastic moduli are concerned, they also pose a di?cult test to ab initio methods. The determination of Cijk is sensitive to errors in energy and stress tensor and requires extremely good convergence of parameters governing the accuracy of computations, which we believe has been reached in our calculations (see Fig. 1). The usage of PAW formalism chosen to solve Kohn-Sham equations seems also not to in?uence the results signi?cantly, since it has been demonstrated that properly performed calculations of the static and dynamical properties for broad range of solids within the PAW, pseudopotential, and LAPW schemes lead to essentially identical results.30 In our opinion, the main problem lies in the approximations to the exchangecorrelation functionals employed in various calculations.
In the present calculations we use GGA-PBE exchangecorrelation functional that is commonly believed to be one of the best in the market. However, even for the second-order elastic constants for GaAs (see Table II), one observes signi?cant di?erences between the calculated and measured values. One possible origin of these discrepancies might be the commonly known tendency of calculations based on GGA functional to underestimate binding strength, and therefore to overestimate lattice constant. Indeed, our calculations predict the equilibrium lattice constant of GaAs to be 5.75 ?, considerably A larger than the experimental value of 5.65 ?.31 This is opA posite to the local density approximation (LDA), which overestimates the binding and leads to lattice constants smaller than experimental. Keeping all the above in mind, we ?nd that the agreement between our computations and measurements for test cases Si and GaAs is reasonably good (see Table II for details). It is also important to note that values of Cijk calculated both from strain-energy and strain-stress relations are consistent with each other. As a cross-check we additionally veri?ed our approach by calculating second-
5 order elastic moduli for GaAs with the aid of the MedeA package.42 It uses its own methodology of calculating cij on the basis of stress computed by the VASP code.32 We obtained values c11 = 99 GPa, c12 = 41 GPa, c44 = 51 GPa, which are in agreement with the results given in Table II. Next interesting issue is to examine for which range of deformations the third-order e?ects really matter. In Fig. 4 we compare energy and stress for the particular deformation ηB in GaN crystal with energy and stress values obtained within linear and nonlinear elasticity theories. One can clearly see that the linear approach is not suf?cient for strains larger than approximately 2.5%. It is also worth noting that for all studied semiconductors and examined range of deformations (i.e., with Lagrangian strains up to 8%) including the terms up to third-order in energy expansion [Eq. (3)] su?ced to obtain good agreement with DFT results. It is also important to note that a quadratic term in ? in the expression for Lagrangian strain η [see Eq. (4)] is usually neglected when the second-order elastic constants are determined. For the third-order elastic constants, such omission is completely unjusti?ed. For example, the approximation η ≈ ? leads to the following third-order wrong wrong elastic constants for Si, C111 = ?256 GPa, C112 = wrong wrong ?375 GPa, C144 = 94 GPa, C155 = ?130 GPa, wrong wrong C123 = ?105 GPa, and C456 = ?2 GPa, which show signi?cant disagreement with the results obtained without the aforementioned simpli?cation (compare results in Table II). As one would expect, the second-order elastic constants remain virtually una?ected by the approximation η ≈ ?, now being c11 = 150 GPa, c12 = 62 GPa, and c44 = 73 GPa.
0.025 fB(ε) nonlinear elasticity linear elasticity
(a)
0.020 Energy [eV/? ]
3
0.015
0.010
0.005
0.000 ?0.08 ?0.06 ?0.04 ?0.02 0 0.02 0.04 0.06 0.08 Strain ε 30 20 10 0 ?10 ?20 ?30 ?0.08 ?0.06 ?0.04 ?0.02 0 Strain ε 0.02 0.04 0.06 0.08
(b)
gB(ε) nonlinear elasticity linear elasticity
FIG. 4: Energy (a) and stress component t3 (b) as a function of linear strain parameter ? for particular deformation ηB for GaN cubic crystal. Full points denote results of DFT computations, solid and dashed lines indicate the curves obtained from nonlinear and linear elasticity theory respectively.
las are given below1 c′ 11 = ?
IV. RELATION TO PRESSURE DEPENDENT ELASTIC CONSTANTS
Lagrangian Stress [GPa]
c′ 12 c′ 44
In the case of materials under large hydrostatic pressure it is useful to describe the nonlinear elastic properties using the concept of pressure dependent elastic constants cij (P ). For many applications, it is su?cient to consider only terms linear in the external hydrostatic pressure cP (P ) ≈ c11 + c′ 11 P, 11 cP (P ) ≈ c12 + c′ 12 P, 12 cP (P ) ≈ c44 + c′ 44 P, 44
2C112 + C111 + 2c12 + 2c11 , 2c12 + c11 C123 + 2C112 ? c12 ? c11 = ? , 2c12 + c11 2C155 + C144 + c44 + 2c12 + c11 = ? . 2c12 + c11
(10)
(9)
with pressure derivatives c′ ij being material parameters. Naturally, the information about c′ ij can be recovered from third-order elastic constants. The necessary formu-
Results for pressure derivatives c′ ij calculated on the basis of our estimates for second- and third-order elastic constants are shown in Tables IV and V. Table IV provides comparison with experimental results for Si38 and GaAs.39 The agreement is very good and shows that the results from the strain-energy relation reproduce the experimental values slightly better than ?ndings based on strain-stress formula. Table V contains values of c′ ij for zinc-blende nitrides and compares the present calculation with our previous work.12 In Ref. 12, the following approach for the determination of pressure dependence of the second-order elastic constants has been used. First, the hydrostatic strain (corresponding to the external pressure P ) has been ap-
6
TABLE II: Comparison of the calculated second- and third-order elastic constants for Si and GaAs with the experimental values and previous calculations. All data are in GPa. Present results strain-energy strain-stress Si c11 c12 c44 C111 C112 C144 C155 C123 C456 C144 + 2C155 GaAs c11 c12 c44 C111 C112 C144 C155 C123 C456
a Reference b Reference c Reference d Reference e Reference f Reference g Reference h Reference i Reference j Reference
Previous calculations
Experiment
153 65 73 -698 -451 74 -253 -112 -57 -430 100 49 52 -561 -337 -14 -244 -83 -22
153 57 75 -687 -439 72 -252 -92 -57 -432 99 41 51 -561 -318 -16 -242 -70 -22
159 61 85 -750 -480
a a a a a
0 -80 -580 126 55 61 -600 -401 10 -305 -94 -43
a a a
-880 -515 74 -385 27 -40 -696
d d d d d d d
-834 -531 -95 -296 -2 -7 -687
e e e e e e e
167 c 65 c 80 c -825 f -451 f 12 f -310 f -64 f -64 f -608 f 113 g 57 g 60 g -620 j -392 j 8j -274 j -62 j -43 j
b b b b b b b b b
-675 -402 -70 -320 -4 -69
h h h h h h
-622 -387 2 -269 -57 -39
i i i i i i
8 (LDA). 9 (LDA). 33 (T = 73K). 34 (T = 4K). 34 (T = 298K). 29 (T = 298K). 31 (extrapolation to T = 0K). 35 (T = 298K). 36 (T = 298K). 37 (T = 298K).
plied to a crystal, and then the crystal has been additionally deformed to determine the pressure dependent elastic constants. The DFT results for the total elastic energy combined with the strain-energy relation have enabled us to determine cij (P ) as well as c′ ij . We would like to stress that the additional nonin?nitesimal strain has not always been trace-free just leading to a spurious hydrostatic component that has modi?ed external hydrostatic pressure. Therefore, we believe that the approach employed in the present paper is not only more direct, but also slightly more accurate. The discrepancies between our present and previous results can also be partly ascribed to the methodological di?erences, such as di?erent exchange-correlation functional used and slightly different calculation parameters (Brillouin zone sampling, energy cuto?s etc.).
¨ V. RELATION TO GRUNEISEN CONSTANTS OF LONG-WAVELENGTH ACOUSTIC MODES
The mode Gr¨ neisen constants constitute a group of u important coe?cients, which characterize anharmonic properties of crystals. These quantities are frequently encountered in theory of phonons and in the description of thermodynamical properties of solids. The mode Gr¨ neisen constants are de?ned as follows: u γ(q, j) = ? ? ln ω(q, j) V ?ω(q, j) =? , ? ln V ω(q, j) ?V (11)
where ω denotes the frequency of phonon with wave vector q and polarization vector j. V stands here for volume of the crystal. On the basis of continuum limit, one may express mode Gr¨ neisen constants for long-wavelength acoustic u modes in terms of second- and third-order elastic con-
7
TABLE III: Theoretical predictions for the third-order elastic constants of zinc-blende nitrides - AlN, GaN, and InN. The second-order elastic constants are included and compared with previous calculations. All data are in GPa. Present results strain-energy AlN c11 c12 c44 C111 C112 C144 C155 C123 C456 GaN c11 c12 c44 C111 C112 C144 C155 C123 C456 InN c11 c12 c44 C111 C112 C144 C155 C123 C456
a Reference
Previous calculations a strain-stress 282 149 179 -1073 -965 57 -757 -61 -9 252 129 147 -1213 -867 -46 -606 -253 -49 159 102 78 -756 -636 13 -271 -310 15 267 141 172
284 167 181 -1070 -1010 63 -751 -78 -11 255 147 148 -1209 -905 -45 -603 -294 -48 160 115 78 -752 -661 16 -268 -357 14
12 (GGA).
252 131 146
149 94 77
stants. The necessary expressions used here have been given by Mayer and Wehner.40 The results for γ(q, j) obtained from our strain-energy estimates of elastic moduli are given in Table VI. Comparison with the experimental data available for Si shows that results calculated by us often di?er signi?cantly from experimental ?ndings. The discrepancy is particularly pronounced for transverse modes (i.e., γ((?, 0, 0), TA) = γ((?, ?, 0), TAz ) and γ((?, ?, 0), TAxy )) for which the magnitudes of Gr¨ neisen constants are u much smaller than for longitudinal modes. In our opinion, this indicates that γ(q, j) are quite sensitive to inaccuracies in Cijk values. Therefore, one has to treat our prediction for mode Gr¨ neisen constants in zinc-blende u nitrides rather as a quite crude approximation. Neverthe-
less, it could be an interesting subject of further studies to compare the above results with ab initio phonon calculations performed via density functional perturbation theory. More detailed experimental studies for a broader range of materials could also shed more light on the value of the presented theoretical predictions.
VI.
CONCLUSIONS
We have presented a detailed ab initio study of thirdorder elastic constants Cijk for selected semiconductors - Si, GaAs, and zinc-blende nitrides AlN, GaN, and InN. Even though third-order e?ects are very subtle, we showed that it is possible to estimate them by means of
8
TABLE IV: Pressure derivatives of second-order elastic constants for Si and GaAs calculated on the basis of Eqs. (10) . For comparison experimental ?ndings are included. Present results strain-energy strain-stress Si c′ 11 c′ 12 c′ 44 GaAs c′ 11 c′ 12 c′ 44
a Reference
Experiment
4.09 4.34 0.27 4.71 4.56 1.27
4.30 4.43 0.34 5.06 4.67 1.48
4.19 4.02 0.80 4.63 4.42 1.10
a a a
for Si and GaAs with available experimental ?ndings. The agreement is reasonable, however, particularly for moduli of smaller magnitude (e.g., for examined cases typically C144 and C456 ) relative di?erences are signi?cant. In our opinion, they can be ascribed to three main factors: shortcomings of GGA-DFT theory, lack of temperature e?ects in our calculations (experimental results for Cijk are usually obtained in room temperature), and
TABLE VI: Gr¨neisen constants γ(q, j) for long-wavelength u acoustic modes. Experimental results for Si were given in Ref. 40 on the basis of ultrasound measurements data from Ref. 29. Theoretical prediction for γ(q, j) were calculated on the basis of strain-energy values for cij and Cijk . Experiment Si q = (?, 0, 0) γ(LA) γ(TA) q = (?, ?, 0) γ(LA) γ(TAxy ) γ(TAz ) q = (?, ?, ?) γ(LA) 1.081 0.973 1.056 1.214 1.173 1.109 -0.049 0.324 0.999 -0.301 0.006 1.066 -0.684 0.423 1.226 -0.613 0.459 1.218 -1.771 -0.055 1.108 0.324 Theory Si 1.098 0.006 Theory GaN 1.279 0.459
b b b
AlN 1.115 0.423
InN 1.415 -0.055
38 (T=4K). b Reference 39 (T=298K).
TABLE V: Prediction of pressure derivatives of second-order elastic constants for zinc-blende nitrides AlN, GaN, and InN calculated on the basis of Eqs. (10). For comparison, results of previous calculations employing di?erent methodology are included. Present results Previous calculations strain-energy strain-stress AlN c′ 11 c′ 12 c′ 44 GaN c′ 11 c′ 12 c′ 44 InN c′ 11 c′ 12 c′ 44
a Reference a
3.53 4.12 1.03 4.03 4.56 1.01 3.89 5.00 0.13
12.
3.68 4.17 1.20 4.28 4.64 1.18 4.15 5.08 0.24
5.21 4.26 1.69 4.17 3.50 1.12 4.58 4.37 0.66
measurement uncertainties. We have also underlined the relation of third-order elastic constants to other anharmonic properties. On the basis of the ab initio results for Cijk , we have computed the pressure derivatives of second-order elastic moduli and provided rough estimations for Gr¨ neisen constants of long-wavelength acousu tic modes. We believe that DFT estimates of third-order elastic constants can be a very useful tool in modeling semiconducting nanostructures, in which nonlinear effects often play an important role.
Acknowledgments
density functional theory on the GGA level. We have used two approaches involving either strain-energy or strain-stress relations, obtaining consistent results from both of them. To benchmark the reliability of the presented method, we have compared our theoretical results
This work was partly supported by the Polish State Committee for Scienti?c Research, Project No. 1P03B03729. One of the authors (M.L.) wishes to acknowledge useful discussion with Alexander Mavromaras from Materials Design.
1 2
3
F. Birch, Phys. Rev. 71, 809 (1947). F. Murnaghan, Finite Deformation of an Elastic Solids (John Willey and Sons, 1951). S. Bhagavantam, Crystal Symmetry and Physical Properties (Academic Press, 1966).
4
5 6 7
R. Thurston and K. Brugger, Phys. Rev. 133, A1604 (1964). K. Brugger, Phys. Rev. 133, A1611 (1964). Y. Hiki, Ann. Rev. Mater. Sci. 11, 51 (1981). P. N. Keating, Phys. Rev. 149, 674 (1966).
9
8
9 10
11
12
13
14
15
16
17
18
19
20
21 22
23
24
O. H. Nielsen and R. M. Martin, Phys. Rev. B 32, 3792 (1985). J. S¨rgel and U. Scherz, Eur. Phys. J. B 5, 45 (1998). o R. Kato and J. Hama, J. Phys.: Condens. Matter 6, 7617 (1994). S. P. Lepkowski and J. A. Majewski, Solid State Commun. 131, 763 (2004). S. P. Lepkowski, J. A. Majewski, and G. Jurczak, Phys. Rev. B 72, 245201 (2005). K. Shimada, T. Sota, K. Suzuki, and H. Okumura, Jpn. J. Appl. Phys. 37, L1421 (1998). G. Bester, X. Wu, D. Vanderbilt, and A. Zunger, Phys. Rev. Lett. 96, 187602 (2006). M. D. Frogley, J. R. Downes, and D. J. Dunstan, Phys. Rev. B 62, 13612 (2000). S. W. Ellaway and D. A. Faux, J. Appl. Phys. 92, 3027 (2002). B. S. Ma, X. D. Wang, F. H. Su, Z. L. Fang, K. Ding, Z. C. Niu, and G. H. Li, J. Appl. Phys. 95, 933 (2004). J. W. Luo, S. S. Li, J. B. Xia, and L. W. Wang, Phys. Rev. B 71, 245315 (2005). G. Vaschenko, C. S. Menoni, D. Patel, C. N. Tom?, e B. Clausen, N. F. Gardner, J. Sun, W. G¨tz, H. M. Ng, o and A. Y. Cho, Phys. Stat. Sol. (b) 235, 238 (2003). S. P. Lepkowski and J. A. Majewski, Phys. Rev. B 74, 035336 (2006). G. Kresse and J. Hafner, Phys. Rev. B 47, 558 (1993). G. Kresse and J. Furthm¨ller, Phys. Rev. B 54, 11169 u (1996). G. Kresse and J. Furthm¨ller, Comput. Mat. Sci. 6, 15 u (1996). P. E. Bl¨chl, Phys. Rev. B 50, 17953 (1994). o
25 26
27
28
29
30
31 32 33 34 35
36
37
38 39
40
41
42
G. Kresse and D. Joubert, Phys. Rev. B 59, 1758 (1999). J. P. Perdew, K. Burke, and M. Ernzerhof, Phys. Rev. Lett. 77, 3865 (1996). J. P. Perdew, K. Burke, and M. Ernzerhof, Phys. Rev. Lett. 78, 1396 (1996). H. J. Monkhorst and J. D. Pack, Phys. Rev. B 13, 5188 (1976). H. J. McSkimin and P. Andreatch, J. Appl. Phys. 35, 3312 (1964). N. A. W. Holzwarth, G. E. Matthews, R. B. Dunning, A. R. Tackett, and Y. Zeng, Phys. Rev. B 55, 2005 (1997). J. S. Blakemore, J. Appl. Phys. 53, R123 (1982). Y. LePage and P. Saxe, Phys. Rev. B 65, 104104 (2002). H. J. McSkimin, J. Appl. Phys. 24, 988 (1953). J. Philip and M. Breazeale, J. Appl. Phys. 52, 3383 (1981). J. Drabble and A. Brammer, Solid State Commun. 4, 467 (1966). H. J. McSkimin and P. Andreatch, J. Appl. Phys. 38, 2610 (1967). Y. Abe and K. Imai, Jpn. J. Appl. Phys. 25 Suppl 25-1, 67 (1986). A. G. Beattie and J. Schirber, Phys. Rev. B 1, 1548 (1970). H. J. McSkimin, A. Jayaraman, and P. Andreatch, J. Appl. Phys. 38, 2362 (1967). A. Mayer and R. Wehner, Phys. Stat. Sol. (b) 126, 91 (1984). In older texts concerning nonlinear elasticity di?erent de?nitions of Cijk may be encountered. In this paper we follow the convention proposed in Ref. 5 which is now a standard approach. See http://www.materialsdesign.com for details about the software.