First-principles study of structural, elastic, electronic and lattice dynamic properties of AsxPyN1−x−yB quaternary alloys

Computational Materials Science 48 (2010) 94–100

Contents lists available at ScienceDirect

Computational Materials Science journal homepage: www.elsevier.com/locate/commatsci

First-principles study of structural, elastic, electronic and lattice dynamic properties of AsxPyN1xyB quaternary alloys B. Ghebouli a,*, M.A. Ghebouli b, M. Fatmi c, S.I. Ahmed d a

Laboratoire d’Etudes de Surfaces et Interfaces des Matériaux Solides (LESIMS), Université Ferhat Abbas, Sétif 19000, Algeria Département de Physique, Centre Universitaire, Bordj-Bou-Arréridj 34000, Algeria c Laboratoire de Physique et Mécanique des Matériaux Métalliques (LP3M), Université Ferhat Abbas, Sétif 19000, Algeria d Département de Physique, Faculté des Sciences, Université Ain Shams, Caire, Egypt b

a r t i c l e

i n f o

Article history: Received 2 August 2009 Received in revised form 8 November 2009 Accepted 6 December 2009 Available online 8 January 2010 Keywords: Electronic structure Optical properties Lattice matched and mismatched alloys Pseudopotential calculations

a b s t r a c t Information on the energy band gaps, the lattice parameters and the lattice matching to available substrates is a prerequisite for many practical applications. A pseudopotential plane-wave method as implemented in the ABINIT code is used to the AsxPyN1xyB quaternary alloys lattice matched to BP substrate to predict their energy band gaps and lattice dynamic properties. The range of compositions for which the alloy is lattice matched to BP is determined. Very good agreement is obtained between the calculated values and the available experimental data. The compositional dependence of direct and indirect band gaps has been investigated. We study the variation of elastic constants, the optical phonon frequencies (xTO and xLO), the high-frequency dielectric coefficient e(1) and the born effective charge Z* with P concentration. Ó 2009 Elsevier B.V. All rights reserved.

1. Introduction Alloying in the group-V Boron has diversified the properties of semiconductor materials, enabling their production with wide band gaps. These semiconductors are used to produce commercially important high-performance electronic and optoelectronic devices and systems, such as light emitting devices covering many regions of the visible spectrum. The availability of the quaternary alloy AsxPyN1xyB/BP structure permits an extra degree of freedom by allowing independent control of the band gap, Eg, and the lattice constant, a0. Tailoring of these compounds could lead to new semiconductor materials with desired band gaps over a continuous broad spectrum of energies [1–9]. These alloys have emerged as interesting materials for device applications because of their energy band structure and lattice parameters. Furthermore, the use of quaternary films enhances the experimental capability for investigating the effects of strain and piezoelectric fields in quantum wells [10]. In order to help understand and control the materials and device properties, we have carried out a theoretical study of the energy band gaps and elastic properties of the quaternary AsxPyN1xyB, lattice matched to BP substrate. The electronic structure of the quaternary alloy system is calculated using the pseudopotential plane-wave method as implemented in the ABINIT code.

* Corresponding author. E-mail address: [email protected] (B. Ghebouli). 0927-0256/$ - see front matter Ó 2009 Elsevier B.V. All rights reserved. doi:10.1016/j.commatsci.2009.12.008

Calculations are performed over the entire composition range of x and y. This paper is organized as follows: A brief introduction was given in Section 1. The computational method adopted for the calculations is described in Section 2. We present our results and compare them to the available experimental data and others previously published theoretical results in Section 3. Then, our work is summarized in Section 4. 2. Computational method The calculations were performed within the GGA to the DFT, using the pseudopotential plane-wave method as implemented in the ABINIT code [11]. ABINIT computer code is a common project of the Université Catholique of Louvain, Corning Incorporated, and other contributors. Only the outermost electrons of each atom were explicitly considered in the calculation. The effect of the inner electrons and the nucleus (the frozen core) was described within a pseudopotential scheme. We used the Hartwigzen–Goedecker– Hutter scheme [12] to generate the norm-conserving nonlocal pseudopotentials, which results in highly transferable and optimally smooth pseudopotentials. A plane-wave basis set was used to solve the Kohn–Sham equations in the pseudopotential implementation of the DFT–GGA. The Brillouin zone integrations were replaced by discrete summations over a special set of k-points, using the standard k-point technique of Monkhorst and Pack [13] where the k-point mesh used is (8  8  8). The plane-wave energy cutoff to expand the wave functions is set to be 90 Hartree

