Linear Scaling ‘Order-N’ Methods in Electronic Structure Theory

1. Linear Scaling ‘Order-N’ Methodsin Electronic Structure Theory Richard M. Martin University of Illinois Acknowledgements:Pablo Ordejon David Drabold Matthew Grumbach Uwe StephanDaniel Sanchez-Portal Satoshi Itoh Thanks to: Jose Soler, Emilio Artacho, Giulia Galli, ... Comp. Mat. Science School 2001

2. Linear Scaling ‘Order-N’ Methodsand Car-Parrinello Simulations • Fundamental Issues of locality in quantum mechanics • Paradigm for view of electronic properties • Practical Algorithms • Results Comp. Mat. Science School 2001

3. Locality in Quantum Mechanics • V. Heine (Sol. St. Phys. Vol. 35, 1980) “Throwing out k-space” Based on ideas of Friedel (1954) , . . . • Many properties of electrons in any region are independent of distant regions • Walter Kohn “Nearsightness” Comp. Mat. Science School 2001

4. Locality in Quantum Mechanics • Which properties of electrons are independent of distant regions? • Total integrated quantitiesDensity, Forceson atoms, . . . • Coulomb Forces are long range but they can be handled in O(N) fashion just as in classical systems Comp. Mat. Science School 2001

5. Non-Locality in Quantum Mechanics • Which properties of electrons arenon-local? • Individual Eigenstates in crystals • Sharp features of the Fermi surface at low T • Electrical Conductivity at T=0Metals vs insulators: distinguished by delocalization of eigenstates at the Fermi energy (metals) vs localization of the entire many-electron system (insulators) • Approach in the Order-N methods:Identify localized and delocalized aspects Comp. Mat. Science School 2001

6. Density Matrix I • Key property that describes the range of the non-locality is the density matrix r(r,r’) • In an insulator r(r,r’) is exponentially localized • In a metalr(r,r’)decays as a power law at T = 0, exponentially for T > 0.(Goedecker, Ismail-Beigi) • For non-interacting Bosons or Fermions, Landau and Lifshitz show that the correlation functiong(r,r’)is uniquely related to thesquare ofr(r,r’) • Thus correlation lengths and the density matrix generally become shorter range at high T Comp. Mat. Science School 2001

7. Density Matrix II • Key property that describes the range of the non-locality is the density matrix r(r,r’) • Definition:r(r,r’) = Sifi *(r)fi (r’) • Can be localized even if each fi *(r)is not! Atom positions r’ r fixed at r =0 Comp. Mat. Science School 2001

8. Toward Working Algorithms I (My own personal view) Heine and Haydock laid the groundwork - but it was applied only to limited Hamiltonians, …. 1985 - Car-Parrinello Methods changed the picture Key quantity isthe total energy E[{fi}] which does not require eigenstates - only traces over the occupied states - the {fi} can be linear combinations of eigenstates Comp. Mat. Science School 2001

9. Toward Working Algorithms II How can we use the advantages of the Car-Parrinello and the local approaches? 1992 - Galli and Parrinello pointed out the key idea - to make a Car-Parrinello algorithm that takes advantage of the locality Require that the states in localized. Note this does not require a localized basis - it may be very convenient, but a localized basis is not essential to construct localized states (example: sum of plane waves can be localized) Comp. Mat. Science School 2001

10. Toward Working Algorithms III What are localized combinations of the eigenfunctions? Wannier Functions (generalized)! Wannier Functions span the same space as the eigenstates - all traces are the same Extended Bloch Eigenfunctions One localized Wannier Ftn centered on each site Wannier Functions Comp. Mat. Science School 2001

11. Toward Working Algorithms IV Can work with either localized Wannier functions wi (r) or localized density matrix r(r,r’) = Sifi *(r)fi (r’) = Si wi *(r)wi (r’) Functions of one variable But not unique Functions of two variables - more complex But unique Comp. Mat. Science School 2001

12. Linear Scaling ‘Order-N’ Methods • Computational complexity ~ N= number of atoms (Current methods scale as N2 or N3) • Intrinsically Parallel • “Divide and Conquer” • Green’s Functions • Fermi Operator Expansion • Density matrix “purification” • Generalized Wannier Functions • Spectral “Telescoping”(Review by S. Goedecker in Rev Mod Phys) Comp. Mat. Science School 2001

13. Divide and Conquer (Yang, 1991) • Divide System into (Overlapping) Spatial Regions.Solve each region in terms only of its neighbors.(Terminate regions suitably) • Use standard methods for each region • Sum charge densities to get total density, Coulomb terms Comp. Mat. Science School 2001

14. Expansion of the Fermi function • Sankey, et al (1994); Goedecker, Colombo (1994); Wang et al (1995) • Explicit T nonzero • Projection into the occupied Subspace • Multiply trial functionf by “Fermi operator”:Ff = [(H - EF)/KBT +1]-1f • Localized f leads to localized projection since the Fermi operator (density matrix) is localized • Accomplish by expandingF in power series in H operator - Comp. Mat. Science School 2001

