Configurational thermodynamics of a 1/2〈111〉 screw dislocation core in Mo–W solid solutions using cluster expansion

In this work, we have developed a methodology to obtain an ab initio cluster expansion of a system containing a dislocation and studied the effect of configurational disorder on the 1/2 〈 111 〉 screw dislocation core structure in disordered Mo 1 − x W x alloys. Dislocation cores control the selection of glide planes, cross slip, and dislocation nucleation. Configurational disorders in alloys can impact the dislocation core structure and affect dislocation mobility. For our calculations, we have used a quadrupolar periodic array of screw dislocation dipoles and obtained the relaxed structures and energies using density functional theory. We have obtained the dislocation core structure as a function of composition and the interaction energies of solutes with the dislocation as a function of position with respect to the core. With these energies, we performed mean-field calculations to assess segregation toward the core. Finally, with the calculated energies of 1848 alloy configurations with different compositions, we performed a first principle cluster expansion of the configurational energetics of Mo 1 − x W x solid solutions containing dislocations.


I. INTRODUCTION
Alloying has played a key role in producing useful materials from pure metals, enabling the tuning of their strength by the introduction of solutes.One of such cases is refractory materials, which are of technological interest due to their high melting point, but they exhibit room temperature brittleness, limiting their resistance to mechanical load.This undesired behavior can be alleviated by alloying.The case of W-Re alloys is a well-known example, in which the addition of Re produces a ductilizing effect. 1,2The mechanisms behind this ductilization have been widely studied but without a clear answer.4][5][6] Re impurities change the dislocation core structure from a symmetric core to an asymmetric core, 4,[6][7][8][9][10] changing the preference of the slip plane from {110} to {112}, which increases the number of available planes from six to 12. On the other hand, in other materials systems, alloying is a typical method to achieve hardening.This underlines the importance of developing simulation tools to treat dislocations with the accuracy needed for alloy problems.
The study of the interaction of point defects with dislocation cores is well suited for density functional theory (DFT) calculations, since it is able to assess the structure of the dislocation core and how it changes with the incorporation of solutes.The majority of first principles studies on dislocation solute interaction have focused on the interaction of a single impurity atom with the dislocation core, 6,[11][12][13] and a limited number of first principles studies on the change of dislocation core properties with varying solute concentration.These higher solute concentrations have been tackled with the virtual crystal approximation (VCA) 14 and with the supercell approach. 8In the VCA, atoms are substituted by one effective type of atom of intermediate nuclear and electronic charge, and in the supercell approach, atoms are randomly substituted in the cell by alloying atoms.The VCA method severely approximates the real electronic and atomic structure and does not consider the local environment effects.Local environments can modify the core structure in different ways along the dislocation line, giving rise to different mobilities for different dislocation line segments enhancing the formation of dislocation jogs, which would be lost in the average picture given by the VCA.On the other hand, the supercell method requires a large number of simulations to account for all different configurations.The dislocation core structure as a function of alloy composition for W-Re, 8 Mo-Ta, 7 Fe-Co 9 solid solutions was investigated by Romaner et al., while the Peierls stress and the interaction of the core with transition metal solutes in W were studied by Samolyuk et al. 10 Medvedeva used the atomic row model with ab initio parametrization of the interatomic interaction to study solid solution softening and hardening in Nb, Ta, Mo, and W. 4 Even though first principles allow us to characterize the core and obtain the Peierls stress, direct estimation of mechanical properties and temperature effects is not attainable, i.e., length scales needed to incorporate effects of extended defects, their mobility, and interactions.The cluster expansion (CE) method has been used to parametrize the configurational dependence of physical properties.The CE method has been employed to study the equilibrium concentration of point defects, [15][16][17] but to our knowledge, it has not been used to study the alloy configurations close to dislocations.The CE permits one to study the effects of alloying on the dislocation properties for different alloy compositions and to obtain dislocation core concentration profiles as a function of temperature, allowing one to analyze if substitutional atoms segregate to the dislocation core, and compare with the ordering tendencies of the bulk.
The aim of this work is to develop a methodology to allow for an ab initio based cluster expansion of the configurational energetics of a system containing screw dislocations.We have chosen the Mo-W system as our model case since it forms a complete solid solution for all the compositional range. 18We have investigated the dislocation core structure and the interaction energies as a function of distance to the dislocation core for different alloy compositions.We used the interaction energy data to perform mean-field calculations to estimate the concentration profile as a function of temperature.In Sec.IV, we present the details of the CE in a Mo-W dislocation system.