B. Ghebouli et al. / Computational Materials Science 48 (2010) 94–100

95

(1 Hartree = 27.211396 eV). For the treatment of the disordered ternary alloy, we used the VCA [14], in which the alloy pseudopotentials are constructed within a first-principles VCA scheme. Elemental ionic pseudopotentials of AsB and BPN are combined to construct the virtual pseudopotential of AsxPyN1xyB. Recently, Marques et al. [15] have reported a linear behavior of the lattice parameter of AlxGayIn1xyN as a function of the composition x, y. Thus, the Vegard’s rule has been assumed for the calculation of the lattice constant of quaternary alloys under study.

aVCA ¼ xaAsB þ yaBP þ ð1  x  yÞaBN 3. Results 3.1. Structural properties The semiconducting XB compounds are crystallized in the zincblende structure. The space group is F-43 m. The B atom is located at the origin and the X atom is located at (1/4, 1/4, 1/4). The equilibrium lattice parameter is computed from the structural optimization, using the Broyden–Fletcher–Goldfarb–Shanno minimization [16–19]. The lattice matching conditions for AsxPyN1xyB quaternary systems on the BP substrate is: y = 1  1.2588x. Fig. 1. The results of lattice parameter a0 are reported for different composition rates (x, y) in Table 1. Fig. 2 shows the lattice parameter plotted versus P fraction (y). The deviation from the linear dependence is distinct. An analytical relation for the compositional dependence of AsxPyN1xyB lattice parameter is given by quadratic fit:

a ¼ 4:7073  0:129y  0:042y2

ð1Þ

3.2. Elastic properties The elastic constants are important parameters that describe the response to an applied macroscopic stress. In Table 1, we show the calculated elastic constants, namely (C11, C12 and C44), bulk modulus B, shear wave modulus Cs and Kleinman parameter f of AsxPyN1xyB/BP structure at zero pressure and for various compositions (y) in the range (0–1). Also given for comparison are the available experimental data and other works. To the best of our knowledge, no experimental data have been reported so far for the elastic constants of AsxPyN1xyB in the 0 < y < 1 composition

Fig. 2. Relaxed lattice parameter of AsxPyN1xyB/BP structure as a function of P our data, quadratic fit to our data and linear fit. composition (y),

range. In Fig. 3, we depict the composition dependence of the elastic constants (C11, C12 and C44) and the bulk modulus B of AsxPyN1xyB. We observe a linear dependence in all curves in the considered range of composition. It is easy to observe that the elastic constants Cij and bulk modulus B increase when the composition (y) is enhanced. However, the enhancement in C12 is weaker. The mechanical stability requires the elastic constants satisfying the well-known born stability criteria [20]:

9 8 1 > = < K ¼ 3 ðC 11 þ 2C 12 þ PÞ > 0 > G ¼ 12 ðC 11  C 12  2PÞ > 0 > > ; : G0 ¼ C 44  P > 0

ð2Þ

In Fig. 4, we show the dependence of stability criteria of AsxPyN1xyB compounds with P concentration (y) at zero pressure. From our calculated Cij, it is known that these compounds are mechanically stable. From the theoretical elastic constants, we computed the elastic wave velocities. The single-crystal elastic wave velocities in different directions are given by the resolution of the Cristoffel equation [21]:

ðC ijkl  nj  nk  qv 2  dil Þul ¼ 0

ð3Þ

where Cijkl is the single-crystal elastic constant tensor, n is the propagation direction, q is the density of material, u is the wave polarization and m is the wave velocity. The solutions of this equation are of two types: a longitudinal wave with polarization parallel to the direction of propagation vl and two shear waves vT1 and vT2 with polarization perpendicular to n. Another important parameter is the elastic anisotropy factor, A, which gives a measure of the anisotropy of the elastic wave velocity in a crystal. In a cubic crystal, the elastic anisotropy factor is given by [22]:

