Periodic Boundary Methods and Applications: Ab-initio Quantum Mechanics for Band Structures CH121b

1. Periodic Boundary Methods and Applications: Ab-initio Quantum Mechanics for Band Structures CH121b Jamil Tahir-Kheli MSC, Caltech May 4, 2011

2. Outlines • What is different about crystalline solids? • Bloch theorem • First Brillouin zone • Reciprocal space sampling • Plane wave, APW, Gaussian basis sets • SeqQuest • Crystal06

3. What is different about solids? H H H H H H H a Infinite repeating pattern of atoms with translational symmetry Even if you have 1 basis function per atom, there is still an infinite number of atoms leading to diagonalization of an infinite matrix! This implies we can never solve crystals By exploiting the translational symmetry of the crystal, we can find a way to break the problem into finite pieces that approximate the solution

4. Bloch Theorem (simplification due to translation symmetry)K-Vectors

5. Bloch Theorem (example: one dimensional hydrogen chain) H H H H H H H a

6. Bloch Theorem (example: one dimensional hydrogen chain)band structure k = 0 k = p/a

7. Density of States

8. Bloch Theorem (example: two dimensional hydrogen surface)

9. The First Brillouin zone The first Brillouin zone contains all possible interactions between two adjacent unit cells.

10. Hartree-Fock-Roothaan Equation in periodic systems Finite diagonalizations

11. We can solve for each k-point, but there are an infinite number of them By evaluating each k-point at the first Brillouin zone and summing them together, we can obtain the properties such as total energy or electron density of the system Impossible !!! In practice, the only computationally feasible approach is to approximate the full BZ integral with summation over a finite set of k-points.

12. Reciprocal Space Sampling (Monkhost-Pack grids)

13. Differences between Molecularand Periodic Codes There is an infinity far away from the molecule where the density decays to zero as an exponential. The exponent is the ionization potential (up to a factor) and can be shown to equal the HOMO eigenvalue. DFT obtains exact density and thus IP.

14. There is no vacuum away from infinite crystal where we can define the zero of the electrostatic potential. No physical significance can be attached to the Kohn-Sham eigenvalues for solid calculations. Empirically, we do it anyway.

15. Ionization potential Fermi level Orbital energies are arbitrary up to a constant. To obtain the work functions, you need to know the surface charge distribution of a finite sample.

16. Reference Density (SeqQuest) Gaussian Orbitals Ewald (CRYSTAL) Ab-Initio Methods Augmented Plane Waves Plane Wave FLAPW, Wien2k VASP “Exact”  GW

17. SIESTA Numerical Basis Sets DMOL3 Green’s Function (GW)

18. Plane Waves Basis functions for each k in Brillouin Zone, where G is a reciprocal lattice vector. Solve for wavefunctions and energies,

19. Practically, to obtain a finite set of states, the basis functions are cutoff, The cutoff is quoted as an energy, or as a cutoff wavelength,

20. Assembling the Fock matrix to diagonalize is easy with Plane waves. Kinetic Nuclear Coulomb + Exchange

21. Problem: cutoff G must be chosen extremely large to capture variation of wavefunction near nuclei. Fock matrix to diagonalize cheap to assemble, but large. Diagonalization becomes time consuming. CASSTEP is a pure plane wave code.

22. Augmented Plane Wave codes try to reduce the number of basis functions of pure plane wave by using atomic orbitals in the vicinity of nuclei that are smoothly joined to plane waves in the interstitial region. Constant potential in interstitial regions Self-Consistent spherical potential inside spheres Wavefunctions in two regions are smoothly joined

23. APW works well for computing band structures, but has three drawbacks: 1.) There are no standard basis functions. This makes it difficult to visualize the wavefunction in terms of atomic orbitals. Mulliken populations are hard to quantify. 2.) Exact exchange is hard to compute. Thus, modern hybrid functionals that include Hartree-Fock exchange are not presently available with this approach. 3.) There is a certain arbitrariness to the choice of sphere radii. Wien2k and FLAPW

24. = + + + = + ….. GW Method Feynman diagram method Poles of propagator are physical excitation energies. Gives good bandgaps and excitations, but computationally very very expensive. Not competitive with DFT.

25. Gaussian Orbitals Trial wavefunctions for crystal momentum k are built up from linear combinations of localized atomic Gaussian orbitals. Atomic Gaussian localized at R

26. Advantages: 1.) Fewer basis functions needed to solve problem. 2.) Intuitive wavefunctions that are easily visulalized. 3.) Mulliken populations 4.) Can do surface problems Disadvantage: 1.) Much harder to calculate elements in Fock matrix.

27. SeqQuest (Sequential QUantum Electronic STructure) Worked out once Varies slowly so solve in Fourier space using Poisson equation, Can obtain linear scaling!!

28. The linear scaling method does not lend itself to an easy way to compute exact Hartree-Fock exchange. HF exchange requires brute force calculation taking the scaling back to O(N^3). In fact, no one has found a fast way to compute exact exchange for periodic systems. If you can, PUBLISH!

29. GaN Quest Input Deck Bohr

31. Online manual for Quest http://www.cs.sandia.gov/~paschul/Quest/

32. CRYSTAL: A Gaussian CodeInput Structure of CRYSTAL Structure Basis set (atomic orbital) Method (HF or DFT) SCF control

33. Input Structure of CRYSTAL (example) Your personal note about this calculation “crystal” “slab” “polymer” “molecule” Space group sequence number Cell parameters Number of non-equivent atoms Atomic coordiantes Basis set

34. Input Structure of CRYSTAL (Basis set) atomic number For example: C: 6 O: 8 Ni: 28 Ni: 228 all electron basis set effective core potential number of shells Si (1s22s22p63s23p2) 14 electrons Si ((function)3s23p2) 4 electrons shell (orbital) type 0: s orbital 1: s+p orbital 2: p orbital 3: d orbital 4: f orbital 1st shell basis set type 0: input by hand 1: STO-nG 2: 3(or 6)-21G scale factor 2nd shell 3th shell 4th shell End of basis set section number of Gaussian functions number of electrons at this shell

35. Crystal06 Input (basis set)http://www.crystal.unito.it/Basis_Sets/Ptable.html

36. Crystal06 Input (SCF control) k-point net for insulator: n n for metal n 2n maximum SCF iterations mixing control 30% P0 + 70% P1 for second step