Gupta potentials for five HCP rare earth metals

Computational Materials Science 112 (2016) 75–79

Contents lists available at ScienceDirect

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

Gupta potentials for five HCP rare earth metals Xiaojie Li, Jie Fu, Ying Qin ⇑, Shengzhi Hao, Jijun Zhao ⇑ Key Laboratory of Materials Modification by Laser, Ion and Electron Beams (Dalian University of Technology), Ministry of Education, Dalian 116024, China

a r t i c l e

i n f o

Article history: Received 22 July 2015 Received in revised form 6 October 2015 Accepted 11 October 2015 Available online 19 November 2015 Keywords: Rare earth Gupta potential Elastic constants Phonon dispersion

a b s t r a c t Empirical Gupta-type potentials based on the second-moment approximation of tight-binding model have been developed for five hexagonal-closed-packed (hcp) rare earth metals: Er, Dy, Gd, Tb, Lu. The potentials can reproduce experimental cohesive energies, lattice constants and elastic constants of Er, Dy, Gd, Tb and Lu. The calculated bulk moduli, shear moduli, sound velocities and Debye temperatures are in reasonable agreement with the measured values. Vacancy formation energies and surface energies of these metals are also calculated and compared with previous theoretical results. Our potentials are able to simulate the phonon dispersions of Tb, Er, Dy, Gd, and Lu solids and reproduce Tb’s experimental data in detail. The present Gupta potentials would be useful for future simulations of elementary metals of Er, Dy, Gd, Tb, Lu and their alloys. Ó 2015 Elsevier B.V. All rights reserved.

1. Introduction Rare earth materials become known to the world by the discovery of ytterbite by Johan Gadolin in 1794. Rare earth metals were once scientifically precious, but now have commercial values with modern methods of separation and applications. With exotic optical, magnetic, electronic and chemical properties, rare earth elements play an important role in material science and industry, for example, improving the performance of steels and cast iron with a small amount of rare earth elements, working as catalysts in chemical engineering; their oxides being used as laser materials, superconductor materials and so on. Under ambient conditions, lutetium (Lu), gadolinium (Gd), terbium (Tb), erbium (Er), dysprosium (Dy) belong to hcp crystals. Lu is the hardest metal among rare earth solids and is mainly used for fundamental research. Gd shows stronger ferromagnetism at 273 K than iron and is widely used in medical science and industry [1]. Dy also shows ferromagnetism at low temperature [1]. Faraday rotator glass [2], which containing Tb, is widely used for manufacture rotator, isolator and circulator in laser technology. Er2O3 is used to make glazes. The properties of rare earth elements and their compounds are interesting and technologically important. Atomistic simulation can help us to understand the underlying mechanism of their attractive performance. Since it is difficult for first-principles methods to deal with large scale systems containing more than thousands atoms [3], semiempirical or empirical methods which can deal with long time and large scale systems have attracted con⇑ Corresponding authors. E-mail addresses: [email protected] (Y. Qin), [email protected] (J. Zhao). http://dx.doi.org/10.1016/j.commatsci.2015.10.014 0927-0256/Ó 2015 Elsevier B.V. All rights reserved.

siderable attentions in the past three decades. The key point is to develop the reliable potential parameters that can describe the atomistic interactions simply and correctly. Pair-wise potentials (e.g. Lennard-Jones potential, Morse potential) own the simplest form but they are unable to describe the many-body effects in metals. Many-body potential, such as tight-binding (TB) potential [4,5], embedded atom method (EAM) potential [6], Finnis–Sinclair (FS) potential [7], Glue potential [8], have been developed and successfully applied to metals and alloys, including hcp metals [9–14]. With a concise form and almost fewest parameters, Gupta potential (or TB potential) [4] has drawn considerable attention among the current many-body empirical potentials. It was first proposed by Gupta [4] to calculate surface relaxation in 1981. Then, Tománek et al. [15] utilized this potential form to consider the electronic contribution and applied the potential to impurity segregation at metal surfaces in 1985. In 1993, Cleri and Rosato [5] developed the TB potentials for transition metals within second-moment approximation and reproduce the fundamental properties of some fcc (Ni, Cu, Pd, Rh, Ag, Ir, Pt, Au, Pb) and hcp (Ti, Zr, Co, Cd, Zn, Mg) metals. Since then, Gupta potential has been widely used in atomistic simulations of different elementary metals and binary alloys [16–29]. So far, several empirical potentials have been developed for rare earth metals. For instance, Hachiya and Ito [29] presented a hybridized nearly-free-electron-tight-binding-bond (NFE-TBB) potential for some rare earth metals and calculated the elastic constants, interatomic spacing and bulk modulus. Baskes and Johnson [9] fitted a modified embedded atom method (MEAM) potential and calculated the defect energies, structural energy and lattice constants for several hcp rare earth metals (Dy, Er, Gd, Tb). Our group has