A¼

Fig. 1. Lattice matching conditions for AsxPyN1xyB quaternary systems on the BP substrate.

2C 44 þ C 12 1 C 11

ð4Þ

which is zero for an isotropic material. The variation of the longitudinal-wave mode speed (VL) and transverse-wave mode speed (VT) in AsxPyN1xyB/BP structure propagating in the [1 0 0], [1 1 0] and [1 1 1] directions at zero pressure for various P composition (y) in the range (0–1) are shown in Fig. 5. It seen that, longitudinal waves are fastest along [1 1 1] and shear waves are slowest along [1 1 0] for AsxPyN1xyB/BP structure which has a positive elastic anisotropy factor. The curves presented in these figures are the best

96

B. Ghebouli et al. / Computational Materials Science 48 (2010) 94–100

Table 1 The calculated lattice constant a0, bulk modulus B and elastic constants C11, C12 and C44 of AsxPyN1xyB compounds at zero pressure.

a b c d e f g h

Parameters

References

x=0 y=1

x = 0.2 y = 0.74

x = 0.4 y = 0.496

x = 0.6 y = 0.244

x = 0.8 y=0

a0 (Å)

This work Experiment Other works

4.5357 4.543a 4.546b

4.5879

4.6328

4.6728

4.7077

C11 (GPa)

This work Experiment Other works

340.89 315d 358.9e

328.82

315.24

302.37

288.34

C12 (GPa)

This work Experiment Other works

71.16 100d 80.6e

70

68.26

66.39

49.07

C44 (GPa)

This work Experiment Other works

204.18 160d 202f

195.26

186.12

177.9

169.63

B (GPa)

This work Experiment Other works

161.07 152a 162c

156.27

150.58

145.05

138.78

Cs

This work Other works

134.86 130.95g

129.41

123.49

117.99

119.63

f

This work Other works

0.359 0.375h

0.364

0.367

0.37

0.321

[28]. [29]. [30]. [31]. [32]. [33]. [34]. [35].

Fig. 3. Composition dependence of the elastic constants and the bulk modulus in AsxPyN1xyB/BP structure.

fit of the calculated data. The analytical expressions of VL and VT are as follows: 2 V 100 L ðm=sÞ ¼ 7852:47 þ 2016:19y þ 827:22y 2 V 100 T ðm=sÞ ¼ 6023:2 þ 1520:33y þ 733:84y 2 V 111 L ðm=sÞ ¼ 8834:55 þ 2227:3y þ 996:75y 2 V 111 T ðm=sÞ ¼ 5299:54 þ 1364:99y þ 616:33y

V 110 L ðm=sÞ V 110 T1 ðm=sÞ V 110 T2 ðm=sÞ

2

¼ 8599:55 þ 1176:37y þ 956:85y ¼ 4897:77 þ 1280:59y þ 549:1y2 2

¼ 6023:2 þ 1520:33y þ 733:84y

ð5Þ

Fig. 4. Stability criteria of AsxPyN1xyB/BP structure as a function of P concentration (y) at zero pressure.

One can note that longitudinal-wave mode speed and transverse-wave mode speed increase monotonously with increasing P composition y in the range (0–1). Once the elastic constants are determined, we would like to compare our results with experiments, or predict what experiment would yield for the elastic constants. For cubic crystal, the isotropic bulk modulus B and the shear constant Cs are given exactly by:

B¼

C 11 þ 2C 12 3

and C s ¼

1 ðC 11  C 12 Þ 2

ð6Þ

Another important parameter is the Kleinman parameter, f, which describes the relative positions of the cation and anion sublattices under volume-conserving strain distortions for which

97

B. Ghebouli et al. / Computational Materials Science 48 (2010) 94–100

since hD may be estimated from the average sound velocity the following equation [27]:

hD ¼

vm by

 1=3 h 3n mm kB 4pV a

ð9Þ

where h is Plank’s constant, kB Boltzmann’s constant and Va the atomic volume. The average sound velocity in the polycrystalline material is given by [21]:

mm ¼

  1=3 1 2 1 þ 3 2 3 ml mt

ð10Þ

