Interplay of force constants in the lattice dynamics of disordered alloys : An ab-initio study
Rajiv K. Chouhan, Aftab Alam, Subhradip Ghosh, Abhijit Mookerjee
IInterplay of force constants in the lattice dynamics of disordered alloys :An ab-initio study
Rajiv K. Chouhan, , Aftab Alam, † , Subhradip Ghosh, and Abhijit Mookerjee ∗ Department of Condensed Matter Physics and Materials Science,S.N. Bose National Center for Basic Sciences, Kolkata 700098, India Theoretical Sciences Unit, Jawaharlal Nehru Centre for Advanced Scientific Research, Jakkur, Bangalore 560064, India Department of Physics, Indian Institute of Technology, Bombay, Powai, Mumbai 400076, India and Department of Physics, Indian Institute of Technology, Guwahati, Assam, India (Dated: April 14, 2014)A reliable prediction of interatomic force constants in disordered alloys is an outstanding problem.This is due to the need for a proper treatment of multisite (at least pair) correlation within a randomenvironment. The situation becomes even more challenging for systems with large difference inatomic size and mass. We propose a systematic density functional theory (DFT) based study topredict the ab-initio force constants in random alloys. The method is based on a marriage betweenspecial quasirandom structures (SQS) and the augmented space recursion (ASR) to calculate phononspectra, density of states (DOS) etc. bcc TaW and fcc NiPt alloys are considered as the two distincttest cases. Ta-Ta (W-W) bond distance in the alloy is predicted to be smaller (larger) than those inpure Ta (W), which, in turn, yields stiffer (softer) force constants for Ta (W). Pt-Pt force constantsin the alloy, however, are predicted to be softer compared to Ni-Ni, due to large bond distance of theformer. Our calculated force constants, phonon spectra and DOS are compared with experimentsand other theoretical results, wherever available. Correct trend of present results for the two alloyspave a path for future studies in more complex alloy systems.
PACS numbers: 63.20.dk, 63.50.-x, 63.50.Lm
I. INTRODUCTION
In spite of years of rigorous research, a reliable ab-initio theoretical model for structural and/or substitu-tional disordered alloys is still lacking. More importantly,the understanding of the problem of lattice dynamics (orphonons) in disordered alloys is still in its adolescence.This is mainly due to the presence of off-diagonal (mul-tisite) disorder arising out of the force constant tensor D µν ( R i − R j ) in the phonon problem. In addition thesum rule D µν ( R i ) = − (cid:80) R j D µν ( R i − R j ) makes the dis-order at a site depend upon its neighborhood or the socalled environmental disorder. As such first principlesmethod based on e.g. coherent potential approximation(CPA) is inapplicable in this problem. Within CPA,the atoms occupancy are assumed in an average senseembedded in a structure-less uniform average medium.This prohibits CPA to include structural relaxations,which contrasts from the experimental observations be-cause bond distance between atomic pairs (e.g. A-A,B-B, and A-B in a binary alloy) are generally different. Various generalizations of CPA have been suggestedover time, each has its own advantages and disadvan-tages. An striking approach, emerged in recent years, isthe so called special quasirandom structure (SQS) pro-posed by Zunger et al. , which carries the signature ofconfiguration correlation with them. In particular, SQSis an ordered supercell which is constructed in such away to mimic the most relevant pair and multisite cor-relation functions of the disordered phase. Unlike CPAand other related approaches, SQS is a local structuralmodel which captures the most relevant microstates of disordered phase.As far as the calculation of force constants for ran-dom alloys are concerned, three approaches are mainlyutilized in the past. The first attempt was to fit anempirical set of force constants to match the available ex-perimental phonon spectra. The second approach was tocompute the force constants from selected ordered struc-tures and then use them for random alloys. This is ofcourse not a proper solution, because dynamical matri-ces are not directly transferable across the environment. In later studies , few SQS methods have been used, butonly ≤ is a powerful techniqueto do the same, and has been described in great detailsin many of our earlier papers. Interested readers are re-ferred to articles 6 and 7 for any details related to formal-ism. In the present rapid communication, we integratea firs principles SQS method with the ASF to demon-strate the interplay of force constants within a disorderenvironment. Unlike previous approaches, a systematiccalculation of the force constants with increasing size ofthe SQS cell is made. Stress on the atomic sites are di-rectly related to the force constant matrix and hence asmall disturbance leads to a large change. To overcomethis effect, we use the SQS cell in conjunction with the small displacement method to construct the dynami-cal matrix D µν . Based on the predicted bond length a r X i v : . [ c ond - m a t . m t r l - s c i ] A p r TABLE I. Dynamical matrices D µν ( | R | ) (Newton/meter) forbcc Ta W (top) and fcc Ni Pt (bottom). N-atom rep-resents the size of the SQS supercell, Other experimental and theoretical data are give for comparison.bcc Ta W directionTa-Ta 27.707 25.103 22.324 22.872 16.983 111 xx W-W 28.934 24.568 21.063 21.069 23.000 111 xx Ta-W 27.734 26.071 23.192 21.783 23.984 111 xx Ta-Ta 11.504 15.727 16.120 17.599 11.201 111 xy W-W 6.655 8.348 9.418 11.156 19.200 111 xy Ta-W 8.585 13.080 13.668 13.984 17.603 111 xy Ta-Ta -0.009 9.114 25.937 17.180 1.182 200 xx W-W 0.016 17.693 36.233 41.915 47.300 200 xx Ta-W 16.448 13.442 30.733 29.367 24.803 200 xx Ta-Ta -0.001 3.105 -2.909 -0.598 1.423 200 yy W-W 0.008 2.755 -1.648 0.244 -0.800 200 yy Ta-W -1.634 3.656 -2.163 -0.620 1.184 200 yy fcc Ni Pt directionNi-Ni 6.289 8.433 9.813 9.107 8.231 110 xx Pt-Pt 13.755 34.576 36.317 26.747 33.494 110 xx Ni-Pt 12.421 17.387 21.210 17.200 17.868 110 xx Ni-Ni 3.791 8.845 11.115 9.712 9.580 110 xy Pt-Pt 9.008 36.546 43.377 31.426 39.655 110 xy Ni-Pt 8.026 18.287 25.091 19.416 20.740 110 xy Ni-Ni 5.394 0.946 -1.845 -0.0423 -0.525 110 zz Pt-Pt 9.487 2.310 -8.351 -4.199 -6.854 110 zz Ni-Pt 8.785 1.035 -4.124 -1.141 -2.820 110 zz distribution and the calculated force constants for eachpairs A-A, B-B, and A-B, it is concluded that a minimumof 32-atom SQS cell is needed to capture the importantdisorder correlations, and hence a reliable phonon disper-sion. Two different alloy systems, bcc TaW and fcc NiPt,are chosen to demonstrate the reliability of the approach.Both the systems have inelastic neutron scattering datato compare our theoretical results. Ta and W belong to5d-metal series with similar size and atomic masses, butquite different force constants. Ni and Pt on the otherhand differ significantly in size ( ∼ φ Pt-Pt is 55% larger than φ Ni-Ni ). As such considerable differ-ences are expected from the lattice dynamical propertiesof the two systems, with different interplay of force con-stant interactions. We compare our calculated dynamicalmatrices, phonon dispersion and DOS with existing ex-perimental and theoretical data, wherever available.To calculate the force constants within a disorder en-vironment, we first develop a structural model based onthe SQS method. SQS is an N -atom periodic structureconstructed in such a way that the associated set of cor-relation function of this structure mimic the ensembleaverage correlation functions of the random alloy. Dif-ferent sized SQS-cells (8-atom, 16-atom, 32-atom and64-atom) are used for both the fcc and bcc systems(See supplementary material [12] for the SQS structures).We use Vienna ab-initio simulation package (VASP) with a pseudo-potential and a projected-augmented-wave(PAW) basis based on the local density approximation(LDA). The cut-off energy for the electronic wavefunc-tions is 500 e V. All the structures are fully relaxed untilthe energy converges to within 10 − e V and the forces oneach atom is less than 0 . e V/˚A. Such a relaxation, ina way, captures the effect of any static displacements thatmay be present in a real system. A Monkhorst-pack Bril-louin zone (BZ) integration with a 8 k -mesh is used for8-atom SQS calculation. Convergence of D µν as a func-tion of k -points is checked, see below. Smaller k -meshesare used for larger supercells. Magnetic (non-magnetic)calculations are done for NiPt (Ta-W) systems. Relaxedlattice constants for 8-atom, 16-atom, 32-atom and 64-atom SQS calculation for fcc Ni Pt are 3 .
72, 3 .
72, 3 . .
72 ˚A respectively, compared to the experimentalvalue of 3 .
785 ˚A. For bcc Ta W , they are 3 .
23 ˚Afor all the structures, compared to 3 .
23 ˚A as observed. To extract the force constant matrices, we use the fullyrelaxed SQS structures and apply the small displace-ment method using PHON package implemented withinVASP. For 32-atom SQS, force fields are constructed byapplying 48 displacements in fcc Ni Pt and 96 dis-placements in bcc Ta W along 3-cartesian axes, eachof amplitude 0 . δ . For example, a particular atom at (1 / / / ± δ ,1 / ± δ , ± δ ) oratoms at (1 / / /
2) and (1,0,0) in a bcc lattice movesto (1 / ± δ ,1 / ± δ ,1 / ± δ ) and (1 ± δ , ± δ , ± δ ) respec-tively. In order to retrieve the desired symmetry of thedynamical matrix, we apply transformation operation onthese average matrices to get the direction specific dy-namical matrices e.g. φ = B φ − B T along one of D µ ν ( N / m ) xx xy xx yy Ta-Ta(black)W-W(red)Ta-W(green)
FIG. 1. (Color online) Convergence of various componentsof dynamical matrix ( D µν ) as a function of the number of k -points (n, along each direction) for 8-atom SQS bcc Ta W .Different colors indicate different pairs of force constants, andvarious components of D distinguished by different symbols. the nearest neighbor direction of bcc lattice, where B isthe transformation matrix. A list of transformation ma-trices along specific directions for both fcc and bcc lat-tices are given in Table S1 of the supplement (See [12]).The calculated dispersion in the occupancy for eachpairs reflect the sensitivity of the bond distances onthe local environment. The calculated nearest neigh-bor (next nearest neighbor) average bond distances forthree pairs d Ta-Ta , d
W-W , and d
Ta-W are 2 .
837 (3 . .
775 (3 . . . W . The nn-bond distance for Ta-Ta inthe alloy is found to be ∼ .
8% smaller compared tothat in pure Ta, while W-W bond distance in the alloyis ∼ .
7% larger than that in pure W. The calculateddynamical matrices (up to 2nd neighbor) for 8-, 16-, 32-and 64-atom SQS Ta W are shown in the upper partof Table I. Note that, experimental force constants forTa-Ta and W-W pairs are not for the alloy, but for pureTa and pure W respectively. Force constants for Ta-W pair, however, are indeed for the alloy. Notably ourcalculated Ta-Ta force constants in the alloy are stiffercompared to those in pure Ta. On the other hand, thecalculated W-W force constants are softer than those inpure W. This prediction actually jibe with the calculatedbond lengths between these two pairs. Alloying shrinks(expands) the Ta-Ta (W-W) bond lengths making thesprings relatively stiffer (softer). As far as the force con-stants for Ta-W pair goes, 64-atom SQS results are ourbest numbers to compare with the experiment. Experi-mental force constants are computed using a polynomialfit to their measured dispersion. Keeping in mind thesensitivity of the estimates, both on the theoretical andexperimental front, the overall agreement between the64-atom SQS results and the experiment for Ta-W forceconstants is fairly well.The calculated average bond lengths for Ni-Ni, Pt-Pt D µ ν ( N / m ) xx xy zz Ni-Ni(black)Pt-Pt(red)Ni-Pt(green)
FIG. 2. (Color online) Same as Fig. 1, but for fcc Ni Pt . and Ni-Pt pairs are 2 . .
692 and 2 . Pt . Ni-Ni (Pt-Pt) bond length inthe alloy is ∼ .
3% larger ( ∼ .
8% smaller) than thosein pure Ni (Pt). As such, Ni-Ni (Pt-Pt) force constants inthe random alloy is expected to get softer (stiffer) com-pared to those in pure Ni (Pt). The calculated ab-initioforce constants for the three pairs in disordered Ni Pt alloy are shown in the lower panel of Table I. As before,results are shown for 8-, 16-, 32- and 64-atom SQS. Be-cause these calculations are done for ferromagnetic NiPt,a one to one comparison with the experiment requiresobserved data at very low temperature, which we wereunable to find in the literature. We have also performednon-magnetic calculations for NiPt, and found similarresults for D µν . In Table I, the force constants underthe column labeled Other are the results from a recentcalculation on Ni Pt alloy by Granas et al. Theseforce constants for each pair (Ni-Ni, Pt-Pt and Ni-Pt)are within the disordered environment, and agree fairlywell with ours within a few percent. Calculated forceconstants for Ni-Ni (Pt-Pt) pairs in the alloy are foundto be softer (stiffer) compared to those in pure Ni (Pt)(See Ref. 21 for the force constants of pure Ni and Pt).This, again, goes in accordance with the bond lengths ofrespective pairs in the alloy vs. those in pure elements.Keeping in mind the sensitivity of dynamical matri-ces to finer details of calculation, we have checked theconvergence of various components of D with respect tothe number of k -points used in BZ integration. This isshown in Fig. 1 and 2 for 8-atom SQS Ta W andNi Pt respectively. Different colors indicate differentpairs of force constants while various components of D are distinguished by different symbols. One can noticethat ALL the components of D for both the systems arewell converged by 6 k -points.Figure 3 shows the calculated phonon dispersion (left)and the configuration averaged phonon DOS (right) forbcc Ta W . Dispersion curves are calculated usingthe force constants of 32-atom SQS, as listed in Table Γ H P Γ N L T L LT L
DOS (Arbit. scale) T T [ζ00] [ζζζ] [ζζζ] [ζζ0] T ν ( T H z ) bcc Ta W
50 50
FIG. 3. (Color online) (Left) Phonon dispersion for bccTa W alloy using the force constants of 32-atom SQS. Land T stands for longitudinal and transverse modes. Errorbars indicate the calculated FWHM’s. Square symbols in-dicate experimental data. (Right) Phonon DOS using theforce constants of 8-, 16-, and 32-atom SQS. I. Phonon DOS, however, are shown with three sets offorce constants, i.e. 8-, 16- and 32-atom SQS, for compar-ison. Error bars in the dispersion curve indicate the fullwidths at half maxima (FWHM). Our calculated phonondispersion and DOS using 64-atom SQS force constantsare very similar to those using 32-atom SQS ones (alsotrue in case of NiPt, see below). As such it is intuitivelyexpected that, with increasing size of the supercell, theforce constant matrix D converges in a collective man-ner. Finally, note that our calculated dispersion com-pares fairly well with the experiment (shown by squaresymbol).Figure 4 (left) shows the phonon dispersion forNi Pt alloy calculated using the force constants of32-atom SQS cell (see Table I). Unlike TaW, NiPt alloyshows interesting split band behavior along each symme-try direction. This is due to the strong disorder in bothmass and force constants, giving rise to resonant modes,and has been evidenced in previous studies. Errorbars with solid circles indicate the calculated FWHM.Error bars with square symbol along [ ζ
00] direction indi-cate the neutron scattering data. The panel on the rightshows the configuration averaged phonon DOS with threesets of calculated force constants. Square symbols indi-cate the generalized phonon DOS derived from inelasticincoherent scattering. Notice that the calculated bandedge increases with increasing the SQS cell size, and com-pare better with experiment. The integral value undereach phonon DOS, however, remain the same. It is im-portant to emphasize that the experimental phonon DOSis only shown for reference. A one to one comparison be-tween our calculated DOS and the experimental DOS isnot feasible. In inelastic neutron scattering, phonon DOScan be represented as N ( ω ) = (cid:80) j ( b j /M j ) n j ( ω ), where b j , M j and n j ( ω ) are the inelastic scattering cross sec-tion, atomic mass, and the partial phonon DOS of atom j respectively. Although the calculated DOS from the force constants of 32-atom (or 64-atom, not shown here) SQScell resembles maximally with the experimental DOS, the ν [ζ00] [ζζ0] [ζζζ] Γ X L Γ L T L LT T T ( T H z ) DOS (Arbit. scale) fcc Ni Pt
50 50
FIG. 4. (Color online) Same as Fig. 3, but for fcc Ni Pt al-loy. Blue square symbols in both left and right panels indicatethe experimental data. calculated band edge is still less than the measured one.This is an inherent problem of LDA-based calculations,which usually underestimates the band edge of the calcu-lated DOS and are also reflected via the bulk modulus.In summary, we propose a systematic first principlescalculation of the interatomic force constants for disor-dered alloys. SQS structures of different cell size are usedto capture the effects of random environment at differ-ent length scales. Two alloy systems with very differentintrinsic properties (i.e. lattice type, masses, force con-stants etc.) are investigated. In bcc TaW alloy, Ta-Taforce constants are predicted to be stiffer compared tothose in pure Ta, however W-W force constants behaveoppositely. In fcc NiPt alloy, Ni-Ni (Pt-Pt) force con-stants within the disordered environment behave softer(stiffer) than those in pure Ni (Pt). Calculated averagebond lengths between each pair of atoms are found toclosely dictate the nature of force constants. For both thealloys, the prediction of bond length distribution and theforce constants are found to improve with increasing sizeof the SQS-cell; in particular, the 32-atom and 64-atomSQS cell for bcc Ta W yields force constants whichagree fairly well with experiment. Ab-initio predictionof force constants in random alloys is one of the key thrustof the present work. Calculated phonon dispersion andthe DOS also compare reasonably well with available ex-perimental data. We propose the present method as apotential solution to the microscopic understanding offorce constants in disordered alloys. Future studies onmore complex alloys, addressing the effects of magneticorder, lattice mismatch etc. on force constants, are on-going which will provide an even more stringent test ofthe present methodology.We thank D. Alf`e from University college London(UK)for helpful discussion on the symmetry of dynamical ma-trices. ∗ †
Corresponding author: [email protected]; [email protected] D. W. Taylor, Phys. Rev. , 1017 (1967). G. Renaud et al. , Phys. Rev.
B 38 , 5944 (1988). D. A. Rowlands et al. , Phys. Rev.
B 72 , 045101 (2005);D. A. Biava et al. , Phys. Rev.
B 72 , 113105 (2005). S. Ghosh et al. , Phys. Rev.
B 66 , 214206 (2002). A. Zunger et al. , Phys. Rev. Lett. , 353 (1990). A. Alam and A. Mookerjee,
Phys. Rev.
B 69 , 024205(2004). Aftab Alam et al. , Phys. Rev.
B 75 , 134202 (2007). A. van de Walle et al. , Rev. Mod. Phys. , 11 (2002). C. K. Gan et al. , Phys. Rev.
B 73 , 235214 (2006); A. vande Walle and G. Ceder, Phys. Rev.
B 61 , 5972 (2000). D. Alf`e, Comp. Phys. Commun. , 2622 (1999). Chao Jiang et al. , Phys. Rev.
B 69 , 214202 (2004); A.van de Walle and G. Ceder, Phys. Rev.
B 61 , 5972(2000); J. von Pezold et al. , Phys. Rev.
B 81 , 094203(2010); Yao-Ping Xie et al. , Comput. Mater. Sci. ,2586 (2011); For other fcc SQS structure, See URL : See Supplemental Material at [URL] for details of the oc-cupation number, transformation matrices and SQS struc-tures for fcc and bcc lattices. G. Kresse and J. Furthm¨uller, Phys. Rev.
B 54 , 11169(1996); Comput. Mat. Sci.
B 6 , 15 (1996). G. Kresse and D. Joubert, Phys. Rev.
B 59 , 1758 (1999). C. Leroux,M.C. Cadeville, V. Pierron-Bohnes, G. Inden,F. Hinz, Journal of Physics F , 2033-2051 (1988). Z.I. Tolmacheva, V.I. Kornilova, Izvestiya Akademii NaukSSSR, Metally B. J. Higuera e t al., Phys. Rev. B 31 , 730 (1985). Oscar Granas e t al., J. Phys.: Condens. Matter , 015402(2012). A. D. B. Woods, Phys. Rev.