Using multibody energy expansions from abinitiocalculations for computation

Using multi-body energy expansions from ab-initiocalculations for computation of alloy phase structures V. Sundararaghavan and Prof. Nicholas Zabaras Materials Process Design and Control Laboratory Sibley School of Mechanical and Aerospace Engineering 188 Frank H. T. Rhodes Hall Cornell University Ithaca, NY 14853 -3801 Email: zabaras@cornell. edu URL: http: //mpdc. mae. cornell. edu Materials Process Design and Control Laboratory

PREDICTION OF STABLE STRUCTURES Ca Cu o. P 12 Stable Pt clusters (Doye and Wales, New J. Chem. , 1998) h. P 6 Stable configurations of adsorbed species Computational techniques -Exhaustive or heuristic search aided by DFT calculations -Cluster expansion Materials Process Design and Control Laboratory

Comparison with CE Cluster expansion • • • Only configurational degrees of freedom Relaxed calculation required but only a few calculations required Periodic lattices, Explores superstructures of parent lattice Multi-body expansion • • • Configurational and positional degrees of freedom Relaxed DFT calculations are not required Periodicity is not required Requires a large number of cluster energy evaluations Convergence issues Materials Process Design and Control Laboratory

Hybrid cluster expansions • Allow positional degrees of freedom in cluster expansions • For periodic lattices Cluster expansion for the fixed lattice Pair potentials for local relaxations Geng, Sluiter et al, Phys Rev B 2006 Materials Process Design and Control Laboratory

Multi-body expansion = ∑ +∑ +… Position and species Total energy Symmetric function JW Martin - Journal of Physics C, 1975, Empirical potentials (3 body): Murrell-Mottram (Mol. Phys 1990) Materials Process Design and Control Laboratory

Multi-body expansion Example of calculation of multibody potentials E 1(X 1) = V (1)(X 1) E 1(X 2) = V (1)(X 2) E 2(X 1, X 2) = V (2)(X 1, X 2) + V (1)(X 1) + V (1)(X 2) Evaluate (ab-initio) energy of several two atom structures to arrive at a functional form of E 2(X 1, X 2) Inversion of potentials V (2)(X 1, X 2) = E 2(X 1, X 2) - (E 1(X 1) + E 1(X 2) ) = Increment in energy due to pair interactions Drautz, Fahnle, Sanchez, J Phys: Condensed matter, 2004 Materials Process Design and Control Laboratory

Multi-body expansion Inversion of potentials Calculation of energies EL is found from ab-initio energy database, L << M Drautz, Fahnle, Sanchez, J Phys: Condensed matter, 2004 Materials Process Design and Control Laboratory

Fitting energy surfaces To calculate the energy of a 3 body structure (E 3), we need to identify E 2 and E 1, values. • Two body energy E 2(X 1, X 2) is the energy of an isolated cluster of 2 atoms at positions X 1 and X 2. • The database may not contain this energy since the energy values have only been obtained for atoms at locations (xi, yi) that are different from (X 1, X 2) • We use interpolation methods for retrieving energy at (X 1, X 2) from the database of energies at (xi, yi). For example, we can use a polynomial interpolation of the form: Interpolation allows us to compute a large number of energies from a well-sampled database Materials Process Design and Control Laboratory

Smolyak algorithm Extensively used in statistical mechanics Provides a way to construct interpolation functions based on minimal number of points Uni-variate interpolation Multi-variate interpolation Smolyak interpolation Accuracy the same as tensor product Within logarithmic constant Increasing the order of interpolation increases the number of points sampled Materials Process Design and Control Laboratory