where vl and vt are the longitudinal and transverse sound velocity of an isotropic aggregate obtained using the shear modulus G and the bulk modulus B from Navier’s equation [26]:

m1 ¼ Fig. 5. The longitudinal-wave mode speed (VL) and transverse-wave mode speed (VT) in AsxPyN1xyB/BP structure propagating in the [1 0 0], [1 1 0] and [1 1 1] directions at zero pressure for various P composition (y).

positions are not fixed by symmetry. It is known that a low value of f implies there is a large resistance against bond bending or bondangle distortion and vice versa [23,24]. We used the following relation [25]:

f¼

C 11 þ 8C 12 7C 11 þ 2C 12

ð7Þ

We also calculated Young’s modulus E and Poisson’s ratio m which are frequently measured for polycrystalline materials when investigating their hardness. These quantities are related to the bulk modulus and the shear modulus by the following equations [26]:

(

9BG E ¼ 3BþG

) ð8Þ

m ¼ 3BE 6B

The shear moduli GR and GV, Young’s modulus E and Poisson’s ratio m for AsxPyN1xyB/BP structure, calculated from the elastic constants are listed in Table 2. 3.3. Calculation of Debye temperature Having calculated Young’s modulus E, bulk modulus B and shear modulus G, one can calculate the Debye temperature, which is an important fundamental parameter closely related to many physical properties such as elastic constants, specific heat and melting temperature. At low temperature, the vibrational excitations arise solely from acoustic mode. Hence, at low temperature the Debye temperature calculated from elastic constants is the same as that determined from specific heat measurements. One of the standard methods to calculate the Debye temperature hD is from elastic data,

 1=2 3B þ 4G 3q

and

mt ¼

 1=2 G

ð11Þ

q

The calculated sound velocities and Debye temperature as well as the density for AsxPyN1xyB/BP structure are given in Table 3. One can note that Debye temperature and bulk modulus increase with increasing P composition (y). The Debye temperature increase when the bulk modulus is enhanced. 3.4. Electronic properties The band structure of AsxPyN1xyB lattice matched to BP was calculated through the high-symmetry points U, X and L in the Brillouin zone. We show for example the band structures of the parents alloys BP and As0.8N0.2B in Fig. 6. It seen that there is an indirect band gaps (C–X) for all composition (y). The direct band gaps (U–U, L–L, X–X) and indirect band gaps (U–L, U–X) of

Table 3 The calculated density q, the longitudinal, transverse and average sound velocity vl, vl and vm calculated from polycrystalline elastic modulus, and the Debye temperatures hD calculated from the average sound velocity of AsxPyN1xyB compounds at zero pressure. References BP As0.2p0.7482N0.0518B As0.4p0.4964N0.1036B As0.6p0.2447N0.1553B As0.8p0.007N0.193B

q

ml

mt

mm

(g cm3)

(m s1)

(m s1)

(m s1)

2.97476 3.41804 3.84783 4.26512 4.6898

11473.69 10502.87 9687.62 9010.12 8393.52

7623.97 6960.07 6405.92 5948.14 5535.65

8338.02 7614.56 7010.22 6510.61 6059.86

Table 2 The calculated shear moduli GR and GV, Young’s modulus E and Poisson’s ratio m of AsxPyN1xyB compounds at zero pressure. Material

Shear modulus, GR (GPa)

Shear modulus, GV (GPa)

Young’s modulus, E (GPa)

Poisson’s ratio, m

BP As0.2p0.7482N0.0518B As0.4p0.4964N0.1036B As0.6p0.2447N0.1553B As0.8p0.007N0.193B

169.36 162.23 161.06 153.93 140.78

176.45 168.92 154.73 147.86 146.64

382.02 367.08 351.01 336.13 320.5

0.147 0.108 0.110 0.113 0.115

Fig. 6. Band structure of parents alloys BP and As0.8N0.2B.

hD (K) 809.81 747.60 741.84 713.02 684.99

98

B. Ghebouli et al. / Computational Materials Science 48 (2010) 94–100

