Integrated Computational Materials Engineering: Calculation of Equation of State Using Density Functional Theory

1. Integrated Computational Materials Engineering EducationCalculation of Equation of State UsingDensity Functional Theory Mark Asta1, Katsuyo Thornton2, Raúl Enrique2, and Larry Aagesen3 1Department of Materials Science and Engineering, University of California, Berkeley 2Department of Materials Science and Engineering, University of Michigan, Ann Arbor 3Idaho National Laboratory

2. Purposes of Density Functional Theory Module • Understand fundamentals of Density Functional Theory (DFT) • Apply DFT to calculate: • Equilibrium lattice constant • Bulk Modulus • Components of elastic constant tensor • Understand how to check for convergence of results

3. Example: Equation of StateA Probe of Interatomic Interactions Schematic Energy vs. Volume Relation Diamond Cubic Structure of Si Energy per atom a a a http://www.e6cvd.com/cvd/page.jsp?pageid=361 Volume per atom (=a3/8 for Si)

4. Equation of StateWhat Properties Can we Learn from It? Pressure versus Volume Relation Equilibrium Volume (or Lattice Constant) Bulk Modulus Recall 1st Law of Thermo: dE = T dS - P dV and consider T = 0 K B related to curvature of E(V) Function

5. Generalize to Non-Hydrostatic DeformationExample of Uniaxial Deformation Lz Lz(1+e) Ly Ly Lx Lx Definition of Deformation In Terms of Strain: (All other strains are zero)

6. Linear-Elasticity for Single Crystals Linear Elasticity: Stress Tensor Elastic Constant Tensor Strain Tensor Voigt Notation: 11 → 1, 22 → 2, 33 → 3, 23 → 4, 13 → 5, 12 → 6 Elastic Energy: In example from previous slide: (All other strains are zero) Note: for cubic crystal C11=C22=C33, C12=C13=C23

7. DFT as a tool • DFT can evaluate energy from atomic configuration • ”Real” materials: accurate physics and chemistry • Use this information to evaluate: • Mechanical properties • Phonon structure • free energy • You can almost use it as a black box BUT…

8. Total Energy in Density Functional Theory Electron Density Electron Wavefunctions Potential Electrons Feel from Nuclei Exchange-Correlation Energy Form depends on whether you use Local Density Approximation (LDA) or Generalized Gradient Approximation (GGA)

9. Kohn-Sham EquationsSchrödinger Equation for Electron Wavefunctions Exchange-Correlation Potential Electron Density Note:fi depends on n(r) which depends on fi  Solution of Kohn-Sham equations must be done iteratively Convergence criterium (e.g. 10-3 eV/atom)

10. Electronic BandstructureExample for Si Bandstructure Brillouin Zone http://en.wikipedia.org/wiki/Brillouin_zone http://de.wikipedia.org/wiki/Datei:Band_structure_Si_schematic.svg Electronic wavefunctions in a crystal can be indexed by point in reciprocal space (k) and a band index (b)

11. Why?Wavefunctions in a Crystal Obey Bloch’s Theorem For a given band b Where is periodic in real space: Translation Vectors: The envelope function represents delocalized distribution of electron density

12. Why?Wavefunctions in a Crystal Obey Bloch’s Theorem For a given band b Each crystal cell interacts with its periodic counterparts: This information must enter the calculation Where is periodic in real space: Translation Vectors: The envelope function represents delocalized distribution of electron density

13. Representation of Electron Density Integral over k-points in first Brillouin zone In practice the integral over the Brillouin zone is replaced with a sum over a finite number of k-points (Nkpt) Band occupation (e.g., the Fermi function) One parameter that needs to be checked for numerical convergence is number of k-points

14. Representation of WavefunctionsPlane-Wave Basis Set For a given band Use Fourier Expansion In practice the Fourier series is truncated to include all G for which: Another parameter that needs to be checked for convergence is the “plane-wave cutoff energy” Ecut

15. Examples of Convergence Checks Effect of Ecut Effect of Number of k Points Note: the different values of kTel corresponds to different choices for occupation function (wj in slide 14) http://www.fhi-berlin.mpg.de/th/Meetings/FHImd2001/pehlke1.pdf

