Calculation of solubility in titanium alloys from first principles

Available online at www.sciencedirect.com

Acta Materialia 57 (2009) 5314–5323 www.elsevier.com/locate/actamat

Calculation of solubility in titanium alloys from first principles Roman V. Chepulskii, Stefano Curtarolo * Department of Mechanical Engineering and Materials Science and Department of Physics, Duke University, Durham, NC 27708, USA Received 1 January 2009; received in revised form 13 July 2009; accepted 18 July 2009 Available online 24 August 2009

Abstract We present an approach to calculate the atomic bulk solubility in binary alloys based on the statistical-thermodynamic theory of dilute lattice gas. The model considers all the appropriate ground states of the alloy and results in a simple Arrhenius-type temperature dependence determined by a ‘‘low-solubility formation enthalpy”. This quantity, directly obtainable from first principles, is defined as the composition derivative of the compound formation enthalpy with respect to nearby ground states. We apply the framework and calculate the solubility of the A solutes in A–Ti alloys (A = Ag, Au, Cd, Co, Cr, Ir, W, Zn). In addition to determination of unknown low-temperature ground states for the eight alloys, we find qualitative agreements with solubility experimental results. The presented formalism, correct in the low-solubility limit, should be considered as an appropriate starting point for estimation of whether a more computationally expensive formalism is needed. Ó 2009 Acta Materialia Inc. Published by Elsevier Ltd. All rights reserved. Keywords: Phase separation; Precipitation; Solubility

1. Introduction High-throughput ab initio methods are capable of predicting materials’ properties and searching for better systems. They are becoming important tools for scientists working on rational materials development. These methods allow researchers to exploit property correlations between different systems and to guide the prediction of new materials [1–6]. The main difference between the ‘‘several-calculations” and the ‘‘high-throughput” philosophies is that the first concentrates on the calculation of a particular property, while the latter focuses on the extraction of property correlations in the search for systems with target characteristics. Several applications of high-throughput ab initio approaches have been proposed in recent years. Examples are the ‘‘data-mining of quantum calculations” method leading to the principle-component analysis of the formation energies of many alloys in several configurations [4,5,7,8], the evolutionary approach for determining lattice *

Corresponding author. Tel.: +1 919 660 5506. E-mail address: [email protected] (S. Curtarolo).

Hamiltonians [9], the ‘‘Pareto-optimal” search for alloys and catalysts [10,11], the prediction of the lithium–boron superconductor [12], the ‘‘high-throughput Kohn-anomalies” search in ternary lithium-borides [13,14], and the ‘‘multi-optimization” techniques used for the study of high-temperature reactions in multicomponent hydrides [15–17]. To overcome the immense quantity of information necessary to build correlations, the high-throughput approach requires a rapid yet approximate estimation of the material’s properties. For example, the ‘‘Automatic-Flow” (AFLOW) Consortium Project, aimed to the parameterization of all the binary intermetallic alloys, contains data for more than 120,000 ab initio calculations [18]. In addition, ad hoc algorithms and appropriate computer software must be developed to analyze and to conclude about correlations. Finally, once the search space has been narrowed, a detailed study is employed to pinpoint the feasible candidate systems. This manuscript focuses on one of the requirements for high-throughput approaches: we present an ad hoc algorithm to calculate rapidly the solubility in binary alloys.

1359-6454/$36.00 Ó 2009 Acta Materialia Inc. Published by Elsevier Ltd. All rights reserved. doi:10.1016/j.actamat.2009.07.037

R.V. Chepulskii, S. Curtarolo / Acta Materialia 57 (2009) 5314–5323

The knowledge of solubility is crucial for designing new alloys with particular physical, chemical and mechanical properties. For example, if we have to enhance an alloy property by adding extra species as solutes, it is necessary to know the equilibrium solubility, to avoid supersaturation-precipitation and subsequent modification of the target property (aging). In superconducting materials research, the problem emerges frequently: often expensive and difficult experiments are undertaken to enhance the critical temperature [19,20]. In catalysis research, experiments and modeling have shown that size-dependent carbon solubility is responsible for thermodynamic instabilities hindering the catalytic activity of very small Fe and Fe:Mo clusters [21,22]. The solubility of Zr in Al has already been addressed with success within the regular solution model [23] fit to ab initio calculations [24]. Solubility can also be extracted from the knowledge of the phase diagram which, in the case of lattice-conserving alloys, can be generated with the cluster expansion method [25,26] to obtain interactions, and with Monte Carlo calculations to determine free energies and phase boundaries [27,28]. However, a straightforward and rapid formalism leading to the estimation of equilibrium solubility for general alloys is still lacking. In the present paper we devised a statistical-thermodynamic approach to tackle the task. Our approach considers all available binary ground states rather than just pure solids. To test the method, we present calculations for a number of binary titanium systems. The paper is organized as follows. In Section 2 we write the equations governing solubility of vacancies and substitutional impurities in binary alloys. Section 3 is devoted to the discussion of capabilities and limits of the formalism. In Section 4, the phase diagrams and solubilities are addressed for A–Ti systems (A = Ag, Au, Cd, Co, Cr, Ir, W, Zn). Conclusions are given in Section 5.