A. Cluster expansion
In the cluster expansion (CE) method, 19 the energy of a crystalline solid is a function of the atomic arrangement of the lattice.The atomic configuration of a binary alloy A 1−x B x can be characterized by assigning to each of the N lattice sites of the system an occupation variable, S i (i ¼ 1, 2, . . ., N) with S i ¼ À1 or S i ¼ þ1 depending if the site i is occupied by an atom type A or B. In this way, any configuration σ can be expressed by the vector of the N site occupational variables.The energies of all the possible configurations can be expressed using a generalized Ising Hamiltonian where V's are the effective cluster interactions (ECIs) for each cluster order.Although this expression is complete, at least for a finite system, in practice, it must be truncated.The cluster expansion method has to be bound by approximations due to the complex nature of real systems, such as non-linear concentration dependence of cluster interactions that has been observed in alloy surfaces. 20In the case of W-Mo, we do not expect any non-linear concentration dependence of the cluster interactions that could not be described by concentration independent cluster expansion.In general, bonding is given by short length scales, and one must expect that long-range interactions would have negligible ECI energies, as such; a finite number of interaction parameters can yield energetics with sufficient accuracy.ECIs can be determined from the energies of a set of configurations through a series of methods but the most straight forward is to follow the approach of Connolly and Williams 21 with a least-square fitting.In this work, we have used the MIT ab initio phase stability (MAPS) 22 code, as implemented in the alloy-theoretic automated toolkit (ATAT), 23 to truncate the expansion and determine ECIs.Out-of-the-box cluster expansion packages struggle to perform expansions of large systems.
The configurational space of the smallest cell in this study is too big to be handled directly by these packages, which often require a complete enumeration of the possible structures.We have opted to generate random configurations for each of the different alloy compositions.For our CE, we have included the energies of pure Mo, pure W, the dilute limit cases, and Mo 1−x W x with x = 0.25, 0.50, and 0.75.For the dilute cases, we have substituted each of the nonequivalent atomic positions by the atom of the opposite species.
Using the symmetry package of Phonopy, 24 we have found 40 unique atomic positions in the dislocation system.For x = 0.25, 0.50, and 0.75, we have used different random configurations and also used the substitution approach, having a total of 1848 unique ab initio energies to use in the CE.

B. Density functional theory (DFT) calculations
The energy calculations required for the construction of the cluster expansion were performed with the Vienna ab initio simulation package (VASP) 25,26 with Perdew-Burke-Ernzerhof generalized gradient approximation to model the exchange correlation effects 27 and the projector augmented wave basis 28 with an energy cutoff of 400 eV, and a 1 × 2 × 16 k-point sampling grid.The structures were relaxed until forces were smaller than 0.01 eV/Å.0][31] The dislocation is placed in the center of gravity of three neighboring 〈111〉 atomic columns.The dislocation dipole is introduced into the simulation cell by applying the displacement field of each dislocation as given by anisotropic linear elasticity as implemented in Babel software. 32Then, the atomic positions are relaxed to minimize the energy of the simulation box.
The dislocation cells created using the Babel software require the lattice and elastic constants as an input.To obtain these parameters for Mo 1−x W x alloy compositions with x = 0.25, 0.50, and 0.75, we used the special quasirandom structure (SQS) approach 33 using 4 × 4 × 4 cubic supercells.The obtained values of the lattice and elastic constants for different alloy compositions are plotted in Fig. 1.

Journal of Applied Physics
The lattice constant closely follows Vegard's law, similar to the reported experimental results, 34 a behavior that is also followed by the elastic constants.
6][37] For most of the calculations, we employed a cell with 1 Burgers vector (b) thickness.This slice comprises three atomic layers and contains 135 atoms.The choice of a cell with 1b along the dislocation line direction is a compromise on the computational cost.This imposes a self-interaction with the periodic image of each site.In the dilute case, for example, we are in fact introducing a column of solute atoms.In order to assess the impact of this approximation, we investigated the intraplanar interaction of solutes, by including calculations of cells with 2 and 3 Burgers vectors in the dilute limit case.

