Fast hybrid density-functional computations using plane-wave basis - PowerPoint PPT Presentation

Fast hybrid density-functional computations using plane-wave basis sets Paolo Giannozzi, Universit` a di Udine, Italy Workshop A fresco of contemporary condensed-matter physics in honour of Pino Grosso – Pisa, 27-28 September 2018 Work mostly done by Ivan Carnimeo, in collaboration with Stefano Baroni (SISSA), in the framework of MaX - Materials at the eXascale Centre of Excellence – Typeset by Foil T EX –

Background and motivations • In many fields of science, atomistic simulations have become very successful ... • ... in particular, first-principle simulations based on density-functional theory (DFT), pseudopotentials (PP) and a plane-wave (PW) basis sets • Hybrid functionals, containing some admixture of exact (Hartee-Fock) exchange, currently offer the best performances in terms of accuracy of results ... • ... but they are so slooow with PWs! The goal of this work: to produce a fast, reliable and robust method to compute the exact-exchange potential in a PW-PP framework

Density-Functional Theory: Kohn-Sham approach Let us introduce the orbitals ψ i for an auxiliary set of non-interacting electrons whose charge density is the same as that of the true system: � | ψ i ( r ) | 2 , � ψ i | ψ j � = δ ij n ( r ) = i Let us rewrite the energy functional in a more manageable way as: � E = T s [ n ( r )] + E H [ n ( r )] + E xc [ n ( r )] + n ( r ) V ( r ) d r where T s [ n ( r )] is the kinetic energy of the non-interacting electrons: � h 2 T s [ n ( r )] = − ¯ � ψ ∗ i ( r ) ∇ 2 ψ i ( r ) d r , 2 m i E H [ n ( r )] is the Hartree energy , due to electrostatic interactions: � n ( r ) n ( r ′ ) E H [ n ( r )] = e 2 | r − r ′ | d r d r ′ , 2

E xc [ n ( r )] is called exchange-correlation energy (a reminiscence from the Hartree-Fock theory) and includes all the remaining (unknown!) energy terms. Minimization of the energy with respect to ψ i yields the Kohn-Sham (KS) equations: � � h 2 − ¯ 2 m ∇ 2 + V ( r ) + V H ( r ) + V xc ( r ) ψ i ( r ) = ǫ i ψ i ( r ) , � �� � H KS where the Hartree and exchange-correlation potentials: � n ( r ′ ) V H ( r ) = δE H [ n ( r )] V xc ( r ) = δE xc [ n ( r )] = e 2 | r − r ′ | d r ′ , δn ( r ) δn ( r ) depend self-consistently upon the ψ i via the charge density. Self-consistency is achieved by looping over the charge density, using some suitable mixing algorithm , such as e.g. simple mixing. At iteration l + 1 : n ( l +1) = βn ( l ) out + (1 − β ) n ( l ) in , 0 < β < 1 in

Hybrid functionals The “standard” functionals, LDA (local density approximation) or GGA (generalized gradient approximation), are functions of the local density n ( r ) and for GGA of gradient |∇ n ( r ) | : � n ( r ) ǫ GGA ( n ( r ) , |∇ n ( r ) | ) d r E xc = Better-performing hybrid functionals , such as B3LYP or PBE0, contain some amount of exact exchange, as in Hartree-Fock theory (spin-restricted): � ψ ∗ � i ( r ) ψ j ( r ) ψ ∗ j ( r ′ ) ψ i ( r ′ ) n ( r ) V ( r ) d r + E H − e 2 � d r d r ′ E = T s + | r − r ′ | 2 i,j � �� � E HF x In hybrid DFT: E hyb xc � � �� � α x E HF + (1 − α x ) E DF T + E DF T E = T s + n ( r ) V ( r ) d r + E H + , x x c with α x = 0 . 2 ÷ 0 . 3 .

Exchange potential The exchange potential V x is, unlike V xc in GGA, nonlocal and depends upon the � density matrix , γ ( r , r ′ ) = ψ j ( r ) ψ ∗ j ( r ′ ) : j � ψ j ( r ) ψ ∗ j ( r ′ ) ψ i ( r ′ ) δE x i ( r ) = − e 2 � d r ′ ( V x ψ i )( r ) ≡ δψ ∗ | r − r ′ | j More generally, it may be defined it via its action on a generic orbital: � ψ j ( r ) ψ ∗ j ( r ′ ) φ i ( r ′ ) ( V x φ i )( r ) ≡ − e 2 � d r ′ | r − r ′ | j Let us assume that our orbitals are expanded into PWs: h 2 1 ¯ � 2 m | k + G | 2 ≤ E cut e i ( k + G ) · r , √ ψ j ( r ) = c j, k + G N Ω G k = Bloch vector, G = reciprocal lattice vectors, E cut = kinetic energy cutoff, Ω = unit cell volume, N Ω = crystal volume.

Exchange potential and plane waves With PW, orbitals are represented in Fourier space as { c j, k + G } , j = 1 , ..., M , where M is the number of orbitals. V x can be easily computed using fast Fourier transforms (FFT) (in the following, simplified case, k = 0 , orbitals are real): F F T • φ i ( G ) − → φ i ( r ) (FFT brings from a discrete grid of G into a discrete grid of r ) • ρ ij ( r ) = ψ j ( r ) · φ i ( r ) F F T • ρ ij ( r ) − → ρ ij ( G ) • v ij ( G ) = ρ ij ( G ) /G 2 F F T • v ij ( G ) − → v ij ( r ) • ( V X φ i )( r ) = � j ψ j ( r ) · v ij ( r ) F F T • ( V X φ i )( r ) − → ( V X φ i )( G ) Each operation is very quick ... but O ( M 2 ) such operations are needed!

Hybrid functionals for plane waves Double-loop algorithm Pure DFT Outer loop: Do n = 1, ... The exchange r ′ ϕ ( n ) r ′ ) ϕ ( n ) N r ′ ) ( � k ( � � integrals have to ( V X ϕ ( n ) ϕ ( n ) � i k )( � r ) = − ( � r ) d � SLOW!! i | � r − � r ′ | be evaluated at i { ϕ ( n ) , ψ ( n , 0) } every iteration, i i even though the Inner loop: Do m = 1, ... { ϕ ( n ) } functions i are not varying r ′ ϕ ( n ) r ′ ) ψ ( n , m ) N r ′ ) ( � ( � � ( V X ψ ( n , m ) ϕ ( n ) � i k in the inner )( � r ) = − ( � r ) d � SLOW!! k i | � r − � r ′ | i ones. Iterative diagonalization Check Convergence on { ψ ( n , m ) } at fixed { ϕ ( n ) } i i End Inner loop Check Convergence on { ϕ ( n ) } i End Outer loop Giannozzi et al., JPCM , 21, 395502 (2009) Ivan Carnimeo (UniUd/SISSA) Exact exchange with plane waves 2017-10-02 4 / 18

Exchange with localized orbitals: inner projection We project the exchange operator over the set of occupied orbitals − 1 � ϕ j | ˆ � � ˆ V X | ϕ i � � ϕ i | ˆ ˆ W X = V X | ϕ j � V X = | ξ i � � ξ i | ij i � ˆ V X | ϕ i � L − T | ξ k � = ik i Adaptively Compressed Exchange (ACE) Inner projection methods. The application of ˆ W X during the SCF is equivalent to ˆ V X within the subspace of the projection � ˆ W X | ψ k � = | ξ i � � ξ i | ψ k � i Lin, JCTC , 12, 2242 (2016); Damle, Lin, Ying, JCTC , 11, 1463 (2015); L owdin, IJQC , 4, 231 (1971). Ivan Carnimeo (UniUd/SISSA) Exact exchange with plane waves 2017-10-02 6 / 18

Hybrid functionals for plane waves Double-loop algorithm Pure DFT Outer loop: Do n = 1, ... With the localization step r ′ ϕ ( n ) r ′ ) ϕ ( n ) N r ′ ) ( � k ( � � the exchange ( V X ϕ ( n ) ϕ ( n ) � i k )( � r ) = − ( � r ) d � SLOW!! i | � r − � r ′ | integrals in the i − 1 � ϕ ( n ) | V X | ξ ( n ) � � ξ ( n ) | = V X | ϕ ( n ) � � ϕ ( n ) | V X | ϕ ( n ) � outer loop now � W X = i involve two W X | ϕ ( n ) | ξ ( n ) � � ξ ( n ) | ϕ ( n ) � k � = − k � localized i i i functions. { ξ ( n ) , ψ ( n , 0) } i i Lin, JCTC , 12, Inner loop: Do m = 1, ... 2242 (2016); W X | ψ ( n , m ) | ξ ( n ) � � ξ ( n ) | ψ ( n , m ) � � = − � Damle, Lin, Ying, k i i k i JCTC , 11, 1463 (2015). Iterative diagonalization Check Convergence on { ψ ( n , m ) } at fixed { ξ ( n ) } i i End Inner loop Check Convergence on { ϕ ( n ) } i End Outer loop Ivan Carnimeo (UniUd/SISSA) Exact exchange with plane waves 2017-10-02 7 / 18

Hybrid functionals for plane waves With the ACE method Wall time along SCF iterations the cost of the inner 1000 iterations is analogous ACE to the cost of a pure Full 800 Exchange step DFT method. Wall time (sec) 600 Example: Ethylene, 1 CPU, 100 Ry, converged to 10 − 7 Ry 400 Pure DFT: 6 200 Outer loop: 5 Inner loop: 2, 3, 3, 2 0 0 2 4 6 8 10 12 14 16 Tot. calls: 32/5 #Iteration Ivan Carnimeo (UniUd/SISSA) Exact exchange with plane waves 2017-10-02 8 / 18

Localizing orbitals The key to achieve better performances is orbital localization . There are various ways to localize orbitals in real space, the most famous being the transformation to Wannier orbitals . The density matrix is naturally localized: for insulators, γ ( r , r ′ ) ∼ exp( − α | r − r ′ | ) asymptotically. This means that w ( r ) = � j c j ψ j ( r ) , with c j = ψ ∗ j ( r ′ ) are localized! The set of w ’s above is redundant: we need a way to select the best manifold among the many possibles, without storing the density matrix. This can be achieved using the selected columns of the density matrix (SCDM) algorithm. This is an algebraic procedure that requires a standard QR decomposition (actually, only the pivoting part): Lin, JCTC , 12, 2242 (2016); Damle, Lin, Ying, JCTC , 11, 1463 (2015). Once orbitals are localized, exchange integrals between non-overlapping orbitals can be discarded.

Hybrid functionals for plane waves N � r ′ ) w k ( � r ′ ) r ′ w i ( � � ( ˆ V X w k )( � r ) = − w i ( � r ) d � | � r − � r ′ | i The products between two canonical orbitals are much more delocalized than the product of two localized orbitals and this can be exploited in order to reduce the cost of the exchange integrals evaluations. Ivan Carnimeo (UniUd/SISSA) Exact exchange with plane waves 2017-10-02 10 / 18