2.1. Enthalpy Let us consider a disordered dilute solution of A atoms and vacancies (V) in a pure B solid with a given Bravais crystal lattice (labeled with ‘‘dis”). Neglecting A–A, A–V and V–V interactions and assuming small A and V concentrations, the approximate enthalpy can be written as: H dis ¼ Edis þ pV dis ¼ H at B N þ H AB N A þ H VB N V ;

ð1Þ

where H at B is the enthalpy of the pure B solid per atom, H AB is the change in enthalpy of the solid upon substitution of one B atom with an A atom, and H VB is the change in enthalpy upon removal of one B atom and the introduction of a vacancy. In addition, N is the total numbers of crystal lattice sites, N a ða ¼ A; B; VÞ is the number of a-type atoms: H AB ¼ E0AB þ pv0AB ;

0 Eat a and va ða ¼ A; BÞ are the energy and volume per atom of the a-pure solid, vAB and v0VB represent the change of volume of the B-pure solid upon introduction of one A atom or one vacancy, and p is the pressure. The framework introduced by Eq. (1) is similar to Wagner–Schottky model of a system of non-interacting particles [29]. The quantities with superscript ‘‘0” are considered to be temperature, pressure and composition independent being calculated at zero temperature and pressure. Hence, the accuracy of the approximation is higher in the case of limited pressures or for systems with very high bulk modulus, where the elastic energy fraction of Eat a is negligible. The quantities in Eq. (2) can be easily approximated as differences of first-principles energies ðEÞ and volumes ðvÞ of large supercell (sc) with and without a defect (substitutional A atom or vacancy):

E0AB ’ Esc ½Bnsc 1 A  Esc ½Bnsc ; v0AB ’ vsc ½Bnsc 1 A  vsc ½Bnsc ; E0VB ’ Esc ½Bnsc 1   Esc ½Bnsc ;

at 0 H at B ¼ EB þ pvB ;

H VB ¼ E0VB þ pv0VB ;

ð2Þ

ð3Þ

v0VB ’ vsc ½Bnsc 1   vsc ½Bnsc : As the size of the supercell increases, the approximate quantities in Eq. (3) approach their exact values. In literature [30] and in this paper, E0aB ; v0aB and H aB ða ¼ A; VÞ are called the ‘‘raw” (composition unpreserving) a-defect formation energy, volume and enthalpy, respectively. The enthalpy per atom is obtained from Eq. (1) as: dis at at H dis at ¼ H =ðN A þ N B Þ ¼ H B þ H AB xA þ ðH VB þ H B ÞxV ;

ð4Þ

where xa ða ¼ A; B; VÞ are the atomic concentrations xa ¼ N a =ðN A þ N B Þ:

ð5Þ

Then, we follow the conventional use of the formation enthalpy DH dis at calculated with respect to the pure A and B solids [31]: dis at at DH dis at ¼ H at  xA H A  ð1  xA ÞH B ;

2. Solubility formalism

N ¼ N A þ N B þ N V;

5315

ð6Þ

at at 0 where H at A ¼ E A þ pvA and H B has been defined in Eq. (2). Combining Eqs. (4) and (6), we get ð7Þ DH dis at ¼ H A xA þ H V xV ;

where the quantities H A and H V , defined as at H A ¼ H AB  H at A þ HB ;

H V ¼ H VB þ H at B;

ð8Þ

are usually called A-defect and V-defect ‘‘true” (composition preserving) formation enthalpies [30]. They can be obtained from [31]:   @DH dis at  ða ¼ A; VÞ: ð9Þ Ha ¼ @xa xa !0 2.2. Equilibrium Gibbs free energy The formation Gibbs free energy, DGdis at , is defined as: dis dis DGdis at ¼ DH at  T DS at ;

ð10Þ

5316

R.V. Chepulskii, S. Curtarolo / Acta Materialia 57 (2009) 5314–5323

where DH dis at is described by Eq. (7). The formation entropy dis DS dis at ¼ S at can be obtained within the mean-field approximation as: X kBN ca ln ca ; ð11Þ DS dis at ¼  N A þ N B a¼A;B;V where T ; k B ; ca ða ¼ A; B; VÞ are the temperature, the Boltzmann constant and the site concentrations of atoms: ca ¼ N a =N :

ð12Þ

By changing variables from site concentrations ca to atomic concentrations xa : xa ¼ N a =ðN A þ N B Þ;

ca ¼ xa =ð1 þ xV Þ;

Δ

ð13Þ

we rewrite Eq. (10) as dis;A þ DGdis;V DGdis at ¼ DGat at ;

Δ

dis;A ¼ H A xA þ k B T ½xA ln xA þ ð1  xA Þ lnð1  xA Þ; DGat dis;V ¼ H V xV þ k B T ½xV ln xV  ð1 þ xV Þ lnð1 þ xV Þ: DGat

ð14Þ In an alloy with fixed atomic composition xA , the equilibrium concentration of vacancies xeq V is determined by minimizing the formation Gibbs free energy:  dis;V  @DGdis @DGat at  ¼ ¼ 0: ð15Þ  @xV xA @xV The manipulation of Eqs. (14) and (15) leads to:    1 HV  1 ¼ exp ; xeq V kBT   HV ; ’ exp  xeq j V k B T EV kBT    HV dis;V eq DGat ; ðxV Þ ¼ k B T ln 1  exp  kBT    HV dis;V eq  : ðxV Þ k T E ’ k B T exp  DGat B V kBT

ð16Þ

