- Email: [email protected]
First-principles study of solute–solute binding in magnesium alloys
Computational Materials Science 103 (2015) 97–104
Contents lists available at ScienceDirect
Computational Materials Science journal homepage: www.elsevier.com/locate/commatsci
First-principles study of solute–solute binding in magnesium alloys Guobao Liu a, Jing Zhang a,b,⇑, Yuchen Dou a a b
College of Materials Science and Engineering, Chongqing University, Chongqing 400044, China National Engineering Research Center for Magnesium Alloys, Chongqing 400044, China
a r t i c l e
i n f o
Article history: Received 23 September 2014 Received in revised form 9 March 2015 Accepted 17 March 2015 Available online 3 April 2015 Keywords: Magnesium alloys Solute–solute binding First-principles calculations
a b s t r a c t Solute–solute interactions play a major role in the properties of materials. In this work, we present an extensive database of solute–solute binding energies that captures the detailed interactions in Mg-based alloys from first-principles calculations based on density functional theory. The effects of solute–solute binding energies on magnesium properties, precipitation hardening responses and stacking fault energies in particular, are inferred and discussed. The results of our calculations regarding bindings between solutes with different chemistries, including Al–Sn, Al–Ca, Ca–Zn, Ca–In, and Sn–Zn, were validated using available experimental investigations. Solute pairs that were predicted to show large positive (e.g., Yb–Bi/ Sn/Pb and Ca–Bi/Sn/Pb) and negative (e.g., Bi–Sn/Pb/Al) values of binding energies exhibited potential in modifying the precipitation sequence and stacking fault energy. Moreover, alloys added with these alloying elements may exhibit unique mechanical properties, which await experimental verification. Finally, the effect of physical features, including atomic radius and electronegativity, on the solute–solute bindings was investigated. Ó 2015 Elsevier B.V. All rights reserved.
1. Introduction Magnesium alloys in wrought forms (e.g., extrusions, forgings, sheet, and plates) have gained considerable worldwide interest in the past decade as potential replacements for heavier steel and aluminum alloys because of their high strength–density ratio [1–3]. However, the vast majority of Mg alloy applications are presently covered by cast products [1,4]. The application of wrought Mg alloys is hindered in part by their poor room temperature formability caused by their strong basal type texture and large anisotropy between basal and non-basal (prismatic and pyramidal) slips. The critical resolved shear stress of a basal slip system in Mg at room temperature is approximately two orders of magnitude lower than those of non-basal slip systems on prismatic and pyramidal planes and somewhat less than that of twinning. Therefore, basal slip and twinning are almost the only deformation mechanisms in polycrystalline Mg alloys at room temperature, which is far from the five independent operating systems required by the von Mises’ criterion for sufficient ductility [5–8]. Meanwhile, the strength for most commercial Mg alloys is still substantially lower than their aluminum alloy counterparts, which is in part attributed
⇑ Corresponding author at: College of Materials Science and Engineering, Chongqing University, Chongqing 400044, China. Tel.: +86 23 65111167; fax: +86 23 65102821. E-mail address: [email protected] (J. Zhang). http://dx.doi.org/10.1016/j.commatsci.2015.03.023 0927-0256/Ó 2015 Elsevier B.V. All rights reserved.
to the low precipitation hardening responses of these alloys [3,9,10]. Therefore, developing new Mg alloys with higher strength and ductility is necessary for Mg alloys to remain competitive with aluminum alloys. In principle, the addition of alloying elements capable of modifying the precipitation sequence and the strengths of various deformation modes, thereby reducing crystallographic anisotropy, holds promise as a means to improve the strength and ductility of Mg alloys. Over the past decades, numerous experimental studies on the influence of alloying additions on precipitation hardening responses of Mg alloys have been conducted [9–17]. The agehardenability of Mg alloys improves significantly when the alloying elements for addition are selected appropriately. The interaction between solutes distinctly affects precipitation hardening responses by forming new types of precipitates or manipulating the density, distributions, sizes, morphologies, and growth orientation of the precipitates. For instance, the ultimate tensile stress of peak-aged Mg–6Gd–2Zn–0.6Zr (wt.%) alloys is nearly two times higher than its Mg–6Gd–0.6Zr (wt.%) counterpart because of the formation of the c0 (Mg70Gd15Zn15, ordered hcp) phase [11]. Small additions of Zn in Mg-rare earth (RE) [11], Mg–Ca [12], and Mg–Sn [13] alloys have been widely reported to promote the formation of fine-scale precipitates and increase the precipitate number density during ageing heat treatment, thereby increasing the peak-aged hardness of alloys. Mendis et al. [14] observed that the number density of Mg2Sn precipitates increased by
98
G. Liu et al. / Computational Materials Science 103 (2015) 97–104
approximately one order of magnitude in Mg–Sn alloys with the collective addition of Li + In, resulting in a hardening increment increase of 150%. Clustering or co-clustering of solute atoms has been proposed to possibly serve as heterogeneities for precipitate nucleation, leading to increased number density of the precipitates. The increased number density of thin prismatic plates in Mg–Ca alloys with added In has also been reported to enhance age-hardening response by three times [15]. It was suggested that the large strain field expected from Ca–In co-clusters changes the habit planes from basal to prismatic. Alloying elements in the form of substitutional solutes also affect the strength and ductility of Mg alloys by modifying certain underlying factors, such as stacking fault energy and solute/dislocation interaction energy [5–7,18]. Stacking fault energy is directly associated with the dissociation behavior of dislocations, further affecting the strength, ductility, and fracture of materials. Nogaret et al. [19] correlated the dislocation structures and Peierls stresses to gamma surfaces (generalized stacking fault surfaces) in pure Mg. The simulation result, in which pyramidal I {1 0 1 1}h1 1 2 3i slip efficiently occurs in comparison with the pyramidal II {1 1 2 2}h1 1 2 3i slip, correlates well with the calculated gamma surfaces. Theoretical calculations of stacking fault energies of several binary Mg–X alloys have been reported in literature [20–23]. The effects of solutes on the basal [5,18], prismatic [6,7], and twinning dislocations [24] by inputting the firstprinciples data (e.g., size misfit, chemical misfit derived from stacking fault energy, and solute/dislocation interaction energy) into strengthening models have been quantitatively calculated in recent years. At present, investigating the multi-solute effects for Mg alloys is severely limited, although these effects have been recently reported for Al alloys [25]. However, the interaction between solute atoms plays an important role in stacking fault energy. For example, the unstable stacking fault energy of basal slip decreased more significantly for systems with simultaneous addition of Y and Zn than those with single Zn or Y additions [26]. The increased separation distance between two Shockley partial dislocations in Mg–Zn–Y alloys observed using transmission electron microscopy also demonstrates that the combined addition of Zn and Y leads to a significant reduction in stacking fault energy. The interaction between solutes plays an important role in determining the underlying factors, such as precipitation process and stacking fault energy, further affecting the mechanical properties of Mg alloys. Binding energy, in terms of solute–solute/vacancy binding, captures the bonding properties of atoms and provides detailed information on the interaction between solute atoms. For the precipitation sequence, solutes with favorable bindings are likely to attract and form new types of precipitates, whereas those with unfavorable bindings tend to repel each other, thereby decreasing the solubility of solutes in the matrix and promoting precipitate density during aging treatments. For stacking fault energy, if two solutes with favorable binding are placed in two different basal planes, extra energy is definitely needed to shear crystals on this plane. Moreover, solute–solute binding energies can serve as input first-principles data to quantitatively predict the multi-solute effects for Mg alloys if a constitutive model is developed. Although solute-vacancy binding energies in Mg have been systematically studied by Shin and Wolverton [27] and Saal and Wolverton [28], calculation for solute–solute binding energies, at least to the best of our knowledge, only touches upon Al, Zn, Y, and Gd [29]. Therefore, it is necessary to take an extensive calculation of the binding energies of solute atoms in Mg to investigate the interaction between solute atoms. In this paper, first-principles density functional calculations were carried out to calculate the binding energies of various solute pairs to construct an extensive first-principle database for predicting solute–solute binding energies in Mg alloy. Alloying elements
with an equilibrium solid solubility larger than 0.5 (at.%) in Mg were selected: Ag, Al, Bi, Ca, Dy, Er, Ga, Gd, Ho, Li, Lu, In, Pb, Nd, Sc, Sm, Sn, Tm, Y, Yb, Zn, and Zr. The solute–solute binding energy trends were then investigated in terms of physical features, such as atomic radius difference and electronegativity difference of solute pairs. The charge density difference in selected planes of Mg–Bi–Yb and Mg–Bi–Sn alloys were computed to reveal the underlying mechanism of solute effects on the binding energies of solute atoms. 2. Computational details In this study, first-principles total-energy DFT calculations were performed using the Vienna Ab-initio Simulation Package. The Perdew–Burke–Ernzerhof (PBE) exchange–correlation functional for the generalized-gradient-approximation (GGA) and projectoraugmented wave method was used in the present work [30–33]. All calculations were fully relaxed with respect to all degrees of freedom. The first-order Methfessel–Paxton with SIGMA = 0.2 eV was used for geometric relaxation until the electronic energy converged to less than 105 eV/cell and the Hellmann–Feynman force on all atomic sites was less than 0.01 eV Å1. Two successive structural optimizations (adapting basis vectors and computational grids to cell parameters) were conducted to ensure that the cell energies and structural parameters fully converged [34]. Total energy calculation was then performed using linear tetrahedron method with Blöchl correction [35,36]. A plane wave cutoff energy, which is 1.3 times the largest cutoff energy associated with the elements of interest recommended by VASP, was automatically performed by setting the PREC = HIGH in INCAR. It is known that the f electrons are not handled well by presently available density functions because of their varying tendency to form localized states. A common solution to this problem is to place the f electrons in the core (so-called ‘‘frozen’’ potentials). The spatial localization of f electrons near the core is a good approximation and works well in Al–RE [37], Mg–RE [38], Mg–Pb, and Mg–Bi [7] alloys. Therefore, frozen potentials are used for energy calculation in systems embedded with f elements (RE, Pb, and Bi in the present work). Table 1 lists the PAW potentials and the corresponding cutoff energies for Mg and 22 solutes Table 1 The PAW valence configuration, energy cutoff for the total energy calculation of pure Mg, Mg–X, and Mg–X–Z, and the atomic radius and electronegativity of Mg and the alloying elements considered in the present work. Alloying element
PAWPP
Cutoff (eV)
Atomic radius (pm)
Electronegativity
Mg Ag Al Bi Ca Dy Er Ga Gd Ho In Li Lu Nd Pb Sc Sm Sn Tm Y Yb Zn Zr
[Ne]3s2 [Kr]4d105s1 [Ne]3s23p1 ([Xe]4f145d10)6s26p3 [Ar]4s2 ([Xe]4f10)6s2 ([Xe]4f12)6s2 ([Ar]3d10)4s24p1 ([Xe]4f7)5d16s2 ([Xe]4f11)6s2 ([Kr]4d10)5s25p1 [He]2s1 ([Xe]4f14)5d16s2 ([Xe]4f4)6s2 ([Xe]4f145d10)6s26p2 [Ar]3d24s1 ([Xe]4f6)6s2 ([Kr]4d10)5s25p2 ([Xe]4f13)6s2 [Kr]4d25s1 ([Xe]4f14)6s2 [Ar] 3d104s2 [Kr]4d35s1
164 325 312 315 155 202 201 175 201 200 125 182 201 237 309 289 230 313 193 263 146 360 299
160 144 143.1 154.7 197 178.1 176.1 135 180.4 176.2 167 173.8 152 181.4 175 162 180.4 151 175.9 180 193.3 134 160
1.31 1.93 1.61 2.02 1 1.22 1.24 1.81 1.2 1.23 1.78 0.98 1.27 1.14 1.87 1.36 1.17 1.96 1.25 1.22 1.1 1.65 1.33
G. Liu et al. / Computational Materials Science 103 (2015) 97–104
considered in this work. Given that the unpaired 4f electrons are primarily responsible for magnetism, magnetic contribution to the total energy is negligible even for elements known to be ferromagnetic at room temperature, such as Gd, when the frozen potentials are used. As a result, non-spin polarized calculation was considered in all calculations in this work. Supercells for solute–solute binding energy calculations were constructed by duplicating the orthorhombic unit cell with four atoms along three lattice vectors (x is aligned with [11–20], y with [1 0 1 0], and z with [0 0 0 1]) as illustrated in Fig. 1. Table 2 summarizes the supercell dimensions and k-point grids used for solute–solute binding energy calculations in the present work. The binding energy of the solute–solute pairs are defined as follows:
EbðXXÞ ¼ EðMgN2 X2 Þ þ EðMgN Þ 2EðMgN1 X1 Þ EbðXZÞ ¼ EðMgN2 X2 Þ þ EðMgN Þ EðMgN1 X1 Þ EðMgN1 Z1 Þ where Eb is the binding energy, E is the total energy of pure Mg or systems embedded with solute atoms, X and Z represent two different solute atoms in the Mg supercell, and N is the number of lattice sites considered. A positive value of Eb corresponds to favorable binding under this definition. Given that alloying solute atoms affect binding energies by disturbing the local electronic environment and strain field, configurations in which two solute atoms are located respectively at the first- and second-nearest neighbor sites (1NN and 2NN) in the hexagonal lattice were considered in the current calculation. In hcp Mg, two nonequivalent 1NN solute pairs existed and distinguished as follows: one in the same basal (0 0 0 1) plane (1NN) and the other in a {1 1 2 2} pyramid plane (1NN0 ), as shown in Fig. 1a. 3. Results 3.1. Effects of supercell size The effects of supercell size on the binding energies of solute– solute pairs were first investigated. The binding energies of solute–solute pairs commonly reported to show remarkable hardenability in experiments, namely, Ag–Ca [10], Ca–In [15], Ca–Sn [39], Ca–Zn [12], and Li–In [14] for basal plane solute–solute pairs (1NN), were calculated in supercells with 48, 64, 96, and 144 atoms. The calculated solute–solute binding energies for different supercell sizes are listed in Table 3. The table shows that 64-atom calculations converged within 0.02 eV and are thus sufficient to provide reliable binding energies of solutes. As a result, the binding energies of all other solute pairs were calculated in 64-atom supercells, as shown in Fig. 1b. 3.2. Solute–solute binding energy Based on our calculated data, we found that the solute–solute binding energies of Mg were not sensitive to the difference between basal plane (1NN) and out-of-basal-plane (1NN0 ) nearest neighbors. The binding energies of 1NN0 solute pairs were therefore not considered in this study. The calculated binding energies for 1NN and 2NN solute pairs are summarized in Table 4. Basic physical parameters, including atom radii and eletronegativity values of various elements considered in the present work, are listed in Table 1[40]. Table 4 shows that 57 kinds of 1NN solute pairs exhibit positive values of binding energies, indicating that the atoms of these solute pairs are attractive and may get together to form pairs or clusters. Interestingly, most of these solute pairs are composed of large atoms (i.e., RE elements and Ca) and relatively small atoms
99
(Bi, Pb, Sn, In, and Ga) located in the IIIA to VA groups in the periodic table of chemical elements. The solute pair Bi–Yb possesses the largest binding energy in all 1NN solute pairs, with a value of 2.4 eV. It is worth to note that solute pairs such as Ag–Ca, Ca–In, Ca–Zn, and Li–In, which are commonly reported in literature to present unique age hardenability in Mg alloys, were predicted to exhibit favorable 1NN bindings. However, Mg–Sn–Zn alloy systems [17], which were also reported to exhibit impressive precipitation hardening responses, showed an extremely large negative binding energy value of –0.27 eV for Sn–Zn solute pairs. For 2NN solute pairs, the binding energies of 34 kinds of solute pairs are positive. Overall, the magnitude of 2NN binding energies was considerably smaller than that of 1NN binding energies. Notably, most of these pairs are composed of atoms with larger atomic radius than Mg. More specifically, combinations of RE and Ca/Yb provide such examples; among which, the Ca–Nd solute pair is the most unique, with the largest binding energy of 0.08 eV. It is interesting to note that Ca and Yb show unique properties among all the alloying elements considered in the present work. Both X–Ca and X–Yb 1NN solute pairs show relatively larger binding energies compared with other solute pairs, whereas solute pairs showing positive binding energies are all composed of Ca– RE or Yb–RE for 2NN. 3.3. Validating solute–solute bindings Direct experimental measurements for solute–solute binding have yet to be reported. However, the precipitation sequence and theoretical investigations for the stacking fault energies of certain Mg alloys are available to validate our calculation results. Mendis et al. [41] observed the precipitation sequence of the Mg–2.4Zn– 0.1Ag–0.1Ca–0.16Zr alloy. Their 3DAP results showed that Ca forms clusters with Zn in the early stage of precipitation, which agrees well with our calculation that Ca and Zn tend to form favorable bindings. Solute clusters, such as Ca–Al [42], Ca–In [15], and Li–In [14], were also experimentally observed after ageing heat treatments. Notably, solute clusters reported in the above works show considerable potential in enhancing precipitation hardening responses, either by forming new types of precipitates or serving as heterogeneities for the nucleation of precipitates, thereby increasing the density of precipitates. Sasaki et al. [17] and Mendis et al. [13] reported that the age-hardening responses of Mg–Sn alloys were substantially improved by small additions of Zn, which is mainly attributed to the increased number density of precipitates. However, no new types of precipitates were observed except for Mg2Sn and MgZn2. The as-calculated negative binding energy of Zn–Sn solute pair in this work indicates that Zn and Sn tend to repel each other, thereby decreasing the solubility of both Zn and Sn in the Mg matrix and promoting the nucleation of precipitates. Yuasa et al. [43] calculated the stacking fault energies of basal and prismatic slip systems for Mg–Zn–Ca alloys, in which the 1NN0 solute pairs were formed. Plastic anisotropy was reduced because of the non-linear nature (the addition of Zn in Mg–Ca alloy increased the unstable stacking fault energy of basal slip while decreasing that of prismatic) of the unstable stacking fault energy of basal slip. The increase in unstable stacking fault energy of basal slip agrees with our calculation, which shows that Ca and Zn tend to form favorable binding, and extra energy is needed to shear the basal plane between these two solutes. 4. Discussions 4.1. Physical contributions to solute–solute binding energies In this section, we attempt to find general trends and physical factors controlling solute binding energies across solute pairs with
100
G. Liu et al. / Computational Materials Science 103 (2015) 97–104
Fig. 1. (a) Schematic of the orthorhombic Mg unit cell (solid red lines) used to generate supercells in the present work. The solute atom, X (the purple sphere), is at the origin of the unit cell, and the solute Z is placed either at the first-nearest neighbor position in basal plane (1NN, the cyan sphere), in pyramid plane (1NN0 , the dark blue sphere), or at the second-nearest neighbor position (2NN, the green sphere). Orange spheres denote Mg atoms. (b) Super-cell with 64 atoms used to calculate binding energies of all solute pairs in the present work. The purple sphere denotes solute X, and the cyan sphere denotes solute Z. Orange spheres are Mg atoms. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
Table 2 Atom number (N) and corresponding solute concentrations, supercell dimensions, and k-point grids used in the calculation. Atom number (N)
Solute concentrations at.%
Supercell dimensions
k-point grids
48 64 96 144
4.17 3.13 2.08 1.39
322 422 423 433
11 11 11 999 777 555
favorable bindings that have been considered in the current work. The interaction between atoms is closely associated with the changing electronic environment and strain field caused by solutes. Yasi et al. [5] computed the so-called solute ‘‘misfits’’, which consist of the following two important components: a size misfit for the change in local volume and a chemical misfit for the change in energy to slip the crystal, to approximate the solute/dislocation interaction energy. For solute–solute interactions, we supposed that the atom radius and electronegativity of solute atoms are key factors that influence solute–solute binding energy. 4.1.1. Atomic radius The image in Fig. 2 shows the calculated 1NN binding energies with positive values as a function of atomic radius difference between solutes. Considering their extremely cluttered nature for clarity, the data were separated into Al–X, Bi–X, Ca–X, Ga–X, In– X, Sn–X, and Yb–X, and the corresponding scatter plots are shown in Fig. 2a–g, respectively. We found that most solute pairs with positive binding energies are composed of a large atom (RE/Ca) and a small atom belonging to IIIA to VA. The solute/dislocation interaction energy via first-principles calculations by Yasi et al. [5] also predicted that the strongest interaction for an Al solute Table 3 Binding energy convergence corresponds to the supercell size. The binding energies for the 1NN solute pairs of Ag–Ca, Ca–In, Ca–Sn, Ca–Zn, and Li–In are computed in supercell with 48, 64, 96, and 144 atoms. Note that positive energies indicate energetically favorable binding. Solute pairs
Ag–Ca Ca–In Ca–Sn Ca–Zn Li–In
Solute–solute binding energy (eV) 48 atoms
64 atoms
96 atoms
144 atoms
0.02 0.09 0.13 0.03 0.04
0.02 0.09 0.13 0.05 0.03
0.03 0.09 0.15 0.06 0.04
0.02 0.10 0.13 0.05 0.03
is at the site of maximum compression near the center of the dislocation partial cores for both edge and screw dislocations. This result is expected, because Al is smaller than Mg and is affected by size changes in the core. In an analogy of the solute/dislocation interaction and in terms of the size misfit component in the present work, atoms with larger atomic radius than Mg is expected to attract atoms with smaller atomic radius than Mg to relax local elastic strain. However, pairings such as In–RE solute pairs are exceptions, because both atoms possess larger atomic radius than Mg. As a result, the empirical principle that large and small solute atoms tend to segregate together to relax the elastic strain cannot always be used to correctly predict co-segregation behaviors. The large discrepancy of solute/dislocation interaction energy by using direct first-principles calculations and size misfit approximation also indicate that some other factors, such as the chemical misfit, are also of considerable importance in determining interaction energy. The image in Fig. 3 shows the binding energies of 2NN solute pairs as a function of atomic radius difference. All solute pairs with positive binding energies are composed of Ca/Yb–RE, both of which possess larger atomic radius than Mg. Interestingly, there is a notable negative correlation between 2NN binding energies and the atomic radius difference of solute pairs.
4.1.2. Electronegativity The influence of electronegativity difference on the binding energies of solute pairs is plotted in Fig. 4. As illustrated in Fig. 4, a positive correlation exists between electronegativity and binding energies for 1NN solute pairs. As electronegativity describes the tendency of an atom to attract electrons towards itself, we reasonably predict that greater difference in electronegativity between two atoms translates to the formation of a stronger attraction bond. Mendis et al. [15] reported that solutes Ca and In tend to segregate together to form a plate-like precipitate, which seems unreasonable according to the size misfit. The atomic radii of Ca, In, and Mg are 197, 167, and 160 pm, respectively. Both Ca and In are expected to result in a compression strain, so how can these solutes segregate together? Electronegativity difference may explain this phenomenon. The electronegativity values are 1.00 for Ca, 1.78 for In, and 1.31 for Mg. The strong electronegativity difference between Ca and In may thus be responsible for the favorable binding between Ca and In. The variation of binding energies of 2NN solute pairs with the electronegativity difference of solute atoms is shown in Fig. 5. In general, binding energy increases with the increase in
Table 4 Binding energies for various NN solute pairs computed in a supercell with 64 atoms. Note that positive energies indicate energetically favorable binding. 1NN
2NN Ag
Al
Bi
Solute–solute binding energy/eV Ag 0.20 0.22 0.25 0.23 Al 0.27 0.21 0.23 0.24 Bi 0.30 0.32 0.22 0.32 Ca 0.02 0.06 0.15
Ca
Dy
Er
Ga
Gd
Ho
0.04
0.15
0.15
0.03
0.16
0.15
0.04
0.14
0.14
0.03
0.15
0.01
0.12
0.11
0.05
0.07
0.07
In
Lu
Nd
Pb
Sc
Sm
Sn
Tm
Y
Yb
Zn
Zr
0.00
0.05
0.14
0.24
0.24
0.24
0.12
0.23
0.13
0.25
0.04
0.22
0.24
0.14
0.00
0.07
0.13
0.23
0.23
0.23
0.21
0.22
0.12
0.24
0.04
0.21
0.24
0.12
0.11
0.00
0.09
0.11
0.20
0.21
0.23
0.18
0.22
0.10
0.21
0.00
0.24
0.24
0.04
0.08
0.07
0.02
0.01
0.06
0.08
0.00
0.04
0.08
0.02
0.07
0.07
0.06
0.05
0.06
0.02
0.08
0.01
0.02
0.02
0.07
0.03
0.01
0.11
0.05
0.01
0.13
0.01
0.01
0.07
0.14
0.02
0.07
0.02
0.02
0.02
0.06
0.04
0.01
0.10
0.02
0.04
0.12
0.01
0.02
0.07
0.14
0.02
0.09
0.08
0.00
0.05
0.07
0.10
0.05
0.06
0.09
0.04
0.07
0.08
0.04
0.02
0.06
0.02
0.03
0.07
0.03
0.00
0.11
0.04
0.00
0.13
0.00
0.01
0.08
0.14
0.01
0.02
0.07
0.03
0.00
0.10
0.04
0.01
0.12
0.01
0.01
0.07
0.14
0.01
0.02
0.02
0.03
0.00
0.02
0.03
0.00
0.02
0.03
0.02
0.00
0.03
0.06
0.07
0.08
0.06
0.07
0.08
0.06
0.07
0.01
0.06
0.05
0.02
0.10
0.05
0.03
0.12
0.02
0.03
0.06
0.13
0.03
0.18
0.12
0.05
0.21
0.00
0.07
0.08
0.23
0.08
0.22
0.17
0.22
0.09
0.20
0.00
0.23
0.23
0.10
0.23
0.03
0.14
0.04
0.23
0.14
0.19
0.00
0.06
0.08
0.21
0.07
0.10
0.22
0.01
0.22
0.24
0.00
0.07
0.12
0.00
0.08
0.23
0.11
0.04
0.06
0.22 0.25 0.17
0.23
Dy
0.12
0.02
0.01
0.06 0.10 0.13
Er
0.11
0.02
0.01
0.13
0.02 0.18 0.17
Ga
0.08
0.06
0.11
0.09
0.07
0.03 0.16 0.07
Gd
0.12
0.02
0.02
0.14
0.18
0.17
0.03 0.06 0.07
Ho
0.11
0.02
0.01
0.13
0.17
0.16
0.07
0.01 0.19 0.18
In
0.05
0.06
0.09
0.09
0.09
0.09
0.04
0.09
0.02 0.16 0.09
Li
0.03
0.06
0.00
0.05
0.17
0.17
0.02
0.17
0.17
0.00 0.07 0.03
Lu
0.11
0.02
0.01
0.12
0.15
0.14
0.06
0.15
0.14
0.08
0.05 0.06 0.17
Nd
0.21
0.06
0.05
0.14
0.19
0.18
0.08
0.19
0.18
0.10
0.17
0.04 0.12 0.16
Pb
0.29
0.29
0.31
0.12
0.01
0.01
0.09
0.01
0.01
0.07
0.01
0.02
0.06 0.25 0.08
Sc
0.23
0.10
0.13
0.13
0.16
0.14
0.05
0.16
0.14
0.07
0.17
0.12
0.25
0.22 0.30 0.15
Sm
0.19
0.04
0.04
0.14
0.19
0.18
0.08
0.19
0.18
0.10
0.17
0.16
0.25
0.07
0.17 0.23 0.24
Sn
0.30
0.29
0.34
0.13
0.04
0.04
0.09
0.05
0.05
0.08
0.02
0.03
0.03
0.32
0.10
0.06 0.25 0.01
Tm
0.09
0.04
0.02
0.13
0.14
0.13
0.07
0.15
0.13
0.09
0.17
0.11
0.15
0.00
0.12
0.15
0.22 0.33 0.06
Y
0.21
0.08
0.08
0.14
0.19
0.18
0.07
0.20
0.18
0.09
0.17
0.16
0.28
0.10
0.26
0.26
0.05
0.01 0.12 0.15
Yb
0.05
0.08
0.20
0.09
0.09
0.09
0.11
0.09
0.09
0.10
0.05
0.09
0.08
0.16
0.11
0.09
0.17
0.09
0.10 0.30 0.10
Zn
0.24
0.24
0.27
0.05
0.04
0.03
0.06
0.03
0.03
0.04
0.04
0.04
0.11
0.26
0.16
0.10
0.27
0.02
0.12
0.06 0.06 0.09
Zr
0.28
0.07
0.19
0.16
0.10
0.08
0.05
0.11
0.08
0.05
0.21
0.05
0.21
0.21
0.15
0.19
0.11
0.05
0.21
0.11
G. Liu et al. / Computational Materials Science 103 (2015) 97–104
Li
0.11 0.01
101
102
G. Liu et al. / Computational Materials Science 103 (2015) 97–104
Fig. 2. Variation of calculated binding energies with atomic radius difference for 1NN solute pairs with favorable bindings.
underlying mechanism of the correlation between the binding energy and electronegativity difference for 2NN solute pairs requires further clarification. 4.2. Electronic properties
Fig. 3. Variation of calculated binding energies with atomic radius difference for 2NN solute pairs with favorable bindings.
To reveal a deeper insight into underlying factors influencing binding energies, the charge density difference contours of Mg– Bi–Yb and Mg–Bi–Sn alloy systems were taken as examples, considering that Bi–Yb and Bi–Sn 1NN solute pairs show the largest positive and negative values of binding energy, respectively. The geometry of the as-built 3 2 2 Mg supercell was initially optimized. Then, this structure was used for the charge density calculation of pure Mg, Mg–Bi, Mg–Yb, Mg–Sn, Mg–Bi–Yb, and Mg–Bi–Sn alloys with lattice parameters and atomic positions fixed at their positions. To accentuate the redistribution of charge density upon the formation of solute pairs, ground state non-interacting charge densities were subtracted from the total charge density. Hence, the point-by-point charge difference is given by:
dqðMg62 X 1 Z 1 Þ ¼ qðMg62 X 1 Z 1 Þ þ qðMg64 Þ qðMg63 X 1 Þ electronegativity difference. Besides, all solutes, except for Sc and Zr, were found to possess lower electronegativity than Mg. The interaction between atoms in 2NN pairs is somewhat complicated, because these solutes are not bonded directly. In fact, four coneighbored Mg atoms are found between 2NN solutes. The
qðMg63 Z 1 Þ The images in Figs. 6 and 7 show charge density difference contours for the (0 0 0 1) plane of the Mg–Bi–Yb and Mg–Bi–Sn alloys, respectively. VESTA 3.0 was used as visualization software [44]. In
G. Liu et al. / Computational Materials Science 103 (2015) 97–104
103
Fig. 4. Variation of calculated binding energies with electronegativity difference for 1NN solute pairs with favorable bindings.
Fig. 5. Variation of calculated binding energies with electronegativity difference for 2NN solute pairs with favorable bindings.
the contour key in Figs. 6 and 7, positive values denote charge accumulation, whereas negative values denote charge depletion. As shown in Fig. 6, charge accumulates between the solute Bi and Yb atoms, as denoted by the red contours, and covalent bonds were formed. The strong interaction between Bi and Yb is thus
Fig. 6. Charge density difference contours showing a (0 0 0 1) basal plane, in which solutes Bi (the purple sphere) and Yb (the cyan sphere) were located in the Mg–Bi– Yb alloy. Orange spheres are Mg atoms. Charge accumulation between Bi and Yb is evident. The color contour key denotes the range of contour values. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
expected to result in a large positive binding energy. Contrastingly in Fig. 7, a dark blue region corresponding to charge depletion appears between the Bi–Sn solute pair, which means that Bi and Sn tend to repel each other.
104
G. Liu et al. / Computational Materials Science 103 (2015) 97–104
References
Fig. 7. Charge density difference contours showing a (0 0 0 1) basal plane, in which solutes Bi (the purple sphere) and Sn (the gray sphere) were located in the Mg–Bi– Sn alloy. Orange spheres denote Mg atoms. The contour region of charge depletion (dark blue) between Bi and Sn is observed. The color contour key denotes the range of contour values. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
5. Conclusions In summary, a large database of solute–solute binding energies in Mg-based alloys from first-principles calculations was presented. The calculated results allow us to phenomenologically predict the precipitation sequence and changes in stacking fault energies. Solute pairs with large positive (e.g. Yb–Bi/Sn/Pb and Ca–Bi/Sn/Pb) and negative (e.g. Bi–Sn/Pb/Al) values of binding energies are predicted to show promise in modifying precipitation sequence and stacking fault energy. In addition, the predicted solute–solute bindings can serve as first-principles input for future work involving quantitative models to quantitatively predict multi-solute effects on mechanical properties, such as strength of Mg alloys. The analysis of physical factors controlling solute–solute bindings shows a slight correlation between the atomic radius difference and binding energy for 1NN solute pairs, whereas a linearly negative correlation between these two factors was found for 2NN solute pairs. In addition, the calculated results indicate that binding energy increases with the increase in electronegativity difference for 1NN solute pairs, while the value decreases for 2NN solute pairs. Favorable binding is caused by charge accumulation and the formation of covalent bonds, whereas unfavorable binding is attributed to charge depletion between solutes. Acknowledgements The authors are grateful for the financial support from the National Natural Science Foundation of China (Nos. 51271207 and 51471038).
[1] V. Savic, L.G. Hector, S. Kim, R. Verma, J. Eng. Mater. Technol. 132 (2010). 0210061-10. [2] M. Marya, L.G. Hector, R. Verma, W. Tong, Mater. Sci. Eng., A 418 (2006) 341– 356. [3] S.R. Agnew, J.F. Nie, Scripta Mater. 63 (2010) 671–673. [4] B.L. Mordike, T. Ebert, Mater. Sci. Eng., A 302 (2011) 37–45. [5] J.A. Yasi, L.G. Hector, D.R. Trinkle, Acta Mater. 58 (2010) 5704–5713. [6] J.A. Yasi, L.G. Hector, D.R. Trinkle, Acta Mater. 59 (2011) 5652–5660. [7] J.A. Yasi, L.G. Hector, D.R. Trinkle, Acta Mater. 60 (2012) 2350–2358. [8] J. Koike, T. Kobayashi, T. Mukai, H. Watanabe, M. Suzuki, K. Maruyama, et al., Acta Mater. 51 (2003) 2055–2065. [9] J.F. Nie, Metall. Mater. Trans. A 43 (2012) 3891–3939. [10] K. Hono, C.L. Mendis, T.T. Sasaki, K. Oh-ishi, Scripta Mater. 63 (2010) 710–715. [11] J.F. Nie, X. Gao, S.M. Zhu, Scripta Mater. 53 (2005) 1049–1053. [12] X. Gao, S.M. Zhu, B.C. Muddle, J.F. Nie, Scripta Mater. 53 (2005) 1321–1326. [13] C.L. Mendis, C.J. Bettles, M.A. Gibson, C.R. Hutchinson, Mater. Sci. Eng., A 435– 436 (2006) 163–171. [14] C.L. Mendis, C.J. Bettles, M.A. Gibson, S. Gorsse, C.R. Hutchinson, Philos. Mag. Lett. 86 (2006) 443–456. [15] C.L. Mendis, K. Oh-ishi, T. Ohkubo, K. Hono, Scripta Mater. 64 (2011) 137–140. [16] T.T. Sasaki, K. Oh-ishi, T. Ohkubo, K. Hono, Scripta Mater. 55 (2006) 251–254. [17] T.T. Sasaki, J.D. Ju, K. Hono, K.S. Shin, Scripta Mater. 61 (2009) 80–83. [18] G.P.M. Leyson, L.G. Hector, W.A. Curtin, Acta Mater. 60 (2012) 5197–5203. [19] T. Nogaret, W.A. Curtin, J.A. Yasi, L.G. Hector, D.R. Trinkle, Acta Mater. 58 (2010) 4332–4343. [20] M. Muzyk, Z. Pakiela, K.J. Kurzydlowski, Scripta Mater. 66 (2012) 219–222. [21] J. Zhang, Y. Dou, G. Liu, Z. Guo, Comput. Mater. Sci. 79 (2013) 564–569. [22] J. Zhang, Y. Dou, H. Dong, Scripta Mater. 89 (2014) 13–16. [23] W.Y. Wang, S.L. Shang, Y. Wang, Z.-G. Mei, K.A. Darling, L.J. Kecskes, et al., Mater. Res. Lett. 2 (2014) 29–36. [24] M. Ghazisaeidi, L.G. Hector, W.A. Curtin, Acta Mater. 80 (2014) 278–287. [25] G.P.M. Leyson, W.A. Curtin, L.G. Hector, C.F. Woodward, Nat. Mater. 9 (2010) 750–755. [26] Q. Zhang, L. Fu, T.W. Fan, B.Y. Tang, L.M. Peng, W.J. Ding, Physica B 416 (2013) 39–44. [27] D. Shin, C. Wolverton, Acta Mater. 58 (2010) 531–540. [28] J.E. Saal, C. Wolverton, Acta Mater. 60 (2012) 5151–5159. [29] H. Kimizuka, M. Fronzi, S. Ogata, Scripta Mater. 69 (2013) 594–597. [30] G. Kresse, J. Hafner, Phys. Rev. B 47 (1993) 558–561. [31] G. Kresse, J. Furthmuller, Comput. Mater. Sci. 6 (1996) 15–50. [32] G. Kresse, J. Furthmuller, Phys. Rev. B 54 (1996) 11169–11186. [33] G. Kresse, D. Joubert, Phys. Rev. B 59 (1999) 1758–1775. [34] S. Ganeshan, L.G. Hector, Z.K. Liu, Comput. Mater. Sci. 50 (2010) 301–307. [35] M. Methfessel, A.T. Paxton, Phys. Rev. B 40 (1989) 3616–3621. [36] P.E. Blöchl, Phys. Rev. B 50 (1994) 17953–17979. [37] M. Gao, A. Rollet, M. Widom, Phys. Rev. B 75 (2007) 174120. [38] L. Huber, I. Elfimov, J. Rottler, M. Militzer, Phys. Rev. B 85 (2012) 144301. [39] T.T. Sasaki, K. Oh-ishi, T. Ohkubo, K. Hono, Mater. Sci. Eng., A 530 (2011) 1–8. [40] J.G. Speigth (Ed.), Lange’s Handbook of Chemistry, 16th ed., McGraw-Hill Companies, New York, 2004. [41] C. Mendis, K. Ohishi, Y. Kawamura, T. Honma, S. Kamado, K. Hono, Acta Mater. 57 (2009) 749–760. [42] J. Jayaraj, C.L. Mendis, T. Ohkubo, K. Oh-ishi, K. Hono, Scripta Mater. 63 (2010) 831–834. [43] M. Yuasa, M. Hayashi, M. Mabuchi, Y. Chino, Acta Mater. 65 (2014) 207–214. [44] K. Momma, F. Izumi, J. Appl. Crystallgr. 44 (2011) 1272–1276.
Contents lists available at ScienceDirect
Computational Materials Science journal homepage: www.elsevier.com/locate/commatsci
First-principles study of solute–solute binding in magnesium alloys Guobao Liu a, Jing Zhang a,b,⇑, Yuchen Dou a a b
College of Materials Science and Engineering, Chongqing University, Chongqing 400044, China National Engineering Research Center for Magnesium Alloys, Chongqing 400044, China
a r t i c l e
i n f o
Article history: Received 23 September 2014 Received in revised form 9 March 2015 Accepted 17 March 2015 Available online 3 April 2015 Keywords: Magnesium alloys Solute–solute binding First-principles calculations
a b s t r a c t Solute–solute interactions play a major role in the properties of materials. In this work, we present an extensive database of solute–solute binding energies that captures the detailed interactions in Mg-based alloys from first-principles calculations based on density functional theory. The effects of solute–solute binding energies on magnesium properties, precipitation hardening responses and stacking fault energies in particular, are inferred and discussed. The results of our calculations regarding bindings between solutes with different chemistries, including Al–Sn, Al–Ca, Ca–Zn, Ca–In, and Sn–Zn, were validated using available experimental investigations. Solute pairs that were predicted to show large positive (e.g., Yb–Bi/ Sn/Pb and Ca–Bi/Sn/Pb) and negative (e.g., Bi–Sn/Pb/Al) values of binding energies exhibited potential in modifying the precipitation sequence and stacking fault energy. Moreover, alloys added with these alloying elements may exhibit unique mechanical properties, which await experimental verification. Finally, the effect of physical features, including atomic radius and electronegativity, on the solute–solute bindings was investigated. Ó 2015 Elsevier B.V. All rights reserved.
1. Introduction Magnesium alloys in wrought forms (e.g., extrusions, forgings, sheet, and plates) have gained considerable worldwide interest in the past decade as potential replacements for heavier steel and aluminum alloys because of their high strength–density ratio [1–3]. However, the vast majority of Mg alloy applications are presently covered by cast products [1,4]. The application of wrought Mg alloys is hindered in part by their poor room temperature formability caused by their strong basal type texture and large anisotropy between basal and non-basal (prismatic and pyramidal) slips. The critical resolved shear stress of a basal slip system in Mg at room temperature is approximately two orders of magnitude lower than those of non-basal slip systems on prismatic and pyramidal planes and somewhat less than that of twinning. Therefore, basal slip and twinning are almost the only deformation mechanisms in polycrystalline Mg alloys at room temperature, which is far from the five independent operating systems required by the von Mises’ criterion for sufficient ductility [5–8]. Meanwhile, the strength for most commercial Mg alloys is still substantially lower than their aluminum alloy counterparts, which is in part attributed
⇑ Corresponding author at: College of Materials Science and Engineering, Chongqing University, Chongqing 400044, China. Tel.: +86 23 65111167; fax: +86 23 65102821. E-mail address: [email protected] (J. Zhang). http://dx.doi.org/10.1016/j.commatsci.2015.03.023 0927-0256/Ó 2015 Elsevier B.V. All rights reserved.
to the low precipitation hardening responses of these alloys [3,9,10]. Therefore, developing new Mg alloys with higher strength and ductility is necessary for Mg alloys to remain competitive with aluminum alloys. In principle, the addition of alloying elements capable of modifying the precipitation sequence and the strengths of various deformation modes, thereby reducing crystallographic anisotropy, holds promise as a means to improve the strength and ductility of Mg alloys. Over the past decades, numerous experimental studies on the influence of alloying additions on precipitation hardening responses of Mg alloys have been conducted [9–17]. The agehardenability of Mg alloys improves significantly when the alloying elements for addition are selected appropriately. The interaction between solutes distinctly affects precipitation hardening responses by forming new types of precipitates or manipulating the density, distributions, sizes, morphologies, and growth orientation of the precipitates. For instance, the ultimate tensile stress of peak-aged Mg–6Gd–2Zn–0.6Zr (wt.%) alloys is nearly two times higher than its Mg–6Gd–0.6Zr (wt.%) counterpart because of the formation of the c0 (Mg70Gd15Zn15, ordered hcp) phase [11]. Small additions of Zn in Mg-rare earth (RE) [11], Mg–Ca [12], and Mg–Sn [13] alloys have been widely reported to promote the formation of fine-scale precipitates and increase the precipitate number density during ageing heat treatment, thereby increasing the peak-aged hardness of alloys. Mendis et al. [14] observed that the number density of Mg2Sn precipitates increased by
98
G. Liu et al. / Computational Materials Science 103 (2015) 97–104
approximately one order of magnitude in Mg–Sn alloys with the collective addition of Li + In, resulting in a hardening increment increase of 150%. Clustering or co-clustering of solute atoms has been proposed to possibly serve as heterogeneities for precipitate nucleation, leading to increased number density of the precipitates. The increased number density of thin prismatic plates in Mg–Ca alloys with added In has also been reported to enhance age-hardening response by three times [15]. It was suggested that the large strain field expected from Ca–In co-clusters changes the habit planes from basal to prismatic. Alloying elements in the form of substitutional solutes also affect the strength and ductility of Mg alloys by modifying certain underlying factors, such as stacking fault energy and solute/dislocation interaction energy [5–7,18]. Stacking fault energy is directly associated with the dissociation behavior of dislocations, further affecting the strength, ductility, and fracture of materials. Nogaret et al. [19] correlated the dislocation structures and Peierls stresses to gamma surfaces (generalized stacking fault surfaces) in pure Mg. The simulation result, in which pyramidal I {1 0 1 1}h1 1 2 3i slip efficiently occurs in comparison with the pyramidal II {1 1 2 2}h1 1 2 3i slip, correlates well with the calculated gamma surfaces. Theoretical calculations of stacking fault energies of several binary Mg–X alloys have been reported in literature [20–23]. The effects of solutes on the basal [5,18], prismatic [6,7], and twinning dislocations [24] by inputting the firstprinciples data (e.g., size misfit, chemical misfit derived from stacking fault energy, and solute/dislocation interaction energy) into strengthening models have been quantitatively calculated in recent years. At present, investigating the multi-solute effects for Mg alloys is severely limited, although these effects have been recently reported for Al alloys [25]. However, the interaction between solute atoms plays an important role in stacking fault energy. For example, the unstable stacking fault energy of basal slip decreased more significantly for systems with simultaneous addition of Y and Zn than those with single Zn or Y additions [26]. The increased separation distance between two Shockley partial dislocations in Mg–Zn–Y alloys observed using transmission electron microscopy also demonstrates that the combined addition of Zn and Y leads to a significant reduction in stacking fault energy. The interaction between solutes plays an important role in determining the underlying factors, such as precipitation process and stacking fault energy, further affecting the mechanical properties of Mg alloys. Binding energy, in terms of solute–solute/vacancy binding, captures the bonding properties of atoms and provides detailed information on the interaction between solute atoms. For the precipitation sequence, solutes with favorable bindings are likely to attract and form new types of precipitates, whereas those with unfavorable bindings tend to repel each other, thereby decreasing the solubility of solutes in the matrix and promoting precipitate density during aging treatments. For stacking fault energy, if two solutes with favorable binding are placed in two different basal planes, extra energy is definitely needed to shear crystals on this plane. Moreover, solute–solute binding energies can serve as input first-principles data to quantitatively predict the multi-solute effects for Mg alloys if a constitutive model is developed. Although solute-vacancy binding energies in Mg have been systematically studied by Shin and Wolverton [27] and Saal and Wolverton [28], calculation for solute–solute binding energies, at least to the best of our knowledge, only touches upon Al, Zn, Y, and Gd [29]. Therefore, it is necessary to take an extensive calculation of the binding energies of solute atoms in Mg to investigate the interaction between solute atoms. In this paper, first-principles density functional calculations were carried out to calculate the binding energies of various solute pairs to construct an extensive first-principle database for predicting solute–solute binding energies in Mg alloy. Alloying elements
with an equilibrium solid solubility larger than 0.5 (at.%) in Mg were selected: Ag, Al, Bi, Ca, Dy, Er, Ga, Gd, Ho, Li, Lu, In, Pb, Nd, Sc, Sm, Sn, Tm, Y, Yb, Zn, and Zr. The solute–solute binding energy trends were then investigated in terms of physical features, such as atomic radius difference and electronegativity difference of solute pairs. The charge density difference in selected planes of Mg–Bi–Yb and Mg–Bi–Sn alloys were computed to reveal the underlying mechanism of solute effects on the binding energies of solute atoms. 2. Computational details In this study, first-principles total-energy DFT calculations were performed using the Vienna Ab-initio Simulation Package. The Perdew–Burke–Ernzerhof (PBE) exchange–correlation functional for the generalized-gradient-approximation (GGA) and projectoraugmented wave method was used in the present work [30–33]. All calculations were fully relaxed with respect to all degrees of freedom. The first-order Methfessel–Paxton with SIGMA = 0.2 eV was used for geometric relaxation until the electronic energy converged to less than 105 eV/cell and the Hellmann–Feynman force on all atomic sites was less than 0.01 eV Å1. Two successive structural optimizations (adapting basis vectors and computational grids to cell parameters) were conducted to ensure that the cell energies and structural parameters fully converged [34]. Total energy calculation was then performed using linear tetrahedron method with Blöchl correction [35,36]. A plane wave cutoff energy, which is 1.3 times the largest cutoff energy associated with the elements of interest recommended by VASP, was automatically performed by setting the PREC = HIGH in INCAR. It is known that the f electrons are not handled well by presently available density functions because of their varying tendency to form localized states. A common solution to this problem is to place the f electrons in the core (so-called ‘‘frozen’’ potentials). The spatial localization of f electrons near the core is a good approximation and works well in Al–RE [37], Mg–RE [38], Mg–Pb, and Mg–Bi [7] alloys. Therefore, frozen potentials are used for energy calculation in systems embedded with f elements (RE, Pb, and Bi in the present work). Table 1 lists the PAW potentials and the corresponding cutoff energies for Mg and 22 solutes Table 1 The PAW valence configuration, energy cutoff for the total energy calculation of pure Mg, Mg–X, and Mg–X–Z, and the atomic radius and electronegativity of Mg and the alloying elements considered in the present work. Alloying element
PAWPP
Cutoff (eV)
Atomic radius (pm)
Electronegativity
Mg Ag Al Bi Ca Dy Er Ga Gd Ho In Li Lu Nd Pb Sc Sm Sn Tm Y Yb Zn Zr
[Ne]3s2 [Kr]4d105s1 [Ne]3s23p1 ([Xe]4f145d10)6s26p3 [Ar]4s2 ([Xe]4f10)6s2 ([Xe]4f12)6s2 ([Ar]3d10)4s24p1 ([Xe]4f7)5d16s2 ([Xe]4f11)6s2 ([Kr]4d10)5s25p1 [He]2s1 ([Xe]4f14)5d16s2 ([Xe]4f4)6s2 ([Xe]4f145d10)6s26p2 [Ar]3d24s1 ([Xe]4f6)6s2 ([Kr]4d10)5s25p2 ([Xe]4f13)6s2 [Kr]4d25s1 ([Xe]4f14)6s2 [Ar] 3d104s2 [Kr]4d35s1
164 325 312 315 155 202 201 175 201 200 125 182 201 237 309 289 230 313 193 263 146 360 299
160 144 143.1 154.7 197 178.1 176.1 135 180.4 176.2 167 173.8 152 181.4 175 162 180.4 151 175.9 180 193.3 134 160
1.31 1.93 1.61 2.02 1 1.22 1.24 1.81 1.2 1.23 1.78 0.98 1.27 1.14 1.87 1.36 1.17 1.96 1.25 1.22 1.1 1.65 1.33
G. Liu et al. / Computational Materials Science 103 (2015) 97–104
considered in this work. Given that the unpaired 4f electrons are primarily responsible for magnetism, magnetic contribution to the total energy is negligible even for elements known to be ferromagnetic at room temperature, such as Gd, when the frozen potentials are used. As a result, non-spin polarized calculation was considered in all calculations in this work. Supercells for solute–solute binding energy calculations were constructed by duplicating the orthorhombic unit cell with four atoms along three lattice vectors (x is aligned with [11–20], y with [1 0 1 0], and z with [0 0 0 1]) as illustrated in Fig. 1. Table 2 summarizes the supercell dimensions and k-point grids used for solute–solute binding energy calculations in the present work. The binding energy of the solute–solute pairs are defined as follows:
EbðXXÞ ¼ EðMgN2 X2 Þ þ EðMgN Þ 2EðMgN1 X1 Þ EbðXZÞ ¼ EðMgN2 X2 Þ þ EðMgN Þ EðMgN1 X1 Þ EðMgN1 Z1 Þ where Eb is the binding energy, E is the total energy of pure Mg or systems embedded with solute atoms, X and Z represent two different solute atoms in the Mg supercell, and N is the number of lattice sites considered. A positive value of Eb corresponds to favorable binding under this definition. Given that alloying solute atoms affect binding energies by disturbing the local electronic environment and strain field, configurations in which two solute atoms are located respectively at the first- and second-nearest neighbor sites (1NN and 2NN) in the hexagonal lattice were considered in the current calculation. In hcp Mg, two nonequivalent 1NN solute pairs existed and distinguished as follows: one in the same basal (0 0 0 1) plane (1NN) and the other in a {1 1 2 2} pyramid plane (1NN0 ), as shown in Fig. 1a. 3. Results 3.1. Effects of supercell size The effects of supercell size on the binding energies of solute– solute pairs were first investigated. The binding energies of solute–solute pairs commonly reported to show remarkable hardenability in experiments, namely, Ag–Ca [10], Ca–In [15], Ca–Sn [39], Ca–Zn [12], and Li–In [14] for basal plane solute–solute pairs (1NN), were calculated in supercells with 48, 64, 96, and 144 atoms. The calculated solute–solute binding energies for different supercell sizes are listed in Table 3. The table shows that 64-atom calculations converged within 0.02 eV and are thus sufficient to provide reliable binding energies of solutes. As a result, the binding energies of all other solute pairs were calculated in 64-atom supercells, as shown in Fig. 1b. 3.2. Solute–solute binding energy Based on our calculated data, we found that the solute–solute binding energies of Mg were not sensitive to the difference between basal plane (1NN) and out-of-basal-plane (1NN0 ) nearest neighbors. The binding energies of 1NN0 solute pairs were therefore not considered in this study. The calculated binding energies for 1NN and 2NN solute pairs are summarized in Table 4. Basic physical parameters, including atom radii and eletronegativity values of various elements considered in the present work, are listed in Table 1[40]. Table 4 shows that 57 kinds of 1NN solute pairs exhibit positive values of binding energies, indicating that the atoms of these solute pairs are attractive and may get together to form pairs or clusters. Interestingly, most of these solute pairs are composed of large atoms (i.e., RE elements and Ca) and relatively small atoms
99
(Bi, Pb, Sn, In, and Ga) located in the IIIA to VA groups in the periodic table of chemical elements. The solute pair Bi–Yb possesses the largest binding energy in all 1NN solute pairs, with a value of 2.4 eV. It is worth to note that solute pairs such as Ag–Ca, Ca–In, Ca–Zn, and Li–In, which are commonly reported in literature to present unique age hardenability in Mg alloys, were predicted to exhibit favorable 1NN bindings. However, Mg–Sn–Zn alloy systems [17], which were also reported to exhibit impressive precipitation hardening responses, showed an extremely large negative binding energy value of –0.27 eV for Sn–Zn solute pairs. For 2NN solute pairs, the binding energies of 34 kinds of solute pairs are positive. Overall, the magnitude of 2NN binding energies was considerably smaller than that of 1NN binding energies. Notably, most of these pairs are composed of atoms with larger atomic radius than Mg. More specifically, combinations of RE and Ca/Yb provide such examples; among which, the Ca–Nd solute pair is the most unique, with the largest binding energy of 0.08 eV. It is interesting to note that Ca and Yb show unique properties among all the alloying elements considered in the present work. Both X–Ca and X–Yb 1NN solute pairs show relatively larger binding energies compared with other solute pairs, whereas solute pairs showing positive binding energies are all composed of Ca– RE or Yb–RE for 2NN. 3.3. Validating solute–solute bindings Direct experimental measurements for solute–solute binding have yet to be reported. However, the precipitation sequence and theoretical investigations for the stacking fault energies of certain Mg alloys are available to validate our calculation results. Mendis et al. [41] observed the precipitation sequence of the Mg–2.4Zn– 0.1Ag–0.1Ca–0.16Zr alloy. Their 3DAP results showed that Ca forms clusters with Zn in the early stage of precipitation, which agrees well with our calculation that Ca and Zn tend to form favorable bindings. Solute clusters, such as Ca–Al [42], Ca–In [15], and Li–In [14], were also experimentally observed after ageing heat treatments. Notably, solute clusters reported in the above works show considerable potential in enhancing precipitation hardening responses, either by forming new types of precipitates or serving as heterogeneities for the nucleation of precipitates, thereby increasing the density of precipitates. Sasaki et al. [17] and Mendis et al. [13] reported that the age-hardening responses of Mg–Sn alloys were substantially improved by small additions of Zn, which is mainly attributed to the increased number density of precipitates. However, no new types of precipitates were observed except for Mg2Sn and MgZn2. The as-calculated negative binding energy of Zn–Sn solute pair in this work indicates that Zn and Sn tend to repel each other, thereby decreasing the solubility of both Zn and Sn in the Mg matrix and promoting the nucleation of precipitates. Yuasa et al. [43] calculated the stacking fault energies of basal and prismatic slip systems for Mg–Zn–Ca alloys, in which the 1NN0 solute pairs were formed. Plastic anisotropy was reduced because of the non-linear nature (the addition of Zn in Mg–Ca alloy increased the unstable stacking fault energy of basal slip while decreasing that of prismatic) of the unstable stacking fault energy of basal slip. The increase in unstable stacking fault energy of basal slip agrees with our calculation, which shows that Ca and Zn tend to form favorable binding, and extra energy is needed to shear the basal plane between these two solutes. 4. Discussions 4.1. Physical contributions to solute–solute binding energies In this section, we attempt to find general trends and physical factors controlling solute binding energies across solute pairs with
100
G. Liu et al. / Computational Materials Science 103 (2015) 97–104
Fig. 1. (a) Schematic of the orthorhombic Mg unit cell (solid red lines) used to generate supercells in the present work. The solute atom, X (the purple sphere), is at the origin of the unit cell, and the solute Z is placed either at the first-nearest neighbor position in basal plane (1NN, the cyan sphere), in pyramid plane (1NN0 , the dark blue sphere), or at the second-nearest neighbor position (2NN, the green sphere). Orange spheres denote Mg atoms. (b) Super-cell with 64 atoms used to calculate binding energies of all solute pairs in the present work. The purple sphere denotes solute X, and the cyan sphere denotes solute Z. Orange spheres are Mg atoms. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
Table 2 Atom number (N) and corresponding solute concentrations, supercell dimensions, and k-point grids used in the calculation. Atom number (N)
Solute concentrations at.%
Supercell dimensions
k-point grids
48 64 96 144
4.17 3.13 2.08 1.39
322 422 423 433
11 11 11 999 777 555
favorable bindings that have been considered in the current work. The interaction between atoms is closely associated with the changing electronic environment and strain field caused by solutes. Yasi et al. [5] computed the so-called solute ‘‘misfits’’, which consist of the following two important components: a size misfit for the change in local volume and a chemical misfit for the change in energy to slip the crystal, to approximate the solute/dislocation interaction energy. For solute–solute interactions, we supposed that the atom radius and electronegativity of solute atoms are key factors that influence solute–solute binding energy. 4.1.1. Atomic radius The image in Fig. 2 shows the calculated 1NN binding energies with positive values as a function of atomic radius difference between solutes. Considering their extremely cluttered nature for clarity, the data were separated into Al–X, Bi–X, Ca–X, Ga–X, In– X, Sn–X, and Yb–X, and the corresponding scatter plots are shown in Fig. 2a–g, respectively. We found that most solute pairs with positive binding energies are composed of a large atom (RE/Ca) and a small atom belonging to IIIA to VA. The solute/dislocation interaction energy via first-principles calculations by Yasi et al. [5] also predicted that the strongest interaction for an Al solute Table 3 Binding energy convergence corresponds to the supercell size. The binding energies for the 1NN solute pairs of Ag–Ca, Ca–In, Ca–Sn, Ca–Zn, and Li–In are computed in supercell with 48, 64, 96, and 144 atoms. Note that positive energies indicate energetically favorable binding. Solute pairs
Ag–Ca Ca–In Ca–Sn Ca–Zn Li–In
Solute–solute binding energy (eV) 48 atoms
64 atoms
96 atoms
144 atoms
0.02 0.09 0.13 0.03 0.04
0.02 0.09 0.13 0.05 0.03
0.03 0.09 0.15 0.06 0.04
0.02 0.10 0.13 0.05 0.03
is at the site of maximum compression near the center of the dislocation partial cores for both edge and screw dislocations. This result is expected, because Al is smaller than Mg and is affected by size changes in the core. In an analogy of the solute/dislocation interaction and in terms of the size misfit component in the present work, atoms with larger atomic radius than Mg is expected to attract atoms with smaller atomic radius than Mg to relax local elastic strain. However, pairings such as In–RE solute pairs are exceptions, because both atoms possess larger atomic radius than Mg. As a result, the empirical principle that large and small solute atoms tend to segregate together to relax the elastic strain cannot always be used to correctly predict co-segregation behaviors. The large discrepancy of solute/dislocation interaction energy by using direct first-principles calculations and size misfit approximation also indicate that some other factors, such as the chemical misfit, are also of considerable importance in determining interaction energy. The image in Fig. 3 shows the binding energies of 2NN solute pairs as a function of atomic radius difference. All solute pairs with positive binding energies are composed of Ca/Yb–RE, both of which possess larger atomic radius than Mg. Interestingly, there is a notable negative correlation between 2NN binding energies and the atomic radius difference of solute pairs.
4.1.2. Electronegativity The influence of electronegativity difference on the binding energies of solute pairs is plotted in Fig. 4. As illustrated in Fig. 4, a positive correlation exists between electronegativity and binding energies for 1NN solute pairs. As electronegativity describes the tendency of an atom to attract electrons towards itself, we reasonably predict that greater difference in electronegativity between two atoms translates to the formation of a stronger attraction bond. Mendis et al. [15] reported that solutes Ca and In tend to segregate together to form a plate-like precipitate, which seems unreasonable according to the size misfit. The atomic radii of Ca, In, and Mg are 197, 167, and 160 pm, respectively. Both Ca and In are expected to result in a compression strain, so how can these solutes segregate together? Electronegativity difference may explain this phenomenon. The electronegativity values are 1.00 for Ca, 1.78 for In, and 1.31 for Mg. The strong electronegativity difference between Ca and In may thus be responsible for the favorable binding between Ca and In. The variation of binding energies of 2NN solute pairs with the electronegativity difference of solute atoms is shown in Fig. 5. In general, binding energy increases with the increase in
Table 4 Binding energies for various NN solute pairs computed in a supercell with 64 atoms. Note that positive energies indicate energetically favorable binding. 1NN
2NN Ag
Al
Bi
Solute–solute binding energy/eV Ag 0.20 0.22 0.25 0.23 Al 0.27 0.21 0.23 0.24 Bi 0.30 0.32 0.22 0.32 Ca 0.02 0.06 0.15
Ca
Dy
Er
Ga
Gd
Ho
0.04
0.15
0.15
0.03
0.16
0.15
0.04
0.14
0.14
0.03
0.15
0.01
0.12
0.11
0.05
0.07
0.07
In
Lu
Nd
Pb
Sc
Sm
Sn
Tm
Y
Yb
Zn
Zr
0.00
0.05
0.14
0.24
0.24
0.24
0.12
0.23
0.13
0.25
0.04
0.22
0.24
0.14
0.00
0.07
0.13
0.23
0.23
0.23
0.21
0.22
0.12
0.24
0.04
0.21
0.24
0.12
0.11
0.00
0.09
0.11
0.20
0.21
0.23
0.18
0.22
0.10
0.21
0.00
0.24
0.24
0.04
0.08
0.07
0.02
0.01
0.06
0.08
0.00
0.04
0.08
0.02
0.07
0.07
0.06
0.05
0.06
0.02
0.08
0.01
0.02
0.02
0.07
0.03
0.01
0.11
0.05
0.01
0.13
0.01
0.01
0.07
0.14
0.02
0.07
0.02
0.02
0.02
0.06
0.04
0.01
0.10
0.02
0.04
0.12
0.01
0.02
0.07
0.14
0.02
0.09
0.08
0.00
0.05
0.07
0.10
0.05
0.06
0.09
0.04
0.07
0.08
0.04
0.02
0.06
0.02
0.03
0.07
0.03
0.00
0.11
0.04
0.00
0.13
0.00
0.01
0.08
0.14
0.01
0.02
0.07
0.03
0.00
0.10
0.04
0.01
0.12
0.01
0.01
0.07
0.14
0.01
0.02
0.02
0.03
0.00
0.02
0.03
0.00
0.02
0.03
0.02
0.00
0.03
0.06
0.07
0.08
0.06
0.07
0.08
0.06
0.07
0.01
0.06
0.05
0.02
0.10
0.05
0.03
0.12
0.02
0.03
0.06
0.13
0.03
0.18
0.12
0.05
0.21
0.00
0.07
0.08
0.23
0.08
0.22
0.17
0.22
0.09
0.20
0.00
0.23
0.23
0.10
0.23
0.03
0.14
0.04
0.23
0.14
0.19
0.00
0.06
0.08
0.21
0.07
0.10
0.22
0.01
0.22
0.24
0.00
0.07
0.12
0.00
0.08
0.23
0.11
0.04
0.06
0.22 0.25 0.17
0.23
Dy
0.12
0.02
0.01
0.06 0.10 0.13
Er
0.11
0.02
0.01
0.13
0.02 0.18 0.17
Ga
0.08
0.06
0.11
0.09
0.07
0.03 0.16 0.07
Gd
0.12
0.02
0.02
0.14
0.18
0.17
0.03 0.06 0.07
Ho
0.11
0.02
0.01
0.13
0.17
0.16
0.07
0.01 0.19 0.18
In
0.05
0.06
0.09
0.09
0.09
0.09
0.04
0.09
0.02 0.16 0.09
Li
0.03
0.06
0.00
0.05
0.17
0.17
0.02
0.17
0.17
0.00 0.07 0.03
Lu
0.11
0.02
0.01
0.12
0.15
0.14
0.06
0.15
0.14
0.08
0.05 0.06 0.17
Nd
0.21
0.06
0.05
0.14
0.19
0.18
0.08
0.19
0.18
0.10
0.17
0.04 0.12 0.16
Pb
0.29
0.29
0.31
0.12
0.01
0.01
0.09
0.01
0.01
0.07
0.01
0.02
0.06 0.25 0.08
Sc
0.23
0.10
0.13
0.13
0.16
0.14
0.05
0.16
0.14
0.07
0.17
0.12
0.25
0.22 0.30 0.15
Sm
0.19
0.04
0.04
0.14
0.19
0.18
0.08
0.19
0.18
0.10
0.17
0.16
0.25
0.07
0.17 0.23 0.24
Sn
0.30
0.29
0.34
0.13
0.04
0.04
0.09
0.05
0.05
0.08
0.02
0.03
0.03
0.32
0.10
0.06 0.25 0.01
Tm
0.09
0.04
0.02
0.13
0.14
0.13
0.07
0.15
0.13
0.09
0.17
0.11
0.15
0.00
0.12
0.15
0.22 0.33 0.06
Y
0.21
0.08
0.08
0.14
0.19
0.18
0.07
0.20
0.18
0.09
0.17
0.16
0.28
0.10
0.26
0.26
0.05
0.01 0.12 0.15
Yb
0.05
0.08
0.20
0.09
0.09
0.09
0.11
0.09
0.09
0.10
0.05
0.09
0.08
0.16
0.11
0.09
0.17
0.09
0.10 0.30 0.10
Zn
0.24
0.24
0.27
0.05
0.04
0.03
0.06
0.03
0.03
0.04
0.04
0.04
0.11
0.26
0.16
0.10
0.27
0.02
0.12
0.06 0.06 0.09
Zr
0.28
0.07
0.19
0.16
0.10
0.08
0.05
0.11
0.08
0.05
0.21
0.05
0.21
0.21
0.15
0.19
0.11
0.05
0.21
0.11
G. Liu et al. / Computational Materials Science 103 (2015) 97–104
Li
0.11 0.01
101
102
G. Liu et al. / Computational Materials Science 103 (2015) 97–104
Fig. 2. Variation of calculated binding energies with atomic radius difference for 1NN solute pairs with favorable bindings.
underlying mechanism of the correlation between the binding energy and electronegativity difference for 2NN solute pairs requires further clarification. 4.2. Electronic properties
Fig. 3. Variation of calculated binding energies with atomic radius difference for 2NN solute pairs with favorable bindings.
To reveal a deeper insight into underlying factors influencing binding energies, the charge density difference contours of Mg– Bi–Yb and Mg–Bi–Sn alloy systems were taken as examples, considering that Bi–Yb and Bi–Sn 1NN solute pairs show the largest positive and negative values of binding energy, respectively. The geometry of the as-built 3 2 2 Mg supercell was initially optimized. Then, this structure was used for the charge density calculation of pure Mg, Mg–Bi, Mg–Yb, Mg–Sn, Mg–Bi–Yb, and Mg–Bi–Sn alloys with lattice parameters and atomic positions fixed at their positions. To accentuate the redistribution of charge density upon the formation of solute pairs, ground state non-interacting charge densities were subtracted from the total charge density. Hence, the point-by-point charge difference is given by:
dqðMg62 X 1 Z 1 Þ ¼ qðMg62 X 1 Z 1 Þ þ qðMg64 Þ qðMg63 X 1 Þ electronegativity difference. Besides, all solutes, except for Sc and Zr, were found to possess lower electronegativity than Mg. The interaction between atoms in 2NN pairs is somewhat complicated, because these solutes are not bonded directly. In fact, four coneighbored Mg atoms are found between 2NN solutes. The
qðMg63 Z 1 Þ The images in Figs. 6 and 7 show charge density difference contours for the (0 0 0 1) plane of the Mg–Bi–Yb and Mg–Bi–Sn alloys, respectively. VESTA 3.0 was used as visualization software [44]. In
G. Liu et al. / Computational Materials Science 103 (2015) 97–104
103
Fig. 4. Variation of calculated binding energies with electronegativity difference for 1NN solute pairs with favorable bindings.
Fig. 5. Variation of calculated binding energies with electronegativity difference for 2NN solute pairs with favorable bindings.
the contour key in Figs. 6 and 7, positive values denote charge accumulation, whereas negative values denote charge depletion. As shown in Fig. 6, charge accumulates between the solute Bi and Yb atoms, as denoted by the red contours, and covalent bonds were formed. The strong interaction between Bi and Yb is thus
Fig. 6. Charge density difference contours showing a (0 0 0 1) basal plane, in which solutes Bi (the purple sphere) and Yb (the cyan sphere) were located in the Mg–Bi– Yb alloy. Orange spheres are Mg atoms. Charge accumulation between Bi and Yb is evident. The color contour key denotes the range of contour values. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
expected to result in a large positive binding energy. Contrastingly in Fig. 7, a dark blue region corresponding to charge depletion appears between the Bi–Sn solute pair, which means that Bi and Sn tend to repel each other.
104
G. Liu et al. / Computational Materials Science 103 (2015) 97–104
References
Fig. 7. Charge density difference contours showing a (0 0 0 1) basal plane, in which solutes Bi (the purple sphere) and Sn (the gray sphere) were located in the Mg–Bi– Sn alloy. Orange spheres denote Mg atoms. The contour region of charge depletion (dark blue) between Bi and Sn is observed. The color contour key denotes the range of contour values. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
5. Conclusions In summary, a large database of solute–solute binding energies in Mg-based alloys from first-principles calculations was presented. The calculated results allow us to phenomenologically predict the precipitation sequence and changes in stacking fault energies. Solute pairs with large positive (e.g. Yb–Bi/Sn/Pb and Ca–Bi/Sn/Pb) and negative (e.g. Bi–Sn/Pb/Al) values of binding energies are predicted to show promise in modifying precipitation sequence and stacking fault energy. In addition, the predicted solute–solute bindings can serve as first-principles input for future work involving quantitative models to quantitatively predict multi-solute effects on mechanical properties, such as strength of Mg alloys. The analysis of physical factors controlling solute–solute bindings shows a slight correlation between the atomic radius difference and binding energy for 1NN solute pairs, whereas a linearly negative correlation between these two factors was found for 2NN solute pairs. In addition, the calculated results indicate that binding energy increases with the increase in electronegativity difference for 1NN solute pairs, while the value decreases for 2NN solute pairs. Favorable binding is caused by charge accumulation and the formation of covalent bonds, whereas unfavorable binding is attributed to charge depletion between solutes. Acknowledgements The authors are grateful for the financial support from the National Natural Science Foundation of China (Nos. 51271207 and 51471038).
[1] V. Savic, L.G. Hector, S. Kim, R. Verma, J. Eng. Mater. Technol. 132 (2010). 0210061-10. [2] M. Marya, L.G. Hector, R. Verma, W. Tong, Mater. Sci. Eng., A 418 (2006) 341– 356. [3] S.R. Agnew, J.F. Nie, Scripta Mater. 63 (2010) 671–673. [4] B.L. Mordike, T. Ebert, Mater. Sci. Eng., A 302 (2011) 37–45. [5] J.A. Yasi, L.G. Hector, D.R. Trinkle, Acta Mater. 58 (2010) 5704–5713. [6] J.A. Yasi, L.G. Hector, D.R. Trinkle, Acta Mater. 59 (2011) 5652–5660. [7] J.A. Yasi, L.G. Hector, D.R. Trinkle, Acta Mater. 60 (2012) 2350–2358. [8] J. Koike, T. Kobayashi, T. Mukai, H. Watanabe, M. Suzuki, K. Maruyama, et al., Acta Mater. 51 (2003) 2055–2065. [9] J.F. Nie, Metall. Mater. Trans. A 43 (2012) 3891–3939. [10] K. Hono, C.L. Mendis, T.T. Sasaki, K. Oh-ishi, Scripta Mater. 63 (2010) 710–715. [11] J.F. Nie, X. Gao, S.M. Zhu, Scripta Mater. 53 (2005) 1049–1053. [12] X. Gao, S.M. Zhu, B.C. Muddle, J.F. Nie, Scripta Mater. 53 (2005) 1321–1326. [13] C.L. Mendis, C.J. Bettles, M.A. Gibson, C.R. Hutchinson, Mater. Sci. Eng., A 435– 436 (2006) 163–171. [14] C.L. Mendis, C.J. Bettles, M.A. Gibson, S. Gorsse, C.R. Hutchinson, Philos. Mag. Lett. 86 (2006) 443–456. [15] C.L. Mendis, K. Oh-ishi, T. Ohkubo, K. Hono, Scripta Mater. 64 (2011) 137–140. [16] T.T. Sasaki, K. Oh-ishi, T. Ohkubo, K. Hono, Scripta Mater. 55 (2006) 251–254. [17] T.T. Sasaki, J.D. Ju, K. Hono, K.S. Shin, Scripta Mater. 61 (2009) 80–83. [18] G.P.M. Leyson, L.G. Hector, W.A. Curtin, Acta Mater. 60 (2012) 5197–5203. [19] T. Nogaret, W.A. Curtin, J.A. Yasi, L.G. Hector, D.R. Trinkle, Acta Mater. 58 (2010) 4332–4343. [20] M. Muzyk, Z. Pakiela, K.J. Kurzydlowski, Scripta Mater. 66 (2012) 219–222. [21] J. Zhang, Y. Dou, G. Liu, Z. Guo, Comput. Mater. Sci. 79 (2013) 564–569. [22] J. Zhang, Y. Dou, H. Dong, Scripta Mater. 89 (2014) 13–16. [23] W.Y. Wang, S.L. Shang, Y. Wang, Z.-G. Mei, K.A. Darling, L.J. Kecskes, et al., Mater. Res. Lett. 2 (2014) 29–36. [24] M. Ghazisaeidi, L.G. Hector, W.A. Curtin, Acta Mater. 80 (2014) 278–287. [25] G.P.M. Leyson, W.A. Curtin, L.G. Hector, C.F. Woodward, Nat. Mater. 9 (2010) 750–755. [26] Q. Zhang, L. Fu, T.W. Fan, B.Y. Tang, L.M. Peng, W.J. Ding, Physica B 416 (2013) 39–44. [27] D. Shin, C. Wolverton, Acta Mater. 58 (2010) 531–540. [28] J.E. Saal, C. Wolverton, Acta Mater. 60 (2012) 5151–5159. [29] H. Kimizuka, M. Fronzi, S. Ogata, Scripta Mater. 69 (2013) 594–597. [30] G. Kresse, J. Hafner, Phys. Rev. B 47 (1993) 558–561. [31] G. Kresse, J. Furthmuller, Comput. Mater. Sci. 6 (1996) 15–50. [32] G. Kresse, J. Furthmuller, Phys. Rev. B 54 (1996) 11169–11186. [33] G. Kresse, D. Joubert, Phys. Rev. B 59 (1999) 1758–1775. [34] S. Ganeshan, L.G. Hector, Z.K. Liu, Comput. Mater. Sci. 50 (2010) 301–307. [35] M. Methfessel, A.T. Paxton, Phys. Rev. B 40 (1989) 3616–3621. [36] P.E. Blöchl, Phys. Rev. B 50 (1994) 17953–17979. [37] M. Gao, A. Rollet, M. Widom, Phys. Rev. B 75 (2007) 174120. [38] L. Huber, I. Elfimov, J. Rottler, M. Militzer, Phys. Rev. B 85 (2012) 144301. [39] T.T. Sasaki, K. Oh-ishi, T. Ohkubo, K. Hono, Mater. Sci. Eng., A 530 (2011) 1–8. [40] J.G. Speigth (Ed.), Lange’s Handbook of Chemistry, 16th ed., McGraw-Hill Companies, New York, 2004. [41] C. Mendis, K. Ohishi, Y. Kawamura, T. Honma, S. Kamado, K. Hono, Acta Mater. 57 (2009) 749–760. [42] J. Jayaraj, C.L. Mendis, T. Ohkubo, K. Oh-ishi, K. Hono, Scripta Mater. 63 (2010) 831–834. [43] M. Yuasa, M. Hayashi, M. Mabuchi, Y. Chino, Acta Mater. 65 (2014) 207–214. [44] K. Momma, F. Izumi, J. Appl. Crystallgr. 44 (2011) 1272–1276.