76

X. Li et al. / Computational Materials Science 112 (2016) 75–79

obtained the parameters of Gupta potential for two fcc rare earth metals: lanthanum and cerium [28]. However, within our knowledge, there is no Gupta potential fitted for hcp rare metals. In this work, we develop Gupta potentials for Er, Dy, Dd, Lu, Tb of the hcp phase, which are able to describe the physical properties (including lattice constants, cohesive energy, elastic modulus, sound velocity, Debye temperature, phonon dispersion, surface energy and vacancy formation energy) of Dy, Er, Gd, Tb and Lu solids reasonably well and would be useful for the future simulations of these rare earth metals and related materials. 2. Theoretical methods and potential parameters

minimize the numerical first derivatives of the sum of square deviations. The experimental lattice constants (a, c), cohesive energy (Ec), elastic constants (Cij) for these five hcp solids [31–33] were employed as inputs and the cutoff distance was chosen as 12 Å. For hexagonal crystals, there are five independent elastic constants (C11, C12, C13, C33 and C44). The bulk modulus (B) and polycrystalline shear modulus (G) were computed using the Hill average scheme [34]

1 ðBR þ Bv Þ 2 1 GH ¼ ðGR þ Gv Þ 2

BH ¼

The Reuss (BR and GR) and Voigt (BV and GV) bounds are given as [35]:

Based on the second-moment approximation of tight-binding theory, the Gupta potential is introduced and the ion–ion interaction is described by an electronic band term and a repulsive term. The general form is given below and the detailed meanings of each term can be found in the original papers [4,5]:

BR ¼

0 1      !12 X U X@X r ij rij A EC ¼ A exp p 1  exp 2q 1 2 i r0 r0 j–i jð–iÞ

and

ð1Þ where U, q, A, p are empirical parameters; r0 is the equilibrium firstneighbor distance in hcp solid; rij represents the distance between atom i and atom j. In this work, the General Utility Lattice Program (GULP) [30] was used to fit the potential parameters and calculate some physics properties. In brief, the potential fitting procedure is aimed to minimize the average deviation between the experimental data and the calculated ones by adjusting four fitting parameters of the Gupta potential. GULP employs the Newton–Raphson algorithm to

ð2Þ

c33 ðc11 þ c12 Þ  2c213 c11 þ c12 þ 2c33  4c13 2ðc11 þ c12 Þ þ 4c13 þ c33 Bv ¼ 9

ð3Þ

5 c44 c66 ðc33 ðc11 þ c12 Þ  2c213 Þ 2 ðc44 þ c66 Þðc33 ðc11 þ c12 Þ  2c213 Þ þ 3Bv c44 c66 12c44 þ 12c66 þ c11 þ c12 þ 2c33  4c13 Gv ¼ 30 GR ¼

ð4Þ

respectively, where c66 ¼ c11  c12 [36]. The longitudinal (v L ) and transversal (v T ) sound velocities can be calculated from the bulk modulus and shear modulus [37],

4 3

qv 2L ¼ B þ G qv 2T ¼ G;

ð5Þ

where q is the density. The Debye temperature is calculated from the following equation: Table 1 Gupta potential parameters for Er, Gd, Tb, Lu and Dy.

Er Gd Tb Lu Dy

HD ¼

U (eV)

A

P

q

r0 (Å)

0.80496 0.89719 1.08718 0.87562 0.82070

0.016840 0.017699 0.038474 0.015867 0.029736

24.63088 23.15789 16.00164 24.94317 21.02317