AsxPyN1xyB at zero pressure are shown in Table 4. The P composition energy variations of the conduction band edges at U, X and L with respect to the top of the valence band are obtained. Results are plotted in Fig. 7, which shows the dependence of direct and indirect band-gap energies, in the lattice matched AsxPyN1xyB/ BP structure, on the P content. Note that both direct and indirect band gaps namely, (L–L) and (U–L) increase monotonically with increasing the composition (y). This is not the case for the direct and indirect band gap (U–U), (X–X) and (C–X) which change nonmonotonically. In order to provide analytical expressions for the energy band gaps of the AsxPyN1xyB quaternary system which could be of practical use and of easy access, we obtained the following analytical expressions by fitting our energy gaps data using a least squares procedure:

AsxPyN1xyB/BP at various compositions (y) are listed in Table 5. Also shown for comparison are the existing experimental and theoretical data in the literature. The agreement between our results and experiment as regards xLO and xTO of BP is better than 3% and 0.3%. No comparison has been made in the 0 < y < 1 composition range, as there are no known data available to date to the best of our knowledge. The variation of the frequencies of the LO and TO phonons in AsxPyN1xyB/BP quaternary alloys as a function of the P concentration (y) is plotted in Fig. 8, one obtains for xLO and xTO:

xLO ¼ 702:97 þ 43:6y þ 56:4y2 xTO ¼ 701:26 þ 21:99y þ 53:16y2

ð13Þ

ð12Þ

These expressions may be used to predict the LO and TO phonon frequencies for any P concentration (y) in the range (0–1) for AsxPyN1xyB/BP structure. Note also that both xLO and xTO change non-linearly as the composition (y) increases. For the same concentration (y), xLO is greater than xTO. A good knowledge of the full electronic structure is an essential feature in order to get the best understanding of the optical properties of semiconductors. Thus, the refractive index which is essential in the design of heterostructure lasers, in optoelectronic

The numerically longitudinal and transversal optical phonon frequencies, referred to as xLO and xTO, respectively of

Table 5 Born effective charge, the longitudinal optical (LO) and transverse optical (TO) phonon frequencies, refractive index, the high-frequency and static dielectric constants of AsxPyN1xyB for different P concentration (y).

ECC

¼ 3:42 þ 0:14y  0:14y

2

X

EC ¼ 1:49 þ 0:44y  0:57y2 ELC ¼ 2:84 þ 0:8y þ 0:1y2 ELL EXX

2

¼ 4:67 þ 0:93y  0:09y

¼ 5:61 þ 0:61y  0:75y2

3.5. Lattice dynamic properties

Material Table 4 Some direct band gaps (U–U, L–L, X–X) and indirect band gaps (U–L, U–X) of AsxPyN1xyB at zero pressure. Elements

U–U (eV)

L–L (eV)

X–X (eV)

U–L (eV)

U–X (eV)

BP This work Experiment Others As0.2p0.7482N0.0518B As0.4p0.4964N0.1036B As0.6p0.2447N0.1553B As0.8p0.007N0.193B

3.42 5b 3.40a 3.45 3.44 3.44 3.43

5.51 8b 5.4a 5.31 5.11 4.9 4.68

5.73 6.9b 5.99c 5.72 5.69 5.59 5.45

3.76

1.37

3.50 3.27 3.04 2.83

1.5 1.59 1.56 1.54

BP This work Experiment Others As0.2p0.7482N0.0518B As0.4p0.4964N0.1036B As0.6p0.2447N0.1553B As0.8p0.007N0.193B a b c d

a b c

[36]. [37]. [27].

Fig. 7. Direct (C–C), (L–L) and (X–X) and indirect (C–L) and (C–X) band-gap energies in AsxPyN1xyB/BP quaternary alloys as a function of P concentration (y).

e

Z*

0.6159 1.34b 0.66a 0.5366 0.4538 0.3882 0.3234

xLO

xTO

(cm1)

(cm1)

803.16 828.9b 837c 766.74 738.88 717.01 702.92

796.6 799b 813c 761.97 735.59 714.69 701.21

n

e(0)

e(1)

3.084 3e

9.673 11d

3.084 3.093 3.114 3.145

9.632 9.656 9.762 9.94

9.516 9.61d 9.17c 9.512 9.57 9.699 9.894