To conclude, the ‘‘true” vacancy formation enthalpy H V determines the equilibrium concentration of vacancies within an Arrhenius-type formalism (see also Ref. [31]). In the next section, we show that H A resolves the solubility in the case of a phase-separating alloy without intermediate ground states. 2.3. Solubility At a given temperature, the solubility of A atoms in a B solid phase, xsol AB , is defined as the maximum homogeneously achievable concentration of A without precipitation of a new phase (Fig. 1a). The accurate calculation of xsol AB requires the knowledge of the ‘‘closest” ground state (labeled as ‘‘gs”). It is followed by the minimization of the formation Gibbs free energy DGmix at ðxÞ for the mixture of (i) the disordered dilute solution of A atoms and vacancies in a B-rich solid phase at composition xAB (the

Fig. 1. (a) Temperature–concentration and (b) formation Gibbs free energy–concentration (at a given temperature T 0 ) graphs illustrating our solubility concepts.

‘‘dis”-phase of the previous section), and (ii) the on- or offstoichiometric ‘‘gs”-phase at composition xAgs . The lever rule gives the fractions of the two phases: xAgs  x x  x AB DGdis DGgs DGmix at ðxÞ ¼ at ðxAB Þ þ at ðxAgs Þ: xAgs  xAB xAgs  xAB ð17Þ The minimization is performed with respect to xAB and xAgs (x is the overall composition of A in the two-phase mixture, xAB < x < xAgs ): @DGmix at ¼0 @xAB

@DGmix at ¼ 0: @xAgs

ð18Þ

Combination of Eqs. (17) and (18) leads to the usual common-tangent rule (‘‘red line” in Fig. 1b): gs @DGgs DGdis @DGdis at ðxAgs Þ at ðxAB Þ  DGat ðxAgs Þ at ðxAB Þ ¼ ¼ : @xAgs xAB  xAgs @xAB

ð19Þ Substituting Eq. (14) into Eq. (19), we obtain 1

xsol AB ¼ ½expðH sol =k B T Þ þ 1 ;

ð20Þ

which approximates in an Arrhenius-type relation at low temperature:   ’ expðH sol =k B T Þ: ð21Þ xsol AB  k B T H sol

R.V. Chepulskii, S. Curtarolo / Acta Materialia 57 (2009) 5314–5323

5317

The quantity H sol is defined as H sol  H A 

gs sol sol DGdis at ðxAB Þ  DGat ðxAgs Þ sol xsol AB  xAgs

:

ð22Þ