0.453437 0.333519 0.496138 0.292646 0.516655

3.56 3.63 3.60 3.50 3.59

 1=3 h 3 vm kB 4pV

ð6Þ

where V is the average atomic volume, h is the Planck’s constant, kB is the Boltzmann’s constant, and v m is the average sound velocity determined from the following equation:

v 3 m ¼

1 3



1

þ 3

2



ð7Þ

v L v 3T

Table 2 Lattice constants (a, c), cohesive energy (EC) and elastic constants (Cij) of Er, Gd, Tb, Lu and Dy crystal from experiment [31–33], NFE-TBB model [29] and this work. a (Å)

c (Å)

EC (eV)

C11 (GPa)

C12 (GPa)

C13 (GPa)

Er

Expt. Gupta NFE-TBB

3.56 3.49

5.59 5.68

3.29 3.31

86.34 77.67 82.35

30.50 32.60 27.02

24.38 23.57 25.38

C33 (GPa) 85.54 91.22 89.31

C44 (GPa) 28.09 21.65 23.51

Gd

Expt. Gupta NFE-TBB

3.63 3.57

5.78 5.85

4.14 4.18

66.67 68.57 70.86

24.90 28.01 24.48

21.32 18.38 18.78

71.91 73.55 74.73

20.69 17.11 18.26

Tb

Expt. Gupta NFE-TBB

3.60 3.59

5.70 5.72

4.05 4.02

67.88 58.29 72.51

24.32 27.51 25.47

22.99 22.46 19.26

72.25 80.26 75.62

21.40 19.26 18.63

Lu

Expt. Gupta NFE-TBB

3.50 3.44

5.55 5.61

4.43 4.44

86.23 79.11 86.72

32.03 32.49 27.77

28.0 23.26 25.24

80.86 91.16 100.36

26.79 21.96 30.11

Dy

Expt. Gupta NFE-TBB

3.59 3.59

5.65 5.73

3.04 3.01

74.66 61.52 69.63

26.16 28.01 16.64

21.01 24.17 29.08

78.71 87.78 110.97

24.27 21.39 29.49

77

X. Li et al. / Computational Materials Science 112 (2016) 75–79

3. Results and discussion The potential parameters for hcp Dy, Er, Gd, Lu and Tb fitted by GULP procedure are summarized in Table 1. The essential physical properties of Dy, Er, Gd, Lu and Tb solids calculated from Gupta potentials are listed in Table 2, compared with experimental results [31–33] and previous theoretical work [29]. One can see that the calculated lattice constants (a, c) and cohesive energies (Ec) of the five hcp crystals coincide with the experimental values [31–33]. The lattice constants (a/c) of Tb and Dy from the present potentials (shown in Table 1) are 3.59 Å/5.72 Å and 3.59 Å/5.73 Å, respectively, which are in excellent consistence with the experimental results (3.60 Å/5.70 Å and 3.59 Å/5.65 Å). The derivations of cohesive energies (Ec) of hcp Tb and Dy between our potentials and experiment are both 0.03 eV. The calculated lattice constants (a/c) of Er, Gd and Lu hcp solids have a derivation of 2% compared to the measured ones and the cohesive energies (Ec) of Er, Gd and Lu almost equal to the experimental values. From the above comparisons, it is reliable that our fitted potentials can describe the lattice constants and cohesive energies very well. From Table 2, one can see that the calculated elastic constants from current Gupta potentials are consistent with the experiment results. The theoretical C11, C12, C13, C33 and C44 of Gd crystal are 68.57 GPa, 28.01 GPa, 18.38 GPa, 73.55 GPa, and 17.11 GPa, respectively, showing better agreement with the experimental ones [31– 33] than the results calculated from NFE-TBB model potential [29]. As for Er, our calculated data underestimate C11, C13 and C44 by 8.67 GPa, 0.81 GPa and 6.44 GPa, respectively, but overestimate C12 and C33 by 2.10 GPa and 5.68 GPa, respectively. In the case of Lu, the largest derivation among the five elastic constants (C11, C12, C13, C33 and C44) is 10.3 GPa belonging to C33. As for C11 of all five hcp crystals, the maximum discrepancy between experimental results and previous predictions is 13.14 GPa. The minimum and maximum deficiency of all five C12 is 0.46 GPa and 3.19 GPa,