III. RESULTS
A. Core structure Atomistic simulations have revealed two types of core configurations in the 1/2〈111〉 screw dislocation distinguished by the symmetry of their spreading into {110} planes of the 〈111〉 zone.The first case, the symmetric or non-degenerate core, the spreading is evenly distributed, whereas in the asymmetric or degenerate core, the spread is not centered on the dislocation core, and it exists in two energetically equivalent symmetric variants.The dislocation core structure can be visualized with the differential displacement maps. 38In these maps, the crystal is projected to the plane perpendicular to the dislocation line and the difference between the vector connecting two neighboring atoms in the relaxed dislocated crystal and the same vector in the perfect crystal is plotted as arrows pointing from one atom to another.These arrows are centered in the midpoint between the atomic positions and have been normalized to b/3.
The nature of the core can be distinguished by the polarization parameter 39,40 given in Eq. ( 2), where d αβ is the relative displacement of neighboring atoms α and β in the [111] direction.A polarization p ¼ 0 corresponds to a symmetric core configuration, and a value of p ¼ j1j to an asymmetric core.The polarization allows a continuous description from a symmetric non-degenerate core to the asymmetric degenerate core configuration.The labels of the positions involved in the polarization expression are shown in Fig. 4, We found no changes in the dislocation core structure as a function of composition.The symmetric core is maintained for all alloy compositions, as shown in Fig. 2.This behavior is expected since both pure metals possess the same core structure, and it is similar to the results obtained by Li et al. for W-Ta 7 and for W-Zr. 10 It is important to note that, in fact, alloying can produce changes in the polarization of the core, as it is the case for the transition from a symmetric to an asymmetric core in W-Re at ∼25% Re 8 and for Fe-Co at ∼50% Co. 9

B. Dislocation energetics
The interaction energy per solute atom between the dislocation and the solute atoms is given by In this expression E disloþsolute is the energy of cell containing a dislocation and a solute atom (E dislo without solute) and E bulkþsolute the energy of a perfect bcc cell in which one solute atom has been introduced (E bulk without solute).Figure 3(a) shows the interaction energies for the dilute case of Mo and W as a function of distance to the dislocation core.Figure 3(b) shows the interaction energy for the positions closest to the dislocation core for all the compositions used in this study.For all cases, there exists an attractive interaction between Mo and the dislocation, preferring to occupy the closest position to the dislocation core.The opposite occurs for W, where it would prefer to sit as far away as possible.The interaction energy

Journal of Applied Physics
value we obtained for dilute Mo in W of −0.05 eV is in close agreement to the one reported by Hu et al. 41 As mentioned before, due to periodic boundary conditions and cell shape, substituting an atom within the cell in fact introduces a column of solutes in the system.In the dilute case, it is rather the interaction of a dislocation with a column of atoms than with an individual solute atom.In the case of higher solute concentrations, we have a periodic repetition of the distribution of atoms along the dislocation line direction.To investigate this size effect, larger cells are needed.2] We have investigated the effect of the cell size using larger systems along the dislocation line direction, with cells with 2 and 3 Burgers vector thickness.The interaction energies for the different sized systems are shown in Table I.We found small values for the solute repulsion energies, given in Table II.These values are small, particularly, when compared with the repulsion between interstitial C atoms in bcc Fe, which have a reported value of 0.21 eV. 12 This is an indicator that there exists little interplanar solute interaction and cells with 1b should be large enough to describe this system.Besides the effect of periodicity, short cells are not able to capture the effect of the local environment on the dislocation.There exists a configurational force acting on the dislocation line due to the local spatial fluctuations in the solute concentrations. 43If the local environment adjacent to a segment of the dislocation line is more favorable, there exists a driving force for that segment to move.In an attempt to assess this effect on the Mo-W system, we relaxed systems with 4b and 5b cell thicknesses for the 50/50 concentration case.To examine the dislocation line behavior, we used the dislocation analysis tool implemented in OVITO 44,45 and found no change in the position of the dislocation lines with the employed cell sizes.

C. Mean-field segregation simulation
We have performed mean-field analysis 46,47 to study the segregation tendencies as a function of temperature using the methodology of Ventelon et al., 12 which studied carbon interstitials in iron.This approximation accounts for configurational entropy but

Journal of Applied Physics
neglects local chemical environment effects, vibrational, and electronic entropy.We consider only nearest neighbor interactions between solute atoms; the solute interaction energy is given by , where V cc is the first neighbor repulsion energy between solute atoms V cc ¼ E int (1b) À E int (2b).The average solute-dislocation interaction is Minimizing the free energy, we get 12 where C d is the concentration on the dislocation core.The segregation energy is given by The solute concentration is connected to nominal concentration C nom by matter conservation, The number of solute sites in the dislocation core is N d ¼ ρVn/b, and the total number of sites N 0 ¼ 2V/a 3 , where ρ is the dislocation density, V the volume of the system, a the lattice constant, b the Burgers vector, and n the number of sites in the core.Using this relationship, Eqs. ( 5) and ( 6) can be solved self-consistently, we can obtain the concentration profile as a function of temperature.In Fig. 4, we show results for Mo 0.2 W 0.8 , Mo 0.8 W 0.2 , and Mo 0.5 W 0.5 compositions for different atomic positions close to the dislocation core.In the case of a W matrix with Mo solutes, all the considered dislocation core positions are completely decorated with Mo atoms at low temperatures.Increasing temperatures will enhance tungsten to occupy core positions, starting with the farthest core position considered here, until the bulk concentration is also obtained at the dislocation core.In the case of a Mo matrix with W solute atoms, the dislocation core positions are occupied by Mo atoms at low temperatures.As temperature increases, tungsten will gradually stop segregating away from the core, with the closest core position being the last one to attain the bulk concentration.