The non-linear problem described by Eqs. (20)–(22) can be linearized within the low-solubility limit (labeled by ‘‘ls”): ( sol dis;V eq sol xAB ’ 0; DGdis at ðxAB Þ ’ DGat ðxV Þ; ls : ð23Þ gs 0 sol DGgs xsol Ags ’ xAgs ; at ðxAgs Þ ’ DH at ; where DH gs at is the formation enthalpy of the ground state ‘‘gs”. In the low-solubility limit, H sol becomes: dis;V eq 0 H lssol ¼ H ls;nv sol þ DGat ðxV Þ=xAgs ;

ð24Þ

where H ls;nv is the non-vacancy contribution (labeled as sol ‘‘nv”): gs 0 H ls;nv sol ¼ H A  DH at =xAgs :

ð25Þ

From the equilibrium vacancy concentration, Eq. (16), the exponential part of Eqs. (20) and (21) becomes:   expðH lssol =k B T Þ ¼ exp H ls;nv sol =k B T   1=x0 Ags HV  1  exp  ; ð26Þ kBT where the two contributions, non-vacancy and vacancy, are factorized. The last expression indicates that the presence of vacancies increases the solubility by decreasing the number of host B atoms in the solution. 3. Interpretation of H ls;nv sol The calculation of low-solubility of non-interacting defects can be implemented through the first-principles calculation of H lssol (Eq. (24)). It requires the knowledge of the 0 enthalpy DH gs at and composition xAgs of the ground-state ‘‘gs”, as well as the ‘‘true” defect formation enthalpies H A and H V (Eq. (8)). To capture the physical meaning of H ls;nv sol , let us consider an arbitrary dilute solution ‘‘sc” at composition xsc A (see Fig. 2). The label ‘‘sc” indicates that the solution is generated as a supercell of the B solid upon insertion of defects randomly distributed. The first part of Eq. (25) can be rewritten (according to Eq. (9)) as: HA ¼

sc @DH sc at ðxA Þ : sc @xA

ð27Þ

The second part of Eq. (25) becomes: ðBþgsÞ

ðBþgsÞ

DH gs DH at ðxsc @DH at ðxsc at AÞ AÞ ¼ ¼ ; sc sc 0 xA @xA xAgs ðBþgsÞ

ð28Þ

ðxsc where DH at A Þ is the formation enthalpy of the mixture of the pure B solid with the ground state ‘‘gs” with overall composition xsc A (point (B + ‘‘gs”) in Fig. 2). Thus, we obtain:

Fig. 2. The enthalpy–concentration graph demonstrating the low-solubility formation enthalpy. A pure, B pure and ‘‘gs” are the ground states forming a convex hull. ‘‘sc” is a compound (supercell) with some intermediate enthalpy and composition. Compounds (B + A) and (B + gs) correspond to appropriate phase mixtures of the same general composition xsc A.

H ls;nv sol ¼

h i ðBþgsÞ sc sc @ DH sc ðxA Þ at ðxA Þ  DH at @xsc A

:

ð29Þ

The comparison of Eq. (27) with Eq. (29) leads to the conclusion that both H A and H ls;nv sol are derivatives of supercell formation energies with respect to A-composition. For H A , the supercell formation enthalpy is determined with respect to pure A and B phases (the distance between points ‘‘sc” and (B + A) in Fig. 2). For H ls;nv sol , the supercell formation enthalpy is determined with respect to B pure and the ground state ‘‘gs” (the distance between points ‘‘sc” and (B + gs) in Fig. 2). To conclude, H A and H ls;nv sol are characterized by the angles a and b between the B–‘‘sc”/B–A and B–‘‘sc”/B–‘‘gs” lines, respectively. In analogy with the H A ls;nv definition in Eq. (8), H sol can be considered as the ‘‘lowsolubility formation enthalpy”. The quantities H A and H ls;nv sol are identical only for phaseseparating alloys without intermediate ground states (i.e. when ‘‘gs”  A). In this case, the low-solubility can be fordis;A (Eq. (14)) with mally determined by minimizing DGat respect to xA (e.g. Ref. [32]). However, in the general case, the existence of ordered ground states must be verified so the appropriate formalism is used. Generally, in ordering ls;nv differ and could even have different alloys, H A and H sol signs. Note that if the ‘‘sc” point is below (B + gs) in Fig. 2 ðH ls;nv sol < 0Þ, the solubility expressions (20) and (21) are not valid and there must exist an undetected ground state (it might be ‘‘sc” itself) with concentration lower than x0Ags . In this case, such undetected ground state should be used for the calculation of solubility [33]. The expression for non-binary low-solubility within the regular solution model derived in Refs. [23,24] coincides with our derivation in the case of binary alloys without vacancies and high-temperature contributions. This is because the regular solution model corresponds to our model for the free energy in case of dilute solution.

R.V. Chepulskii, S. Curtarolo / Acta Materialia 57 (2009) 5314–5323

0.01 0 FCC

HCP

−0.01 −0.02 eV/atom

As an example of our formalism application, we calculate the solubility of a set of metals in titanium. First, we explore the possible ground states of the A–Ti systems (A = Ag, Au, Cd, Co, Cr, Ir, W, Zn) and then we apply the construction described in the previous section. The low-temperature stability of A–Ti is studied by using our high-throughput quantum calculations framework [5,4,12,14] based on first-principles energies obtained with the VASP software [34]. We use projector-augmented waves (PAW) pseudopotentials [35] and exchange-correlation functionals as parameterized by Perdew–Burke–Ernzerhof [36] for the generalized gradient approximation (GGA). Simulations are carried out without spin polarization except Ti–Co, at zero temperature and pressure [37], and without zero-point motion. All structures are fully relaxed (shape and volume of the cell and internal positions of the atoms). The effect of lattice vibrations is omitted. Numerical convergence to within about 1 meV atom1 is ensured by enforcing a high energy cut-off (357 eV) and dense 6000 k-point meshes. For each system, the energies of 200 crystal structures were calculated. In addition to the 176 configurations described in Ref. [5], the structures included all the symmetrically distinct hexagonal-close packed (hcp), body-centered cubic (bcc) and face-centered cubic (fcc) based superstructures with up to four atoms per cell, and the prototypes A5, A6, A7, A8, A9, A11, B20, C36, D519 ; Al2 Zr4 ; Al3 Zr2 ; Ca7 Ge, CdTi, Cu3 Ti2 ; Ga2 Hf, Ga4 Ni, Ga3 Pt5 ; Ga4 Ti5 ; Hg2 Pt, ITl, InTh, LiB-MS1/2 [12,38], NbNi8 ; NiTi2 SeTl, and V4 Zn5 . The prototypes were considered because they are common or related to Ti alloys [40,39]. This protocol gives reasonable results. In Ref. [5], it was shown that the probability of reproducing the correct ground state, if well defined and not ambiguous, is gc  96:7% (‘‘reliability of the method”, Eq. (3)). We did not consider hcp superstructures with more than four atoms per cell due to the fact that their number increases considerably and that many hcp–hcp systems do not form hcp fields but instead fcc or bcc (i.e. Co–Ti and Ti–Zn in Figs. 6–10, or off-lattice (i.e. Hf–Ti, Hf-Zr [8]). Therefore, even if it is impossible to rule out the existence of an undetected ground state, the protocol is expected to give a reasonable balance between highthroughput speed and scientific accuracy. The whole process is performed in an automatic fashion through the software AFLOW, which generates the prototypes, optimizes the parameters, performs the calculations, corrects possible errors and calculates the phase diagrams [5,18]. Our results of ground state calculation are presented in Figs. 3–10 and Table 1. The correspondence of ab initio vs. experimental results is very good and typical for this type of calculation [5]. We propose the degenerate results for: (1) Co2 Ti: experimentally reported as C15, but with ab initio formation

energies of 311.4 and 304:0 meV atom1 for C14 and C15, respectively; (2) CoTi2 : experimentally reported as NiTi2 , but with ab initio formation energies of 291.2, 287.3 and 285:7 meV atom1 for C37, CuZr2 , and NiTi2 , respectively; (3) TiZn: experimentally reported as B2, but with ab initio formation energies of 195.4 and 193:4 meV atom1 for L10 and B2, respectively. We propose novel results for: (1) Ir7 Ti: experimentally reported as a two-phase region above 500 °C, but with a pos-

−0.03 −0.04 −0.05 −0.06 C11

−0.07

b

B11

−0.08 0 Ag

20

40

60

80

Atomic Percent Titanium

100 Ti

Fig. 3. Ag–Ti ground state convex hull.

0 FCC

HCP

−0.1 eV/atom

4. Ground states of selected titanium alloys

−0.2

−0.3

D1

a

A15

−0.4 C11b

−0.5 0 Au

20

Cu4Ti3

B11

40 60 Atomic Percent Titanium

80

100 Ti

Fig. 4. Au–Ti ground state convex hull.

0.01 HCP

0 HCP

−0.01 −0.02 eV/atom

5318

−0.03 −0.04 −0.05 −0.06 B11

−0.07

C11b

−0.08 0 Cd

20

40 60 80 Atomic Percent Titanium

Fig. 5. Cd–Ti ground state convex hull.

100 Ti

R.V. Chepulskii, S. Curtarolo / Acta Materialia 57 (2009) 5314–5323

5319

0 HCP

0

HCP

−0.05

BCC

HCP

−0.1

−0.2

eV/atom

eV/atom

−0.02 −0.15

−0.25 L12

−0.3

C37/CuZr /NiTi 2

−0.04 −0.06

2

C14/C36

−0.35

−0.08

−0.4

BCC [211] AB2

B2

−0.45 0 Co

20

40 60 80 Atomic Percent Titanium

100 Ti

0

20

Ti

40

60

D1

a

80

Atomic Percent Tungstem

Fig. 6. Co–Ti ground state convex hull.

100 W

Fig. 9. Ti–W ground state convex hull.

0.02 0

0 BCC

HCP

HCP

HCP

−0.02 eV/atom

eV/atom

−0.05 −0.04 −0.06 −0.08

−0.1 A15

−0.15 −0.1 CuZr2

−0.12

−0.2

L10/B2

C15

−0.14 0 Cr

20

40 60 80 Atomic Percent Titanium

100 Ti

Fig. 7. Cr–Ti ground state convex hull.

HCP

−0.1 −0.2 eV/atom

−0.3 −0.4

Ca7Ge

−0.5 −0.6

A15

−0.7 −0.8

0

b

2

−0.9

Ir

C11

L1

L1

0

20

40 60 80 Atomic Percent Titanium

20

40 60 Atomic Percent Zinc

80

100 Zn

Fig. 10. Ti–Zn ground state convex hull.

0

FCC

0 Ti

C14 L12

100 Ti

Fig. 8. Ir–Ti ground state convex hull.

sible ab initio low-temperature ground state Ca7 Ge with formation energy of 373:8 meV atom1 . (2) Ir2 Ti: experimentally reported as a two-phase region above 500 °C, but with a possible ab initio low-temperature ground state C11b , with formation energy of 716:0 meV atom1 . (3) TiW2: experimentally reported as a two-phase region above 500 °C, but with a possible ab initio low-temperature ground state ½211 BCCAB2 , with formation energy of 82:7 meV atom1 (see the notation for prototype in Ref. [5]). (4) TiW4 : experimen-

tally reported as a two-phase region above 500 °C, but with a possible ab initio low-temperature ground state D1a , with formation energy of 83:8 meV atom1 . (5) Ti3 Zn: experimentally not explored, but with a possible ab initio low-temperature ground state A15, with formation energy of 120:0 meV atom1 . In particular, the results indicate that in the Ir–Ti system the low-temperature Ir-rich part of the known phase diagram is not complete, and that the Ti–W alloy has an ordering tendency at low temperature, in contrast to common belief [39,40]. In addition, in Ti–Zn there must exist a Ti-rich compound with Ti composition higher than the reported Ti2 Zn [39,40]. 5. Results and discussion of solubility in Ti alloys The ‘‘raw” formation enthalpies of a substitutional atom or a vacancy in Ti were obtained through Eq. (3) considering a 3  3  3 supercell of hcp Ti. The supercell’s dimensions were chosen to limit the defect–defect interactions by having their distance at least three times larger than the nearest neighbor Ti–Ti bond (see also Fig. 13 and the corresponding discussion below). The results are presented in Table 2. The temperature dependency of solubilities is presented in Fig. 11, while Fig. 12 illustrates a

5320

R.V. Chepulskii, S. Curtarolo / Acta Materialia 57 (2009) 5314–5323

comparison of experimental and theoretical data at T = 700 °C. Table 2 shows negative ‘‘true” formation enthalpies H A for A = Zn, Au, and Ir. The data indicates that calculation of Zn, Au and Ir solubilities in Ti is not feasible without considering nearby intermetallic ground states as performed in our model. Figs. 11 and 12 show that the highest theoretical solubility in Ti occurs for Zn, Cd, Ag and Au as a consequence of their low solubility formation enthalpies (see Table 2). The high solubility has also been observed experimentally for Cd, Ag and Au, whereas, to the best of our knowledge, the solubility of Zn has not been studied. The high solubility of late transition metals Zn, Cd, Ag and Au in the early transition metal Ti is due to the substantial localized stability provided by the filling tendency of the d-

band of Ti [42]. Because of the very high formation enthalpy of vacancy in Ti (reported in Table 2), the vacancy equilibrium concentration was found to be very 6 at low at all considered temperatures (e.g. xeq V < 10 T < 1300 C). Thus, the effect of vacancies on the solubility in Ti alloys is negligible. Although theoretical and experimental results follow similar trends, theoretical solubilities are considerably smaller for most of the considered alloys. A similar discrepancy was also observed for Al–Zr in Ref. [24]. The discrepancy could be due to shortcomings of theory and/or experiment. The main approximations of our model consist of (i) neglecting the defect interaction, (ii) neglecting the spatial defect correlation and (iii) assuming a low concentration of defects. Since these approximations are somehow

Table 1 Low-temperature phase comparison chart for Ag–Ti, Au–Ti, Cd–Ti, Co–Ti, Cr–Ti, Ir–Ti, Ti–V and Ti–Zn. Experimental data correspond to the lowest available temperature [5]. 1 DEgs at ðmeV atom Þ

Space group [41]

Pearson symbol

67.6 63.3a

P4/nmm I4/mmm

tP4 tI6

283.2 430.4 430.6 429.8 356.1a

I4/m I4/mmm I4/mmm P4/nmm Pm3n

tI10 tI6 tI14 tP4 cP8

B11 C11b

62.0 73.7a

P4/nmm I4/mmm

tP4 tI6

C15

120.3a

Fd3m

cF24

L12 C14 C15 B2 C37 CuZr2 NiTi2

254.8 311.4 304.0 386.4 291.2a 287.3 285.7

 Pm3m P63 /mm Fd3m Pm3m Pnma I4/mmm Fd3m

cP4 hP12 cF24 cP2 oP12 tI6 cF96

Ir–Ti Two-phase region above 500 °C Ir3 Ti–L12 IrTi-d Two-phase region above 500 °C IrTi3 –A15

Ir7 Ti–Ca7 Ge L12 L10 Ir2 Ti–C11b A15/tie

373.8 720.3 851.8 716.0 570.3a

Fm3m Pm3m P4/mmm I4/mmm Pm3n

cF32 cP4 tP4 tI6 cP8

Ti–W Two-phase region above 500 °C Two-phase region above 500 °C

TiW2 –BCCAB2 TiW4 –D1a

P3m1 I4/m

hP3 tI10

Pm3n I4/mmm P4/mmm Pm3m P63 =mm Pm3m

cP8 tI6 tP4 cP2 hP12 cP4

Experimental (Refs. [40,39])

Ab initio

Ag–Ti AgTi–B11 AgTi2 –C11b

B11 C11b

Au–Ti Au4 Ti–D1a Au2 Ti–C11b Two-phase region above 500 °C AuTi–B11 AuTi3 –A15

D1a C11b Au4 Ti3 –Cu4 Ti3 /tie B11 A15

Cd–Ti CdTi–B11 CdTi2 –C11b Cr–Ti Cr2 Ti–C15 Co–Ti Co3 Ti–L12 Co2 Ti–C15 CoTi–B2 CoTi2 –NiTi2

Ti–Zn Non explored Ti2 Zn–CuZr2 TiZn–B2 TiZn2 –C14 TiZn3 –L12 a

½211

Ti3 Zn–A15 CuZr2 L10 B2 C14/tie L12

The energy values used in solubility calculations.

82.7a 83.8 120.0a 158.0 195.4 193.4 197.2 198.0

R.V. Chepulskii, S. Curtarolo / Acta Materialia 57 (2009) 5314–5323 Table 2 ‘‘Raw” formation enthalpies H ATi ; H VTi (from Eqs. (2) and (3) at zero pressure), ‘‘true” formation enthalpies H A ; H V (from Eq. (8) at zero pressure) and low-solubility (non-vacancy) formation enthalpy H ls;nv sol (from Eq. (25) at zero pressure) of substitutional A defects (A = Ag, Au, Cd, Co, Cr) and vacancies (V) in Ti. The elements are ordered from low to high ‘‘low-solubility formation enthalpy”. All quantities are in eV units. A

H ATi

HA

ls;nv H sol

Zn Cd Ag Au W Ir Co Cr H VTi ¼ 10:002

6.405 7.267 5.412 3.882 4.303 2.059 1.138 0.674

0.262 0.24 0.304 0.781 0.714 1.14 0.316 1.024 H V ¼ 2:068

0.218 0.461 0.494 0.643 0.838 1.141 1.19 1.205

5321

In order to validate our model of non-interacting solutes, we estimated the convergence of the ‘‘raw” formation enthalpies with respect to the size and composition of the supercell (see Fig. 13). The smaller the supercell, the higher the solute composition, the stronger solute– solute interaction, and correspondingly, the larger deviation of the formation enthalpy from the largest supercell ð3  3  3Þ. Fig. 13 illustrates the trend. For a solute compositions up to 16%, there is a 10% maximum deviation in the formation enthalpies, which does not greatly affect the solubilities at high temperatures. Thus, 16% could be considered as a solute composition threshold. Below this composition, our model is valid. If the criterion is violated, a more precise solute interaction and correlation parameterizations are necessary (requiring, for instance, cluster expansion and Monte Carlo simulations for lattice-based alloys). To conclude, our formalism should be considered as an appropriate starting point to determine if more computationally expensive formalisms are needed. For all the alloys we found that at intermediate temperatures the calculated solubilities are small and within our approximation, but smaller than in experiments, as mentioned before. Ti–Zn has very low formation enthalpy (see Table 2) and, correspondingly, the highest theoretical solubility among the considered alloys (still within the 16% threshold). In Ref. [24], the authors added the defect interactions through the cluster expansion method but did not achieve a significant increase in solubility. In fact, our formalism includes a solute–solvent ordering tendency by considering the real intermetallic ground state other then the pure Ti solid. Thus, we conclude that our approximations

Fig. 11. The solubilities of eight considered transition metals in titanium as functions of temperature (calculated through Eqs. (20) and (24)). Panel (b) is a magnified version of panel (a), to visualize the values for Ir, Co and Cr.

related, calculations leading to low solubility values validate the assumptions, self-consistently. The mean-field approximation neglects the interatomic positional correlations but it is expected to work well when the deviation from complete stoichiometric (pure Ti solid in our case) is small—and so is the solubility in our case—see Section 19 in Ref. [43].

Fig. 12. Comparison of experimental [39,40] and theoretical solubilities of eight considered transition metals in Ti at T = 700 °C. From left to right, the elements are ordered from low to high low-solubility formation enthalpy (and correspondingly theoretical solubility). The error bars characterize the scattering of data measured in different experiments. We could not find the experimental data for solubilities in Ti–Zn.

5322

R.V. Chepulskii, S. Curtarolo / Acta Materialia 57 (2009) 5314–5323

Fig. 13. Relative ‘‘raw” formation enthalpies H ATi ; H VTi (from Eqs. (2) and (3) at zero pressure) as functions of supercell composition and size. We considered ð3  3  3Þ; ð3  3  2Þ; ð2  2  2Þ; ð2  2  1Þ and ð1  1  1Þ supercells of hcp Ti. Enthalpies are presented relative to the largest supercell ð3  3  3Þ, considered as reference. The area between horizontal lines corresponds to 10% error.

(or those of the cluster expansion method) are not responsible for the theory–experiment solubility discrepancy. The other approximation is the assumed independence of the parameters with respect to temperature. The dependence can be caused by atom vibrations, anharmonicity and magnetic ordering. However, theory–experiment solubility discrepancy is observed even at low temperatures (e.g. T 6 HD ðTiÞ ¼ 374–385 K [44]), where the vibrational contribution to the free energy is not important. In fact, vibrational entropy is smaller than the configurational contribution [45], so its inclusion into theory cannot modify the solubility results much. Furthermore, even in Ref. [24], the authors added phonon contribution without increasing the solubility appreciably. Furthermore, in the case of dilute magnetic impurities (such as considered Co), the effect of magnetic ordering on solubility should also be small as the exchange interactions are short-ranged. On the experimental side, the equilibrium solubility tends to be overestimated. In fact, the formation of metastable and/or unstable states, which are subsequently frozen at low temperatures, can make solubility measurements very challenging. In such scenarios, the measured solubility could correspond to spinodal composition instead of the binodal concentration, or simply characterize the undercooled equilibrium solution remaining from the initial specimen prepared at higher temperature. In addition, the segregation of defects into grain boundaries, especially in multicrystalline samples prepared through non-optimal cooling, can dramatically affect the amount of frozen defects. 6. Conclusions Based on the statistical-thermodynamic theory of dilute lattice gas, we have developed an approach to calculate the

atomic bulk solubility in alloys. The approach considers all available binary ground states rather than just pure solids. It is shown that low-solubility follows simple Arrheniustype temperature dependence determined by a ‘‘low-solubility formation enthalpy”. This quantity is defined as the composition derivative of the compound formation enthalpy with respect to nearby ground states. ‘‘Low-solubility formation enthalpy” coincides with the usual defect formation enthalpy only in the case of a phase-separating alloy without intermediate ground states and vacancies. The key quantities of our model can be directly obtained from first principles or by fitting the experimental temperature dependence of solubility. Generalization of our model to intermediate phases and/or to multicomponent, multi-sublattice, interstitial-substitutional alloys is straightforward and it has been performed for the Mg–B–X systems [33]. As example, we applied the developed framework for a set of eight Ti alloys A–Ti (A = Ag, Au, Cd, Co, Cr, Ir, W, Zn). We have found that the highest solubility for Zn, Cd, Ag and Au is in qualitative agreement with available experimental data and band structure considerations. The quantitative differences between the theory and experiment observed in the present and other similar studies are discussed. The formalism is expected to be valid in the limit of lowsolubility, and should be used to estimate if more computationally expensive formalisms are otherwise needed. Acknowledgements We acknowledge Wahyu Setyawan, Mike Mehl, Gus Hart and Ohad Levy for fruitful discussions. This research was supported by ONR (N00014-07-1-0878, N00014-07-11085, N00014-09-1-0921) and NSF (DMR-0639822). We thank the Teragrid Partnership (Texas Advanced Computing Center, TACC) for computational support (MCA07S005). References [1] Ceder G, Chiang YM, Sadoway DR, Aydinol MK, Jang YI, Huang B. Nature 1998;392:694. [2] Johannesson GH, Bligaard T, Ruban AV, Skriver HL, Jacobsen KW, Nørskov JK. Phys Rev Lett 2002;88:255506. [3] Stucke DP, Crespi VH. Nano Lett 2003;3:1183. [4] (a) Curtarolo S, Morgan D, Persson K, Rodgers J, Ceder G. Phys. Rev. Lett. 2003;91:135503; (b) Morgan D, Ceder G, Curtarolo S. Meas Sci Technol 2005;16: 296. [5] Curtarolo S, Morgan D, Ceder G. Calphad 2005;29:163. [6] Fischer C, Tibbetts K, Morgan D, Ceder G. Nat Mater 2006;5:641. [7] Levy O, Chepulskii RV, Hart GLW, Curtarolo S. New face of Rhodium alloys: revealing ordered structures from first principles, submitted for publication. [8] Levy O, Hart GLW, Curtarolo S. Chasing exotic compounds: the necessary synergy of cluster expansion and high-throughput methods, submitted for publication. [9] Hart GL, Blum V, Walorski MJ, Zunger A. Nat Mater 2005;4:391.

R.V. Chepulskii, S. Curtarolo / Acta Materialia 57 (2009) 5314–5323 [10] Bligaard T, Johannesson GH, Ruban AV, Skriver HL, Jacobsen KW, Nørskov JK. Appl Phys Lett 2003;83:4527. [11] Persson M, Bligaard T, Kustov A, Larsen KE, Greeley J, Johannessen T, et al. J Catal 2006;239:501. [12] Kolmogorov AN, Curtarolo S. Phys Rev B 2006;73:180501. R. [13] Calandra M, Kolmogorov AN, Curtarolo S. Phys Rev B 2007;75:144506. [14] Kolmogorov AN, Calandra M, Curtarolo S. Phys Rev B 2008;78:094520. [15] Wolverton C, Siegel DJ, Akbarzadeh AR, Ozolin ß sˇ V. J Phys: Condens Matter 2008;20:064228. [16] Siegel DJ, Wolverton C, Ozolin ß sˇ V. Phys Rev B 2007;76:134102. [17] Akbarzadeh AR, Ozolin ß sˇ V, Wolverton C. Adv Mater 2007;19:3233. [18] Curtarolo S, Hart GLW, Setyawan W, Chepulskii RV, Levy O, Morgan D. Aflow: a computational software to calculate properties of materials in a high-throughput fashion. . [19] Cava RJ, Zandbergen HW, Inumaru K. Physica C 2003;385:8. [20] Singh RK, Shen Y, Gandikota R, Wright D, Carvalho C, Rowell JM, et al. Supercond Sci Technol 2008;21:025012. [21] Harutyunyan AR, Awasthi N, Mora E, Tokune T, Jiang A, Setyawan W, et al. Phys Rev Lett 2008;100:195502. [22] Curtarolo S, Awasthi N, Setyawan W, Jiang A, Bolton K, Harutyunyan AR. Phys Rev B 2008;78:054105. [23] Sigli C, Maenner L, Sztur C, Shahni R. In: Sato T, et al., editors. Proceedings of the 6th international conference on aluminum alloys (ICAA-6), vol. 1. Japan Institute of Light Metals; 1998. p. 87. [24] Clouet E, Sanchez JM, Sigli C. Phys Rev B 2002;65:094105. [25] de Fontaine D. In: Ehrenreich H, Turnbull D, editors. Solid state physics, vol. 47. New York: Academic Press; 1994. p. 33–176.

[26] [27] [28] [29] [30] [31] [32] [33] [34] [35] [36] [37] [38] [39]

[40] [41] [42] [43] [44] [45]

5323

Ceder G. Comput Mater Sci 1993;1:144. Xu Q, Van der Ven A. Phys Rev B 2007;76:064207. Van der Ven A, Ceder G. Phys Rev B 2005;71:054102. Wagner C, Schottky W. Z Phys Chem B 1930;11:163. Mishin Y, Farkas D. Philos Mag A 1997;75:169; Mishin Y, Herzig C. Acta Mater 2000;48:589. Korzhavyi PA, Ruban AV, Lozovoi AY, Vekilov YK, Abrikosov IA, Johansson B. Phys Rev B 2000;61:6003. Matysina ZA, Milyan MI, Zaginaichenko SY. J Phys Chem Solids 1988;49:737. Chepulskii RV, Curtarolo S. Phys Rev B 2009;79:134203. Kresse G, Hafner J. Phys Rev B 1993;47:558. Blochl PE. Phys Rev B 1994;50:17953. Perdew JP, Burke K, Ernzerhof M. Phys Rev Lett 1996;77: 3865. So that by calling the enthalpy we actually mean the energy in this section. Kolmogorov AN, Curtarolo S. Phys Rev B 2006;74:224507. Villars P, Cenzual K, Daams JLC, Hulliger F, Massalski TB, Okamoto H, et al. Inorganic materials database and design system. Binaries Edition. Metals Park, OH: ASM International; 2003. Massalski TB, editor. Binary alloy phase diagrams. Metals Park, OH: ASM International; 1992. Hahn T, editor. International tables for crystallography, vol. A. Dordrecht: D. Reidel; 1983. Friedel J. In: Ziman JM, editor. Physics of metals. Cambridge: Cambridge University Press; 1969. p. 495. Krivoglaz MA, Smirnov AA. The theory of order–disorder in alloys. London: Macdonald; 1964. Chen Q, Sundman B. Acta Mater 2001;49:947. van de Walle A, Ceder G. Rev Mod Phys 2002;74:11.