respectively. The theoretical results for C13, C33 and C44 of all five hcp crystals are reasonably consistent with the experimental ones. The calculated elastic constants of hcp Dy and Tb solids show a little bit larger derivation compared to the measured ones. These derivations may be mainly attributed to the simple form of Gupta potential with only four adjustable parameters and the fact that Gupta potential is based on the central interactions [5]. To further assess the reliability of the current Gupta potential, the bulk modulus, shear modulus, sound speed, and Debye temperature of the five hcp crystals were calculated and compared with the values from experiments [31–33] and NFE-TBB model [29]. From Table 3, we can see that the theoretical bulk modulus for Er, Gd, Tb, Lu and Dy solids is 45.10 GPa, 37.79 GPa, 37.77 GPa, 45.26 GPa and 40.06 GPa, respectively, which shows excellent agreement with the experimental values with the maximum derivation of 4.9% for Lu. The computed results of the shear modulus for all five hcp crystals show small underestimations when compared to the measured ones [31–33] and the maximum derivation is 6 GPa for Er. But the shear modulus of Gd calculated from Eq. (2) is 26.86 GPa, in good agreement with the experimental value of 28.58 GPa. Using Eq. (5), the transverse and longitudinal sound velocities were computed from bulk modulus and shear modulus. The theoretical sound velocities are compared with the ones calculated using the experimental elastic constants [31–33] and mass density since there are no direct experimental data. In the case of Er, both transverse and longitudinal velocities are underestimated by about 9.5% and 5.7%, respectively. As for Gd, Tb, Lu and Dy, the sound velocities are also underestimated with regard to the experimental ones, with the maximum deficiency of 8.9% of transverse velocity for Dy. Most of the calculated Debye temperatures are smaller than the experimental ones with the maximum derivation of 8.27% for Er. The only exception case is Dy, whose Debye temperature is overestimated by about 11.76%.

Table 3 Physical properties of hcp solids of Er, Gd, Tb, Lu and Dy from experiments [31–33,43], NFE-TBB model [29] and this work.

a b c d

Gupta

NFE-TBB

Experiment

Er

Bulk modulus B (GPa) Shear modulus G (GPa)a Transverse sound velocity vT (m s1)b Longitudinal sound velocity vL (m s1)b Debye temperature H (K)c

45.10 31.75 1849.2 3068.7 247.2

45.47 35.30 – – –

45.47 37.75 2042.0 3253.0 269.5

Gd

Bulk modulus B (GPa) Shear modulus G (GPa)a Transverse sound velocity vT (m s1)b Longitudinal sound velocity vL (m s1)b Debye temperature H (K)c

37.79 26.86 1821.6 3015.3 237.4

37.83 29.09 – – –

37.83 28.58 1900.1 3097.0 245.4, 200d

Tb

Bulk modulus B (GPa) Shear modulus G (GPa)a Transverse sound velocity vT (m s1)b Longitudinal sound velocity vL (m s1)b Debye temperature H (K)c

37.77 24.87 1737.5 2934.3 227.4

38.72 29.58 – – –

38.78 29.33 1885.0 3071.6 246.0

Lu

Bulk modulus B (GPa) Shear modulus G (GPa)a Transverse sound velocity vT (m s1)b Longitudinal sound velocity vL (m s1)b Debye temperature H (K)c

45.26 32.41 1793.4 2963.2 242.6

47.59 40.52 – – –

47.59 36.04 1910.8 3112.7 256.3, 210d

Dy

Bulk modulus B (GPa) Shear modulus G (GPa)a Transverse sound velocity vT (m s1)b Longitudinal sound velocity vL (m s1)b Debye temperature H (K)c

40.06 27.17 1791.3 3001.5 234.7

41.07 37.60 – – –

41.07 33.05 1965.3 3154.2 210d

Values calculated from the elastic constants using the Eq. (2). Sound velocity evaluated from elastic constants and mass density using Eq. (5). Debye temperature calculated from Eq. (6). From Ref. [43].

78

X. Li et al. / Computational Materials Science 112 (2016) 75–79