D. Cluster expansion
From the 1848 ab initio calculated energies, 70% conformed our training set to perform the fit of V ij .We evaluated the reliability of the fit based on the ability to reproduce the energies of the 555 remaining configurations that were not included in the fit.For the dislocation CE, we included on-site and first nearest neighbor pairs.ECIs for the on-site clusters are shown in Fig. 5(a) and pair interactions in Fig. 5(b).In Fig. 5(b), we have also included the pair interactions for a perfect Mo-W system (without a dislocation).In both dislocation and perfect cases, it can be observed that the ECI values are small, as is expected from a material that according to the experimental phase diagram forms a solid solution for the full composition range and without strong ordering tendencies.In the dislocation system, the on-site ECIs are stronger closer to the dislocation core, as is expected from the dislocation influence and, in general, have a larger value than the pair interactions.The cluster expansion is dominated by the on-site clusters.In fact, there are four on-site ECIs that dominate the CE, which correspond to the closest positions to the core.The rest of the on-site ECIs are of comparable magnitude to the pair ECIs.From Fig. 5(b), it can be noted that the presence of the dislocation introduces a distinctive influence in the pair ECI, as seen by comparing the values between ECIs in the bulk and with those in the dislocation system.It shows that changes in pair interactions due to a dislocation core should not be disregarded as possible in future studies of alloys in the presence of dislocations.
To measure the quality of the fit, we used the cross-validation score.The cross-validation score provides a measure of the predictive power of the cluster expansion, it is analogous to the root mean square error but designed to estimate the error of the structures not included in the least squares fit. 22It should be noted that a fit that only includes on-site terms gives a cross-validation score of 0.0003 eV/atom, and the CV decreases when more pairs are

Journal of Applied Physics
included, with the fit that includes all nearest neighbor pairs (148 pairs) the one with the lowest CV score.We have plotted the value of the ECI for different numbers of clusters included in the fit in Fig. 6, where we have included fits with only on-site clusters, 80 pairs, and all nearest neighbor pairs.The value of the ECI seems to converge with the inclusion of more pair clusters in the fit.The self-interaction image of each atomic position is part of the nearest neighbor set in the direction of the Burgers vector, as such the selfinteraction pair is included, and it has a value of 0 as it can be seen in Fig. 6.ECIs with the largest values in Fig. 6 correspond to the positions closest to the dislocation core (Fig. 5).

IV. CONCLUSIONS
In this work, we applied the cluster expansion method to a system containing a screw dislocation.We chose the Mo-W system because it forms a complete solid solution for the whole concentration range and thus, in the absence of ideal bulk effects, would highlight the impact of dislocation on the configurational energetics.First, we obtained the energetics of Mo and W solutes for different alloy compositions.We found that Mo atoms prefer to be closer to the dislocation core, while W prefers to be far away from the dislocation.We investigated the cell size effect on the interaction energies and found weak solute-solute interaction, indicating that a 1b cell size is a reasonable approximation to describe this system.With the interaction energies, we employed a mean-field analysis to estimate segregation profiles as a function of temperature for different global compositions.In the case of equiatomic Mo 0.50 W 0.50 at low temperatures, the dislocation core is fully decorated with Mo, with tungsten segregating first to the positions farthest from the dislocation core, with the closest position attaining a 0.37 W concentration at 1500 K, until reaching the nominal composition at high temperatures.We obtained a CE that properly predicts the calculated DFT energies of the dislocation system for structures from both the training set and the test set.The ECI values for both on-site and pair interactions are small, with the on-site term interactions dominating the ECI values and being stronger close to the dislocation.In a future study, we will use the obtained ECIs to perform Monte Carlo simulations to obtain the concentration profiles of the dislocation core as a function of temperature and compare it with the mean-field approximation performed in this study.This methodology to perform a cluster expansion in a supercell containing a dislocation dipole allows the identification of the dislocation effect on relevant on-site and pair interactions and, in principle, can be applied to more complicated alloy systems, e.g., systems with strong ordering tendencies and magnetic interactions.

FIG. 1 .
FIG. 1.(a) Lattice and (b) elastic constants as a function of composition for Mo 1−x W x alloys.

FIG. 2 .
FIG. 2. Differential displacement maps and polarization values of the 1/2〈111〉 screw dislocation core for different compositions of Mo 1−x W x solid solutions.

TABLE I .
Interaction energies as a function of cell size for the dilute case of Mo and W.

TABLE II .
Repulsion energy between first nearest neighbor solutes.