[38]. [39]. [40]. [41]. [42].

Fig. 8. Frequencies of the longitudinal optical (LO) and transversal optical (TO) phonons in AsxPyN1xyB/BP quaternary alloys as a function of P concentration (y).

B. Ghebouli et al. / Computational Materials Science 48 (2010) 94–100

devices as well as in solar cell applications has been calculated. Our results for various concentrations (y) are listed in Table 5. The knowledge of the refractive index makes it possible to proceed with the calculations of the dielectric constants. The high-frequency dielectric constant has been calculated and the obtained data are depicted in Table 5. The calculations in the present work were extended to include the static dielectric constant. We have calculated the values of e(0) at various P concentrations (y), for the lattice matched AsxPyN1xyB/BP structure. Our results are presented in Table 5. For the fact that both experimental and theoretical data regarding n, e(1), and e(0) for AsxPyN1xyB lattice matched to BP were not available, our results are predictions and may serve for reference. The variations of both static and high-frequency dielectric constants as a function of composition (y) are plotted in Fig. 9. The high-frequency dielectric and static dielectric constants decrease rapidly for P concentrations varying from 0% to about 0.84% then show a weaker enhancement. In order to provide analytical expressions for both e(1) and e(0) of the material of

99

interest, the data obtained by our calculations was found to fit best by a least squares procedure giving the following relations:

eð0Þ ¼ 9:93  0:86y þ 0:6y2 eð1Þ ¼ 9:89  0:91y þ 0:54y2

ð14Þ

These expressions may be useful for obtaining e(1) and e(0) for any concentration (y) in AsxPyN1xyB alloy lattice matched to BP. In Fig. 10, the absolute values of the calculated Born effective charge Z* are plotted as a function of the P concentration (y). The value of Z* increases monotonously with increasing P concentration (y). A quadratic fit to our data gives:

Z  ¼ 0:323  0:244y  0:048y2

ð15Þ

4. Conclusion In conclusion, we studied the composition dependence of energy band gaps and dielectric constants of AsxPyN1xyB quaternary alloys, lattice matched to BP substrate, by means of the pseudopotential plane-wave method as implemented in the ABINIT code. Based on the calculated lattice constants, we predicted the range of compositions for which AsPNB is lattice matched to BP. Very good agreement was obtained between the band gap predictions and available experimental data. Our findings indicate that the high-frequency and the static dielectric constants decrease monotonically with increasing P concentration (y) from 0 to 1. The expressions derived for direct and indirect energy band gaps, static and high-frequency dielectric constants as a function of P concentrations (y) are useful for tailoring electric and opto-electric devices based on cubic AsxPyN1xyB quaternary alloys. Due to the lack of experimental and theoretical data regarding the refractive index and high-frequency and static dielectric constants, our results are predictions and may serve as reference for future experimental work. Acknowledgment This work is supported by LESIMS Laboratory, University of Setif.

Fig. 9. Static and high-frequency dielectric constants as a function of P concentration (y) for the lattice matched AsxPyN1xyB/BP structure.

Fig. 10. Dependence of the born effective charge with P composition (y) for AsxPyN1xyB/BP structure. Solid line is quadratic fit to our data.