The phonon dispersion curves, as a crucial examination for the performance of a potential, were calculated using the present sets of parameters. Here the path of phonon dispersion was set along  0, ½0 1 1  0 and ½0 0 0 1 directions in the Brillouin zone ½1 0 0 0, ½1 12 of the hexagonal lattice. The calculated phonon dispersion curves for all five hcp metals, together with available experimental data for Tb [38] are shown in Fig. 1. One can see that the calculated trends of phonon dispersion for Tb basically agree with the measured ones, which indicates that the present potential can describe the lattice dynamics reasonably. Since there are no available measurements for hcp Dy, Er, Lu and Gd, the present theoretical results can provide vital guidance for future experiments. In addition to the perfect hcp crystals of five rare earth elements, we also examined some properties of imperfect crystals, including surface energy and vacancy formation energy. Since there are no vacancy formation energies for the five hcp elements in experiments, we compare the monovacancy formation energies calculated from present Gupta potentials, MEAM potential [9] and EAM potential [38]. The monovacancy formation energies calculated for Er and Dy from present potentials are 1.55 eV and 1.50 eV, respectively, which are close to the results from Baskes’s MEAM (i.e., 1.45 eV and 1.34 eV) [9], but larger than the values calculated by Hu’s EAM potential [38]. The calculated monovacancy formation energies for hcp solids of Gd and Tb (1.96 eV and 2.96 eV, respectively) are larger than those calculated from Hu’s and Baskes’s. We also calculated the vacancy formation energy for Lu (2.16 eV), which is 48.1% of cohesive energy of the perfect crystal. Divacancy formation energies of in-plane and out-ofplane modes are also calculated and compared with MEAM potential [9] in Table 4. The formation energies of in-plane divacancy for Er, Gd, Tb, Lu and Dy are 3.05 eV, 3.35 eV, 4.79 eV, 4.12 eV, and 2.90 eV, respectively, and those of out-of-plane ones for Er, Gd, Tb, Lu and Dy are 3.04 eV, 3.76 eV, 4.82 eV, 4.18 eV, 2.90 eV, respectively. It is known that vacancies reflect the kinetic and thermodynamic properties of metals and intermetallic compounds [39]. Gupta potential with only four adjustable parameters can only roughly reproduce this effect but cannot describe it very accurately. In Table 4, the calculated (0 0 0 1) surface energy using Gupta potential was also compared to the results from previous MEAM potential [9], Miedema’s theory [40], full charge density (FCD) method [41] and EAM potential [38]. Previously, Tyson [42] estimated the surface energies for the five rare earth solids and his results are also shown in Table 4. One can see that the calculated (0 0 0 1) surface energies for Er and Dy are smaller than those calculated from MEAM potential [9], but larger than those from Miedema’s theory [40] and EAM potential [38]. The theoretical (0 0 0 1) surface energies for Gd and Tb from Gupta potentials are 2.14 J/ m2 and 1.90 J/m2, respectively, which are larger than the three sets of other theoretical results. Due to different atomic arrangement of  0Þ surfaces, the two ð1 0 1  0Þ surface energies second layer of ð1 0 1 are calculated and compared with MEAM potential [9] and EAM potential [38]. From Table 4, we can see that the difference of  0Þ surfaces is small. The the energies between the two ð1 0 1  0Þ surface energies for Er and Dy are 1.80/1.79 J/m2 and ð1 0 1 1.65/1.69 J/m2, respectively, which lie between Hu’s results (0.58 J/m2 and 0.52 J/m2, respectively) and Baskes’s values (2.25 J/  0Þ plane m2 and 2.24 J/m2, respectively). Surface energies of ð1 0 1 for Gd and Tb are 2.23/2.22 J/m2 and 2.19/2.15 J/m2, which are the largest values among the three sets of theoretical results. The surface energies for hcp Lu solid are 2.71 J/m2 of (0 0 0 1) surface,  0Þ surface, respectively, which are larger 2.70/2.71 J/m2 of ð1 0 1 than the results from full charge density (FCD) method [41]. Due to the lack of experimental values, we cannot draw any definitive

Fig. 1. Phonon dispersions of Dy, Er, Gd, Lu and Tb solids with hcp structure. For Tb, the experimental data [38] are given in open circles.