16. DFT Module • Problem 1: Calculate equilibrium volume and bulk modulus of diamond cubic Si using Quantum Espresso on Nanohub (http://nanohub.org/) • Outcome 1: Understand effect of numerical parameters on calculated results by testing convergence with respect to number of k-points and plane-wave cutoff • Outcome 2: Understand the effect of the theoretical model for exchange-correlation potential on the accuracy of the calculations by comparing results from Local Density Approximation (LDA) and Generalized Gradient Approximation (GGA) with experimental measurements

17. DFT Module • Problem 2: Calculate the single-crystal elastic constants C11 and C12 • Outcome 1: Understand how to impose homogeneous elastic deformations in a DFT calculation • Outcome 2: Understand the effect of the theoretical model for exchange-correlation potential on the accuracy of the calculations by comparing results from Local Density Approximation (LDA) and Generalized Gradient Approximation (GGA) with experimental measurements

18. DFT Module • For problem 1 you will make use of the unit cell for diamond-cubic Si shown below. You will vary only the lattice constant a. a Unit Cell Vectors a1 = a (-1/2, 1/2 , 0) a2 = a (-1/2, 0, 1/2) a3 = a (0, 1/2, 1/2) a a Basis Atom Positions (Fractional Coordinates) 0 0 0 ¼ ¼ ¼

19. DFT Module • For problem 2 you will impose a homogeneous tensile strain (e) along the [001] axis (see slide 4) • Such a strain results in the x3 coordinate of all atoms changing to x3*(1+e) • This homogeneous deformation can be represented by changing the unit cell vectors as follows: a a Unit Cell Vectors a1 = a (-1/2, 1/2 , 0) a2 = a (-1/2, 0, (1+e)/2) a3 = a (0, 1/2, (1+e)/2) a Basis Atom Positions (Fractional Coordinates) 0 0 0 ¼ ¼ ¼

20. Extra Slides

21. Equation of StateHow to Calculate from Density Functional Theory Formulation: for a given arrangement of nuclei defined by the lattice constant, crystal structure, and non-hydrostatic strains, compute the total energy corresponding to the optimal arrangement of the electron density Theoretical Framework: Quantum mechanical calculation of energy of electrons and nuclei interacting through Coulomb potential Practical Implementation: Density functional theory

22. Implementation of DFT for a Single Crystal Crystal Structure Defined by Unit Cell Vectors and Positions of Basis Atoms Example: Diamond Cubic Structure of Si a Unit Cell Vectors a1 = a (-1/2, 1/2 , 0) a2 = a (-1/2, 0, 1/2) a3 = a (0, 1/2, 1/2) a a Basis Atom Positions 0 0 0 ¼ ¼ ¼ All atoms in the crystal can be obtained by adding integer multiples of unit cell vectors to basis atom positions

23. Representation of WavefunctionsFourier-Expansion as Series of Plane Waves For a given band: Recall that is periodic in real space: can be written as a 3D Fourier Series: where the are primitive reciprocal lattice vectors

24. Recall Properties of Fourier Series Black line = (exact) triangular wave Colored lines = Fourier series truncated at different orders http://mathworld.wolfram.com/FourierSeriesTriangleWave.html General Form of Fourier Series: For Triangular Wave: Typically we expect the accuracy of a truncated Fourier series to improve as we increase the number of terms

25. Electron Density in Crystal Lattice Unit-Cell Vectors a1 = a (-1/2, 1/2 , 0) a2 = a (-1/2, 0, 1/2) a3 = a (0, 1/2, 1/2) a a a Electron density is periodic with periodicity given by Translation Vectors:

26. Self-Consistent Solution to DFT Equations Set up atom positions Make initial guess of “input” charge density (often overlapping atomic charge densities) Solve Kohn-Sham equations with this input charge density Compute “output” charge density from resulting wavefunctions If energy from input and output densities differ by amount greater than a chosen threshold, mix output and input density and go to step 2 Quit when energy from input and output densities agree to within prescribed tolerance (e.g., 10-5 eV) Input Positions of Atoms for a Given Unit Cell and Lattice Constant guess charge density compute effective potential compute Kohn-Sham orbitals and density different compare output and input charge densities same Energy for Given Lattice Constant Note: In this module the positions of atoms are dictated by symmetry. If this is not the case another loop must be added to minimize energy with respect to atomic positions.