References [1] J.W. Orton, C.T. Foxon, Rep. Prog. Phys. 61 (1998) 1. and references therein. [2] S.C. Jain, M. Willander, J. Narayan, R. Van Overstraelen, J. Appl. Phys. 87 (2000) 965. and references therein. [3] K. Kassali, N. Bouarissa, Microelectron. Eng. 54 (2000) 277. [4] I. Vurgaftman, J.R. Meyer, J. Appl. Phys. 94 (2003) 3675. and references therein. [5] S. Nakamura, G. Fasol, Springer, Berlin, 1997. [6] S. Nakamura, Semiconduct. Sci. Technol. 14 (1999) R27. [7] P. Kung, M. Razegui, Opt. Rev. 8 (2000) 201. [8] M. Asif Khan, J.N. Kusnia, D.T. Olson, W.J. Schaff, J.W. Burm, M.S. Shur, Appl. Phys. Lett. 65 (1994) 1121. [9] C.I.H. Ashby, C.C. Willan, J. Han, N.A. Missert, P.P. Provencio, D.M. Follstaedt, G.M. Peake, L. Griego, Appl. Phys. Lett. 77 (2000) 3233. [10] F.G. McIntosh, K.S. Boutros, J.C. Roberts, S.M. Bedair, E.L. Piner, N.A. El-Masry, Appl. Phys. Lett. 68 (1996) 40. [11] The ABINIT Computer Code is a Common Project of the Universite Catholique of Louvain, Corning Incorporated and Other Contributors. . [12] C. Hartwigsen, S. Goedecker, J. Hutter, Phys. Rev. B 58 (1998) 3641. [13] H.J. Monkhorst, J.D. Pack, Phys. Rev. B 13 (1976) 5189. [14] L. Nordheim, Ann. Phys. (Leipzig) 9 (1931) 607. [15] M. Marques, L.K. Teles, L.M.R. Scolfaro, J.R. Leite, J. Furthmuller, F. Bechstedt, Appl. Phys. Lett. 83 (2003) 890. [16] C.G. Broyden, J. Inst. Math. Appl. 6 (1970) 222. [17] R. Fletcher, Comput. J. 13 (1970) 317. [18] D. Goldfarb, Math. Comput. 24 (1970) 23. [19] D.F. Shanno, Math. Comput. 24 (1970) 647. [20] G.V. Sin’ko, A. Smirnov, J. Phys.: Condens. Mat. 14 (2002) 6989. [21] O.L. Anderson, J. Phys. Chem. Solids 24 (1963) 909. [22] S.E. Lofland, J.D. Hettinger, A. Bryan, P. Finkel, S. Gupta, M.W. Barsoum, G. Hug, Phys. Rev. B 74 (2006) 74501.

100

B. Ghebouli et al. / Computational Materials Science 48 (2010) 94–100

[23] K. Kim, W.R.L. Lambrecht, B. Segal, Phys. Rev. B 50 (1994) 1502. [24] L. Kleinman, Phys. Rev. B 12 (1962) 2614. [25] W.A. Harrison, Electronic Structure and Properties of Solids, Dover, New York, 1989. [26] M.W. Barsoum, T. El-Raghi, W.D. Porter, H. Wang, S. Chakraborty, J. Appl. Phys. 88 (2000) 6313. [27] P. Wachter, M. Filzmoser, J. Rebisant, J. Physica B293 (2001) 199. [28] H. Xia, A.L. Ruoff, J. Appl. Phys. 74 (1993) 1660. [29] A. Zaoui, F. El-Haj Hassan, J. Phys.: Condens. Mat. 13 (2001) 253. [30] R. Mohammad, S. Katırcıoglu, J. Alloy. Compd. 478 (2009) 531. [31] M. Grimsditch, E.S. Zouboulis, J. Appl. Phys. 76 (2) (1994) 832. [32] S.Q. Wang, H.Q. Ye, Phys. Status Solidi B 240 (1) (2003) 45.

[33] P.R. Hernandez, M.G. Diaz, A. Munoz, Phys. Rev. B 51 (1994) 14705. [34] A. Bouhemadou, R. Khenata, M. Kharoubi, T. Seddik, Ali H. Reshak, Y. Al-Douri, Comp. Mater. Sci. 45 (2009) 474. [35] S.Q. Wang, H.Q. Ye, Phys. Status Solidi B 240 (2003) 45. [36] M.P. Surh, S.G. Louie, M.L. Cohen, Phys. Rev. B 43 (1991) 9126. [37] C.C. Wang, M. Cardona, A.G. Fischer, RCA Rev. 25 (1964) 159. [38] Kh. Bouamama, P. Djemia, N. Lebga, K. Kassali, High Pressure Res. 27 (2007) 1. [39] J.A. Sanjurjo et al., Phys. Rev. B 28 (1983) 4579. [40] D. Touat, M. Ferhat, A. Zaoui, J. Phys.: Condens. Mat. 18 (2006) 3647. [41] W. Wettling, J. Windscheif, Solid State Commun. 50 (1983) 33. [42] D.J. Moss, E. Ghahramani, J.E. Sipe, Phys. Status Solidi B 164 (1991) 578.