15. Density Matrix “Purification” • Li, Nunes, Vanderbilt (1993); Daw (1993)Hernandez, Gillan (1995) • Idea: A density matrix rat T=0 has eigenvalues = occupation = 0 or 1 • Suppose we have an approximater that does not have this property • The relationrn+1 = 3 (rn )2 - 2 (rn )3 always produces a new matrix with eigenvalues closer to 0 or 1. Comp. Mat. Science School 2001

16. Density Matrix “Purification” • The relationrn+1 = 3 (rn )2 - 2 (rn )3 always produces a new matrix with eigenvalues closer to 0 or 1. 3 x2 - 2 x3 1 instability x 1 Comp. Mat. Science School 2001

17. Generalized Wannier Functions • Divide System into (Overlapping) Spatial Regions. • Require each Wannier function to be non-zero only in a given region • Solve for the functions in each region requiring each to be orthogonal to the neighboring functions • New functional invented to allow direct minimization without explicitly requiring orthogonalization • Mauri, et al.; Ordejon, et al; 1993; Stechel, et al 1994; Kim et al 1995 Comp. Mat. Science School 2001

18. Generalized Wannier Functions • Factorization of the density matrixr(r,r’) = Siwi* (r ) wi(r’) • Can chose localized Wannier functions (really linear combinations of Wannier functions) • Minimize functional: E = Tr [ (2 - S) H] • Since this is a variational functional, the Car-Parrinello method can be used to use one calculation as the input to the next • Mauri, et al.; Ordejon, et al; 1993; Stechel, et al 1994; Kim et al 1995 Overlap matrix Comp. Mat. Science School 2001

19. Functional (2-S)(H - EF) • Minimization leads to orthonormal filled orbitals focres empty orbitals to have zero amplitude • Each matrix element (S and H) contains two factors of the wavefunction - amplitude ~ x. • For occupied states (eigenvalues below EF) - ( 2 x2 - x4 ) x 1 Minimum at zero for empty states above EF Minimum for normalized wavefunction (x = 1) 1 Comp. Mat. Science School 2001

20. Example of Our workPrediction of Shapes of Giant Fullerenes S. Itoh, P. Ordejon, D. A. Drabold and R. M. Martin, Phys Rev B 53, 2132 (1996).See also C. Xu and G. Scuceria, Chem. Phys. Lett. 262, 219 (1996). Comp. Mat. Science School 2001

21. U. Stephan Wannier Function in a-Si Comp. Mat. Science School 2001

22. Combination of O(N) Methods Comp. Mat. Science School 2001

23. Collision of C60 Buckyballs on Diamond Galli and Mauri, PRL 73, 3471 (1994) Comp. Mat. Science School 2001

24. Deposition of C28 Buckyballs on Diamond • Simulations with ~ 5000 atoms, TB Hamiltonian from Xu, et al. ( A. Canning, G.~Galli and J .Kim, Phys.Rev.Lett. 78, 4442 (1997). Comp. Mat. Science School 2001

25. Daniel Sanchez-Portal(Phys. Rev. Lett. 1999) Simulation of a gold nanowire pulled between two gold tips Full DFT simulation Explanation for very puzzling experiment! Thermal motion of the atoms makes some appear sharp, others weak in electron microscope Example of DFT Simulation (not order N) Comp. Mat. Science School 2001

26. Simulations of DNA with the SIESTA code • Machado, Ordejon, Artacho, Sanchez-Portal, Soler (preprint) • Self-Consistent Local Orbital O(N) Code • Relaxation - ~15-60 min/step (~ 1 day with diagonalization) Iso-density surfaces Comp. Mat. Science School 2001

27. HOMO and LUMO in DNA (SIESTA code) • Eigenstates found by N3 method after relaxation • Could be O(N) for each state Comp. Mat. Science School 2001

28. O(N) Simulation of Magnets at T > 0 • Collaboration at ORNL, Ames, Brookhaven • Snapshot of magnetic order in a finite temperature simulation of paramagnetic Fe. These calculations represent significant progress towards the goal of full implementation of a first principles theory of the finite temperature and non-equilibrium properties of magnetic materials. • Record setting performance for large unit cell models (up to 1024-atoms) led to the award of the 1998 Gordon Bell prize. • The calculations that were the basis for the award were performed using the locally self-consistent multiple scattering method, which is an O(N)Density Functional method • Web Site:http://oldpc.ms.ornl.gov/~gms/MShome.html Comp. Mat. Science School 2001

29. FUTURE! ---- Biological Systems • Examples of oriented pigment molecules that today are being simulated by empirical potentials Comp. Mat. Science School 2001

30. Conclusions • It is possible to treat many thousands of atoms in a full simulation - on a workstation with approximate methods - intrinsically parallel for a supercomputer • Why treat many thousands of atoms? • Large scale structures in materials - defects, boundaries, …. • Biological molecules • The ideas are also relevant to understanding even small systems Comp. Mat. Science School 2001