Smolyak algorithm: reduction in points For 2 D interpolation using Chebyshev nodes Left: Full tensor product interpolation uses 256 points Right: Sparse grid collocation used 45 points to generate interpolant with comparable accuracy Results in multiple orders of magnitude reduction in the number of points to sample For multi-atom systems, sample all combinations of atoms (eg. E(A-A-A), E(A-A-B), E(A-B-B), E(B-B-B) and construct interpolants. Materials Process Design and Control Laboratory

CLUSTER REPRESENTATION Specification of clusters of various order by position variables 4 4 5 a 1 2 a 3 b 1 2 A point in 6 dimensional space b 3 5 • Convex hull technique to represent all atoms in the positive z-direction • Use independent coordinates to represent the cluster geometry Materials Process Design and Control Laboratory

CLUSTER ENERGY COMPUTATIONS • Executables – Cluster coordinates – Energy interpolation – Batch input for PWSCF – Read energies from PWSCF – Energy calculation Computations were performed in parallel on a 64 node quad-processor LINUX cluster • Plane-wave electronic density functional program ‘quantum espresso’ (http: //www. pwscf. org) calculations are used to compute energies given the atomic coordinates and lattice parameters. • These calculations employ LDA and use ultra-soft pseudopotentials. • Single k-point calculations were used for isolated clusters, the cell size was selected so that the effect of periodic neighbors are negligible. • For multi-component systems, a constant energy cutoff equal to cutoff for the hardest atomic potential (e. g. B in B-Fe-Y-Zr) is used. MP smearing (ismear=1, sigma=0. 2) is used for the metallic systems. Materials Process Design and Control Laboratory

LINKING THE MULTIBODY EXPANSION TO OTHER SOFTWARE Multi Body Expansion (MBE) The multibody expansion software written in C++ Two parts: potential generation & energy computation Energy computation part is the Hamiltonian Molecular dynamics- LAMMPS Large-scale Atomic/Molecular Massively Parallel Simulator (LAMMPS) is a classical molecular dynamics (MD) code developed by S. Plimpton et. al (Sandia national lab) Directly linked energy computation part in LAMMPS with MBE Useful for molecular dynamics and energy minimization Monte Carlo- MCCCS Towhee Monte Carlo for Complex Chemical Systems (MCCCS) developed by M. G. Martin, J. I. Siepmann et. al. Available at http: //towhee. sourceforge. net/ Fortran based code. Linked Towhee and MBE using a library Performs a variety of calculations in all ensembles Materials Process Design and Control Laboratory

ENERGY SURFACES FOR ISOLATED CLUSTERS Y X Platinum isolated cluster energies computed using multi-body potentials a = (1. 5*X+0. 5)*7. 5 Bohr b = (1. 5*Y+0. 5)*7. 5 Bohr b a Materials Process Design and Control Laboratory

COMPUTATION OF CLUSTER ENERGIES Distance between atoms 1 -3 (Bohr) The complete potential surface for a 3 Pt cluster. Figure (a) shows computed Platinum threeatom cluster energies, while (b) shows extension of energies using pair potential terms beyond the cutoff. Distance between atoms 1 -2 (Bohr) (a) Distance between atoms 1 -2 (Bohr) (b) Materials Process Design and Control Laboratory

Convergence results for different energy functions Energy of the system scales as n 2 where n is the number of atoms Energy of the system scales as n 1/2 where n is the number of atoms Order of interactions necessary for full convergence: 2 Order of interactions necessary for full convergence: 4 5 body term contribution =0 3 body term contribution =0 Materials Process Design and Control Laboratory

Convergence results for different energy functions Using pair potentials: Lennard Jones for Helium atoms Order of interactions necessary for full convergence is 2 as expected (since it is a pair potential) 3 body term contribution =0 Materials Process Design and Control Laboratory

Oscillations in MB energy for complex energy functionals EAM potentials: Platinum system Energies (En) calculated from an n-body expansion -Energies oscillate around the true energy -Approach: Low pass filtering (convolution operation) that cuts off high frequency oscillations. correct energy -Compute the energy at the minima using self consistent field calculation Materials Process Design and Control Laboratory

Computation of MBE energy filters + + +. . Weighted MBE Is the total energy correlated with structural energies of clusters ? Materials Process Design and Control Laboratory

Weighted MB energy a 1 a 2 True energy a 3 Materials Process Design and Control Laboratory

Extrapolatory tests on weighted MBE 16 atom Au-Cu FCC cluster Au. Cu 3 4 unit cell, 4 at/cell True energy Weighted MBE energies once built for a small set of configurations provide accurate energy fit for various different inter-atomic distances within that configuration. MBE 4 th order Materials Process Design and Control Laboratory

Selection of order of expansion Test various MBE orders in extrapolatory modes Weighted 4 th order MBE Cohesive energy (Ryd) Weighted order MBE True energies Weighted 3 rd order MBE Weighted MBE expansion coefficients are fitted using 12 atom cluster energies and the results are presented for a 16 atom cluster. Cohesive energy (Ryd) 2 nd True energies Energies may differ but the weighted MBE captures the energy minima within 4 th order expansion. Materials Process Design and Control Laboratory

Platinum clusters 16 atom FCC cluster 0 Weighted MBE 4 th order Cohesive energy (Ryd) -1 -2 4 Number of isolated cluster calculations 120 4 560 4 1820 Depth of interpolation Actual energy + -3 -4 + -5 -6 6 6. 5 7 7. 5 8 Lattice parameter (Bohr) 8. 5 9 Energy minima • Coefficients obtained using an 12 atom cluster energies at different lattice parameters Materials Process Design and Control Laboratory

Platinum clusters 24 atom FCC cluster 1 0 Cohesive energy (Ryd) -1 Weighted MBE 4 th order -2 -3 Actual energy 4 Number of isolated cluster calculations 276 4 2024 4 10626 Depth of interpolation -4 + -5 -6 -7 -8 -9 + 6 6. 5 7 7. 5 8 Lattice parameter (Bohr) 8. 5 9 Energy minima Materials Process Design and Control Laboratory

Platinum clusters A random 24 atom configuration -1 Weighted MBE 4 th order Cohesive energy (Ryd) -2 -3 Actual energy 4 Number of isolated cluster calculations 276 4 2024 4 10626 Depth of interpolation -4 + -5 -6 -7 + -8 -9 6 6. 5 7 7. 5 8 Lattice parameter (bohr) 8. 5 9 Materials Process Design and Control Laboratory

Stable phase structures of Au-Cu alloy Super-cell approach For computing stable structures of periodic lattices, a 4 x 4 x 4 supercell (216 atoms) is used as an approximation. Weighted MBE is several orders of magnitude faster than a relaxed DFT calculation. Useful for amorphous structures Small cluster calculations are used to compute the weights in the weighted MBE expansion FCC structures are considered here for Au-Cu. Materials Process Design and Control Laboratory

Stable phase structures of Au-Cu alloy Au. Cu 3 cell relaxation a = 6. 62 bohr = 3. 50 A Au 3 Cu cell relaxation a = 7. 3 bohr = 3. 86 A 3 x 3 x 3 supercell Au. Cu 3 lattice parameter: 3. 76 A Au 3 Cu lattice parameter: 4. 04 A Materials Process Design and Control Laboratory

Stable phase structures of Au-Cu alloy Au. Cu 3 cell relaxation a = 6. 71 bohr = 3. 55 A Au 3 Cu cell relaxation a = 7. 4 bohr = 3. 92 A 4 x 4 x 4 supercell Au. Cu 3 lattice parameter: 3. 76 A Au 3 Cu lattice parameter: 4. 04 A Materials Process Design and Control Laboratory

APPLICATION TO SURFACE PHENOMENA Minimum energy surface of h on Pt(111) Plot of minimum energy in z direction for the primitive cell Highly anharmonic potential energy surface FCC->HCP (55 mev), FCC->TOP (160 mev) H confined to FCC-HCP-FCC valleys (Baskar and Zabaras, 2007) FCC site G. Kallen, G. Wahnstrom, Quantum treatment of H on a Pt(111) surface, Phys Rev B, 65 (2001) Materials Process Design and Control Laboratory

Conclusions • MB expansion provides atom position dependent potentials that are used to identify stable structures. • Ab-initio database of cluster energies are created and interpolation for various cluster positions are generated using efficient sparse grid interpolation algorithms. • Weighted MBE is fast and captures the energy minima within a small order of expansion. • Technique is applicable to study stability of amorphous systems, molecules and clusters. Materials Process Design and Control Laboratory