conclusion about the performance of our potentials on the surface energy. Nevertheless, according to the above discussions and comparisons about vacancy formation energies and surface energies, the present Gupta potentials are still in an acceptable manner for the description of defective properties of Er, Gd, Tb, Lu and Dy solids.

79

X. Li et al. / Computational Materials Science 112 (2016) 75–79

Table 4 Surface energy and vacancy formation energy of hcp solids of Er, Gd, Tb, Lu and Dy from present Gupta potential, MEAM potential [9], AMEAM potential [38], Miedema’s theory [40], full charge density (FCD) method [41], and Tyson’s results [42].

Surface energy (0 0 0 1) (J/m2)

This work Previous work

 0ÞA (J/m2) Surface energy ð1 0 1  0ÞB (J/m2) Surface energy ð1 0 1

This work This work Previous work

Monovacancy formation energy (eV)

This work Previous work

Divacancy formation energy (eV)

In-plane Previous workb Out-of-plane Previous workb

a b c d e

From From From From From

Ref. Ref. Ref. Ref. Ref.

Gd

Tb

Dy

Lu

1.82 0.58a 2.32b 1.17c 0.84e 1.80 1.79 0.58a 2.25b 1.55 1.31a 1.45b 3.05 2.46 3.04 2.45

2.14 0.49a 1.33b 1.11c 1.00e 2.23 2.22 0.48a 1.29b 1.96 1.14a 1.25b 3.35 2.12 3.76 2.12

1.90 0.52a 1.52b 1.13c 1.00e 2.19 2.15 0.51a 1.49b 2.96 1.18a 1.29b 4.79 2.21 4.82 2.20

1.39 0.52a 2.33b 1.14c 0.77e 1.65 1.69 0.52a 2.24b 1.50 1.22a 1.34b 2.90 2.27 2.90 2.26

2.71 1.60d 1.23c 1.15e 2.71 2.70 1.42d 1.62d 2.16 – – 4.12 – 4.18 –

[38]. [9]. [40]. [41]. [42].

4. Conclusion A series of empirical Gupta-type potentials have been developed for Er, Dy, Gd, Tb and Lu with hcp solids. The experimental cohesive energies, lattice constants and elastic constants are reproduced. The bulk modulus, shear modulus, sound velocities, Debye temperature for each hcp crystal are calculated using the present potential parameters, and show satisfactory agreement with experimental results. Since there are no experimental data available, vacancy formation energies and surface energies are calculated and compared to previous theoretical values. The phonon dispersion curves for hcp Er, Dy, Gd, Tb and Lu metals are also calculated, and the theoretical trend of different branches for Tb is in line with the measured data. Unfortunately, it is difficult to describe all properties of rare earth elements with f electrons using such a simple empirical potential. Even with some deficiencies, the fitted Gupta potential, with only four adjustable parameters, can describe most fundamental properties of Er, Dy, Gd, Tb and Lu solids. Hence, we expect that the present version of Gupta potentials would be useful for the future atomic simulations of Er-, Dy-, Gd-, Tb- and Lu-based materials. Acknowledgement This work was supported by the National Natural Science Foundation of China (11134005, 11574040, 51471043). References [1] [2] [3] [4] [5] [6] [7] [8]

Er

R. Elliott, Magnetic Properties of Rare Earth Metals, 2013. Z. Peiming, W. Biao, M. Hanfen, J. Yasi, Chinese J. Lasers 23 (1996) 857–860. A. Carlsson, Solid State Phys. 43 (1990) 1–91. R.P. Gupta, Phys. Rev. B 23 (1981) 6265–6270. F. Cleri, V. Rosato, Phys. Rev. B 48 (1993) 22–33. M.S. Daw, M.I. Baskes, Phys. Rev. Lett. 50 (1983) 1285–1288. M. Finnis, J. Sinclair, Philos. Mag. A 50 (1984) 45–55. F. Ercolessi, E. Tosatti, M. Parrinello, Phys. Rev. Lett. 57 (1986) 719–722.

[9] M.I. Baskes, R.A. Johnson, Modelling Simulat. Mater. Sci. Eng. 2 (1994) 147– 163. [10] M. Igarashi, M. Khantha, V. Vitek, Philos. Mag. B 63 (1991) 603–627. [11] R. Johnson, D. Oh, J. Mater. Res. 4 (1989) 1195–1201. [12] R. Johnson, Philos. Mag. A 63 (1991) 865–872. [13] G.J. Ackland, Philos. Mag. A 66 (1992) 917–932. [14] G. Ackland, S. Wooding, D. Bacon, Philos. Mag. A 71 (1995) 553–565. [15] D. Tománek, A.A. Aligia, C.A. Balseiro, Phys. Rev. B 32 (1985) 5051–5056. [16] J. Guevara, A.M. Llois, M. Weissmann, Phys. Rev. B 52 (1995) 11509–11516. [17] G.C. Kallinteris, N.I. Papanicolaou, G.A. Evangelakis, D.A. Papaconstantopoulos, Phys. Rev. B 55 (1997) 2150–2156. [18] M.M.G. Alemany, O. Diéguez, C. Rey, L.J. Gallego, Phys. Rev. B 60 (1999) 9208– 9211. [19] L. Gómez, A. Dobry, H.T. Diep, Phys. Rev. B 55 (1997) 6265–6271. [20] F. Willaime, C. Massobrio, Phys. Rev. Lett. 63 (1989) 2244–2247. [21] W.S. Lai, B.X. Liu, Phys. Rev. B 58 (1998) 6063–6073. [22] C. Massobrio, V. Pontikis, G. Martin, Phys. Rev. B 41 (1990) 10486–10497. [23] F. Willaime, C. Massobrio, Phys. Rev. B 43 (1991) 11653–11665. [24] K. Michaelian, N. Rendón, I.L. Garzón, Phys. Rev. B 60 (1999) 2000–2010. [25] C. Rey, L.J. Gallego, J. García-Rodeja, J.A. Alonso, M.P. Iñiguez, Phys. Rev. B 48 (1993) 8253–8262. [26] X. Li, Y. Qin, J. Fu, J. Zhao, Comput. Mater. Sci. 98 (2015) 328–332. [27] J. Fu, J. Zhao, Modell. Simul. Mater. Sci. Eng. 21 (2013) 065003. [28] J. Fu, C. Zhang, J. Zhao, Comput. Mater. Sci. 85 (2014) 142–146. [29] K. Hachiya, Y. Ito, J. Phys.: Condens. Matter 11 (1999) 6543. [30] J.D. Gale, A.L. Rohl, Mol. Simul. 29 (2003) 291–341. [31] K. Salama, F. Brotzen, P. Donoho, J. Appl. Phys. 43 (1972) 3254–3258. [32] J. Tonnies, K. Gschneidner Jr., F. Spedding, J. Appl. Phys. 42 (1971) 3275–3283. [33] G. Simmons, H. Wang, Single Crystal Elastic Constants and Calculated Aggregate Properties, 1971. [34] R. Hill, Proc. Phys. Soc. A 65 (1952) 349. [35] L. Vitos, Computational Quantum Mechanics for Materials Engineers: The Emto Method and Applications, 2007. [36] R. Pasianot, E. Savino, Phys. Status Solidi (b) 176 (1993) 327–334. [37] H. Ledbetter, J. Appl. Phys. 44 (1973) 1451–1454. [38] W. Hu, H. Deng, X. Yuan, M. Fukumoto, Eur. Phys. J. B – Condens. Matter Complex Syst. 34 (2003) 429–440. [39] P.A. Korzhavyi, I.A. Abrikosov, B. Johansson, A. Ruban, H.L. Skriver, Phys. Rev. B 59 (1999) 11693. [40] F.R. De Boer, W. Mattens, R. Boom, A. Miedema, A. Niessen, Cohesion in Metals: Transition Metal Alloys, ed., Amsterdam, 1988. [41] L. Vitos, A. Ruban, H.L. Skriver, J. Kollar, Surf. Sci. 411 (1998) 186–202. [42] W.R. Tyson, Can. Metall. Quart. 14 (4) (1975) 307–314. [43] K. Charles, M. Paul, M. Paul, Introduction to Solid State Physics, 1976.