E S DFT: GGA

baconossifiedMechanics

Oct 29, 2013 (4 years and 2 months ago)

138 views

E
LECTRONIC

S
TRUCTURE

WITH

DFT:

GGA
AND

BEYOND

Jorge
Kohanoff

Queen’s University Belfast

United Kingdom


j.kohanoff@qub.ac.uk



Kohn
-
Sham equations:













This

PDE

must

be

solved

self
-
consistently
,

as

the

KS

potential

depends

on

the

density,

which

is

constructed

with

the

solutions

of

the

KS

equations
.

D
ENSITY

FUNCTIONAL

THEORY

(
DFT
):

KOHN
-
SHAM

EQUATIONS



The inhomogeneous electron gas is considered as locally
homogeneous:










LDA XC hole centred at r
, interacts with the electron also at
r
.
The exact XC hole is centred at r’


This is partially compensated by multiplying the pair correlation
function with the density ratio

(r)/

(r’)









E
XCHANGE

AND

CORRELATION

IN

DFT
:

THE

LOCAL

DENSITY

APPROXIMATION

(
LDA
)

XC energy density of the HEG


Location of the XC hole
(Jones and Gunnarsson, 1982)








E
XCHANGE

AND

CORRELATION

IN

DFT
:

THE

LOCAL

DENSITY

APPROXIMATION

(
LDA
)


Favors more homogeneous electron densities


Overbinds molecules and solids (Hartree
-
Fock underbinds)


Geometries, bond lengths and angles, vibrational frequencies
reproduced within 2
-
3%


Dielectric constants overestimated by about 10%


Bond lengths too short for weakly bound systems (H
-
bonds, VDW)


Correct chemical trends, e.g. ionization energies



Atoms (core electrons) poorly described (HF is much better)


XC potential decays exponentially into vacuum regions. It should
decay as

e
2
/r.
Hence, it is poor for dissociation and ionization


Poor for metallic surfaces and physisorption


Very poor for negatively charged ions (
self
-
interaction error
)


Poor for weakly bound systems: H
-
bonds (

), VDW (
non
-
local
)


Band gap in semiconductors too small (~40%)


Poor for strong on
-
site correlations (
d

and
f

systems, oxides, UO
2
)









LDA
-
LSDA
:
TRENDS

AND

LIMITATIONS


Inhomogeneities in the density


Self
-
interaction cancellation


Non
-
locality in exchange and correlation


Strong local correlations



Gradient expansions


Weighted density approximation


Exact exchange in DFT (OEP local vs HF non
-
local)


DFT
-
HF hybrids


Self
-
interaction correction


Van der Waals

and RPA functionals


LSDA+U


Multi
-
reference Kohn
-
Sham


GW approximation (Many
-
body)









B
EYOND

THE

LDA


E
XC

expanded in gradients of the density




where


is the spin polarization

s
=|

|/2k
F


is the density gradient

And
F
XC

is the enhancement factor



First
-
order term is fine, but higher
-
order terms diverge. Only by
some re
-
summation to ∞
-
order the expansion converges.


GGA:
F
XC

is designed to fulfil a number of exactly known
properties
, e.g.
Perdew
-
Burke
-
Ernzerhof

(PBE)

1.
Exchange
: uniform scaling, LSDA limit, spin
-
scaling
relationship, LSDA linear response,
Lieb
-
Oxford bound

2.
Correlation
: second
-
order expansion, hole sum rule, vanishes
for rapidly varying densities, cancels singularity at high densities

G
RADIENT

EXPANSIONS
:

GENERALIZED

GRADIENT

APPROXIMATION

r+
d
r


⡲+
d
r)


Improves atomization and surface energies


Favors density inhomogeneities


Increases lattice parameters of metals


Favors non
-
spherical distortions


Improves bond lengths


Improves energies and geometries of H
-
bonded systems


There is error cancellation between X and C at short range



XC potential still decays exponentially into vacuum regions


Some improvement in band gaps in semiconductors


What was correct in LDA is worsened in GGA


Still incorrect dissociation limit. Fractionally charged fragments


Inter
-
configurational errors in I
P

and E
A


Error cancellation between X and C is not complete at long
-
range.
X hole is more long
-
ranged than XC hole










P
ROPERTIES

OF

THE

GGA


Combine GGA local exchange with
Hartree
-
Fock

non
-
local
exchange:




Parameter


fitted to experimental data for molecules (~0.75), or
determined from known properties.


PBE0, B3LYP, HSE06



Properties
:

1.
Q
uite accurate in many respects, e.g. energies and geometries

2.
Improve on the self
-
interaction error,
but not fully SI
-
free

3.
Improve on band gaps

4.
Improve on electron affinities

5.
Better quality than MP2

6.
Fitted hybrids unsatisfactory from the theoretical point of view


H
YBRID

FUNCTIONALS


Self
-
interaction can be removed at the level of classical
electrostatics
:








Potential is
state
-
dependent
. Hence it is not an eigenvalue
problem anymore, but a
system of coupled PDEs


Orthogonality of SIC orbitals not guaranteed
, but it can be
imposed

(Suraud)


Similar to HF, but the Slater determinant of SIC orbitals is
not invariant against orbital transformations


The result depends on the choice of orbitals (
localization
)








S
ELF
-
INTERACTION

CORRECTION

(
SIC
)

Perdew
-
Zunger 1982

Mauri, Sprik, Suraud


Van der Waals (dispersion) interactions
: are a dynamical
non
-
local correlation effect



Dipole
-
induced dipole interaction due to quantum density
fluctuations in spatially separated fragments




Functional (Dion et al 2004):





Expensive double integral


Efficient implementations (Roman
-
Perez and
Soler

2009)


Good approximations based on dynamical response theory


Beyond VDW: Random Phase Approximation (
Furche
)















D
YNAMICAL

CORRELATION
:
VDW


1
(
r
,t
)


2
(
r
,t
)


1
(
t
)

E
1
(
t).

2
(
t
)



=
VDW kernel fully non
-
local.
Depends on

(
r
)
and

(
r’
)


Empirical approaches:





With
f
a function that removes smoothly the singularity at
R=0, and interferes very little with GGA (local) correlation.



Grimme

(2006):
C
6

parameters from atomic calculations.
Extensive parameterization: DFT
-
D.


Tkatchenko

and
Scheffler

(2009):
C
6

parameters
dependent
on the density.



Random Phase Approximation (RPA):
captures VDW and
beyond. Can be safely combined wit exact exchange (SIC).
Infinite order perturbation (like Coupled clusters in QC).

Furche

(2008);
Paier

et al
(2010);
Hesselmann

and
Görling

(2011)















D
YNAMICAL

CORRELATION
:
VDW

AND

BEYOND


Strong onsite Coulomb correlations are ot captured by LDA/GGA


These are important for localized (
d

and
f
) electronic bands, where
many electrons share the same spatial region:
self
-
interaction
problem



Semi
-
empirical solution: separate occupied and empty state by an
additional energy
U

as in Hubbard’s model:





This induces a splitting in the KS eigenvalues:








S
TRONG

STATIC

CORRELATION
:
LSDA
+
U

f
i
=orbital occupations









S
UMMARY

OF

DFT

APPROXIMATIONS

ELECTRONIC

STRUCTURE

OF

UO
2


Using the quantum
-
espresso package
(
http://www.quantum
-
espresso.org/
)



Pseudopotentials


Plane wave basis set


P
ROPERTIES


fluorite structure


fcc
, 3 atoms un unit cell


Lattice constant = 10.26 Bohr


Electronic insulator.
E
g
=2.1
eV


Electronic configuration of U:
[
Rn
]7s
2
6d
1
5f
3


U
4+
:
f
2


5f
-
band partially occupied
(2/7)


UO
2
:
splitted

by crystal field:

t
1u
(3)+t
2u
(3)+
a
g
(1)


Still partially occupied

(2/3)


Jahn
-
Teller distortion opens gap.

P
SEUDOPOTENTIAL

C
ONVERGENCE

WITH

ENERGY

CUTOFF

E
NERGY
-
VOLUME

CURVE

GGA(PBE)
DENSITY

OF

STATES

GGA+U
DENSITY

OF

STATES

GGA+U
DENSITY

OF

STATES
:
DISTORTED


An

important

family

of
Room

Temperature

Ionic

Liquids

(
Green
solvents
)


Competing

electrostatic

vs
dispersion

interactions


Large

systems

studied

with

DFT,
within

LDA
or

GGA

[Del
Popolo
,
Lynden
-
Bell and
Kohanoff
, JPCB
109
, 5895 (2005)]


Force

fields

fitted

to

DFT
-
GGA
calculations

[
Youngs
, Del
Popolo

and
Kohanoff
, JPCB
110
, 5697 (2006)]


Electrostatics

well

described

in DFT (LDA
or

GGA)


Dispersion

(van
der

Waals
)
interactions

are
absent

in
both
, LDA and GGA



V
AN

DER

WAALS

FOR

IMIDAZOLIUM

SALTS

DFT
IMIDAZOLIUM

SALTS

R
ESULTS

FOR

SIMPLE

DIMERS

Ar
2

and Kr
2

(C
6
H
6
)
2

Bond lengths:
5
-
10%

too long

Binding energies:
50
-
100%

too large

M. Dion et al, PRL
92
, 246401 (2004)

R
ESULTS

FOR

SOLIDS

Polyethylene

Silicon



Reasonable results for molecular systems



Keeps GGA accuracy for covalent systems




General purpose functional

S
OLID

FCC

A
RGON

(E.
A
RTACHO
)

Some overbinding, and lattice constant still 5% too large


… but much better than PBE

(massive underbinding and lattice constant 14% too large)

T
HE

DOUBLE

INTEGRAL

PROBLEM



(
q
1
,
q
2
,
r
12
) decays as
r
12
-
6



E
c
nl
= (1/2)




d
3
r
1

d
3
r
2


(
r
1
)


(
r
2
)


(
q
1
,q
2
,r
12
) can be truncated for
r
12

>
r
c

~ 15Å



In principle O(
N
) calculation for systems larger than 2
r
c
~ 30Å



But... with

x

~ 0.15Å (
E
c
=120Ry) there are ~(2

10
6
)
2

= 4

10
12

integration points



Consequently, direct evaluation of vdW functional is much more
expensive than LDA/GGA

F
ACTORING


(
Q
1
,
Q
2
,
R
12
)

G. Román
-
P
ė
rez

and J. M. Soler, Phys. Rev. Lett. 2009

Expand


in a basis set of functions
p

(
q
)

FT

I
NTERPOLATION

AS

AN

EXPANSION

20 grid points are sufficient

Smoothening of


required at small q

Basis functions
p

(
q
)

interpolate between grid points




Lagrange polynomials: grid given by zeros of orthogonal polynomials



Cubic splines
: grid points defined on a
logarithmic mesh


O(
N
LOG
N
)
ALGORITHM

do, for each grid point
i


find

i

and


i



find
q
i
=q(

i
,


i
)


find


i
=

i

p

(q
i

)



end do

Fourier
-
transform



i





k



do, for each reciprocal vector
k


find
u

k

=





(k)


k




end do

Inverse
-
Fourier
-
transform
u

k



u

i



do, for each grid point
i


find

i

,


i

, and
q
i



find


i
,



i
/

i

, and



i
/




i





find

v
i

end do

Implemented into
SIESTA
, but
not SIESTA
-
specific:

Input:


i

on a regular grid

Output:
E
xc

, v
i
xc

on the grid

No need for supercells in solids

No cutoff radius of interaction

A
LGORITHM

EFFICIENCY

Message

If you can simulate a system with LDA/GGA,
you can also simulate it with vdW
-
DFT


System


Atoms

CPU time in
GGA
-
XC

CPU time in
vdW
-
XC

vdW/GGA
overhead

Ar
2

2

0.75 s (44%)

7.5 s (89%)

400%

MMX polymer

124

1.9 s (2%)

10.6 s (16%)

17%

DWCN

168

11.9 s (0.6%)

109 s (5.2%)

4%

I
MIDAZOLIUM

CRYSTALS
:
V
OLUMES

LDA

PBE

VDW

EXP

[
emim
][PF6]

877.4

1088.3

1059.1

1023.9

[
bmim
][PF6]

513.4

636.6

620.0

605.0

[
ddmim
][PF6]

1710.9

2258.0

2095.5

2000.6

[
mmim
][
Cl
]

607.9

750.5

728.6

687.6

[
bmim
][
Cl
]
-
o

843.8

1039.9

1019.8

961.1

[
bmim
][
Cl
]
-
m

836.3

1051.9

1024.6

966.7

Error

-
13.5%

+8.5%

+5%

I
MIDAZOLIUM

CRYSTALS
:
H
EXAFLUOROPHOSPATES

bmimpf6

emimpf6

ddmimpf6

[bmim][PF6]

[ddmim][PF6]

C
4
-
C
4

C
1

-
C
1


C
12
-
C
12


C
4

-
C
4


C
2
-
C
4


C
4

-
C
1


P
-
P

C
2
-
P

C
1

-
P

C
4

-
P

C
12

-
P

Error

LDA

5.18

4.28

4.83

4.52

4.23

6.00

6.66

3.55

3.95

4.27

4.35

6.0%

PBE

5.62

4.24

5.44

5.12

5.10

6.56

6.89

4.00

4.75

4.96

4.78

7.0%

VDW

5.47

4.29

5.23

4.90

4.79

6.44

6.80

3.82

4.22

4.68

4.65

3.3%

EXP

5.47

4.36

5.01

4.85

4.61

6.25

6.45

4.23

4.01

4.62

4.63

C5
-
C5

C1”
-
C1”

&?¶
-
C4’

C2
-
C4’

C4’
-
C1”

P
-
P

C2
-
P

C1”
-
P

C4’
-
P

Error

LDA

4.32

4.24


3.73

4.61

3.64

5.05

3.70

3.88

4.75

7.0%

PBE

4.64

4.85

3.99

5.07

3.97

5.69

4.16

4.21

4.95

2.6%

VDWE

4.60

4.55

4.05

4.80

4.07

5.64

4.21

4.20

4.91

2.6%

VDW

4.54

4.60

4.04

4.80

3.95

5.62

3.90

4.11

4.87

2.0%

EXP

4.56

4.52

4.05

5.08

3.86

5.54

3.95

4.29

4.94

I
MIDAZOLIUM

CRYSTALS
:
H
EXAFLUOROPHOSPATES

I
MIDAZOLIUM

CRYSTALS
:
C
HLORIDES

[bmim][Cl]
-
m

C5
-
C5

C1’
-
C1’


-
C1’

C2
-
Cl

C1’
-
Cl

LDA

3.07

3.54

3.37

3.34

3.81

PBE

3.94

3.86

3.63

3.48

3.74

VDW

3.34

3.82

3.56

3.56

3.77

C5
-
C5

C1”
-
C1”

C4’
-
C4’

C2
-
C4’

C4’
-
C1”

C2
-
Cl

C1”
-


C4’
-


LDA

4.87

4.50

5.17

3.25

3.83

3.28

3.51

3.70

PBE

5.26

5.14

5.72

3.63

3.81

3.38

3.86

3.89

VDW

5.21

4.97

5.42

3.53

3.84

3.46

3.84

3.94

EXP

5.19

4.96

5.33

3.48

3.62

3.39

3.82

3.80

[mmim][Cl]

[
BMIM
][
CL
]:
POLYMORPHISM

[bmim][Cl]
-
Monoclinic

[bmim][Cl]
-
Orthorhombic

Energy difference per neutral ion pair

LDA

PBE

VDW

FF(
CLaP
)


䔨
-


[
eV
]

-
0.05

-
0.06

-
0.01

0.08

I
MIDAZOLIUM

CLUSTERS
:
T
RIFLATES

[bmim][Tf]

[bmim][Tf]
2

LDA

PBE

VDW

FF(
CLaP
)




[
eV
]

-
1.411

-
1.053

-
1.453

-
1.457

Dimer association energy

LDA

PBE

VDW

FF(
CLaP
)

Error


7.2%

2.3%

1%

--

Geometry

C
ONCLUSIONS

1.
The

non
-
empirical

van

der

Waals

functional

of

Dion

et

al
.

(DRSLL)

improves

significantly

the

description

of

the

geometry

of

imidazolium

salts
.


2.
Volumes

are

improved

respect

to

PBE,

but

still

overestimated

by

5
%
.


3.
Energetics

is

also

improved
.

It

is

similar

to

that

of

empirical

force

fields

such

as

CLaP
.


4.
The

cost

of

calculating

the

van

der

Waals

correlation

correction

is

10

times

that

of

PBE
.

However,

in

a

self
-
consistent

calculation

for

100

atoms

the

overhead

is

only

20
%
.

T
HE

THEORETICAL

LANDSCAPE


Size (number of atoms)

1000

10000

1000000

100000

100

High

Low


Treat relevant part of the system quantum
-
mechanically, and the
rest classically.



The problem is how to match the



two regions. Easy for non
-
bonded



interactions, more difficult for



chemical bonds



One can also treat part of the



system as a polarizable



continuum, or reaction field (RF)







QM/MM

Q
UANTUM

M
ECHANICS

IN

A

LOCAL

BASIS

T
IGHT

BINDING

T
IGHT

BINDING

T
IGHT

BINDING

MODELS

FOR

WATER

Ground
-
up philosophy


Water molecule


1.
Minimal basis. On
-
site energies to reproduce
band structure

a)

1s orbital for H:

Hs

b)

2s & 2p
orbitals

for O:

Os
,

Op


2.
O
-
H hopping integrals:

a)

Values at equilibrium length to reproduce
HOMO
-
LUMO
gap:
t
ss

, t
sp


b)

GSP functional form
. Cut
-
off between first and second neighbours


3.
Charge transfer: Hubbard terms fitted to reproduce
dipole moment
:
U
O
, U
H


4.
O
-
H pair potential:
GSP

form, fitted to reproduce
bond length
and
symmetric stretching

force constant


5.
Crystal field parameter

spp

selected to reproduce
polarizability

A. T. Paxton and J.
Kohanoff
, J. Chem. Phys.
134
, 044130 (2011)

T
IGHT

BINDING

MODELS

FOR

WATER

Ground
-
up philosophy


Water
dimer


6.
O
-
O hopping integrals:
t
ss

, t
sp


,
t
pp

,
t
pp




GSP

form. Cut
-
off after first neighbours


7.
O
-
O pair potential


Various forms (
GSP
, quadratic) to


reproduce
binding energy curve


8.
Fitting procedure

a)

By hand
(intuitive)

b)

Genetic

algorithm


This is the end of the fitting


All the rest are predictions

I
CE
-
XI



DFT ice sinks in water!



Polarizable model marginal



Point charge model is fine

T
IGHT
-
BINDING

LIQUID

WATER

A. T. P
AXTON

AND

J.
K
OHANOFF
, J. C
HEM
. P
HYS
.
134
, 044130 (2011)

T
IGHT

BINDING

MODEL

FOR

T
I
O
2

Band structure: A. Y.
Lozovoi
, A. T. Paxton and J.
Kohanoff

Rutile

Anatase

DFT

TB

O on
-
site energies

Os
,

Op
and O
-
O hopping integrals:
t
ss

, t
sp


, t
pp

, t
pp



W
ATER
/T
I
O
2

INTERFACES

(S
ASHA

L
OZOVOI
)

Single water adsorption: dissociation

W
ATER
/T
I
O
2

INTERFACES

(S
ASHA

L
OZOVOI
)

A water layer on TiO
2

W
ATER
/T
I
O
2

INTERFACES

(S
ASHA

L
OZOVOI
)

Bulk water on TiO
2


4 TiO
2
layers
(480 atoms)



184 water
molecules
(552 atoms)



1032 atoms



TBMD 1
ps

(one week)

W
ATER
/T
I
O
2

INTERFACES

(S
ASHA

L
OZOVOI
)

z
-
density profile of bulk water between TiO
2
surfaces

P
ROTON

TRANSFER

IN

WATER

Chemical bonds broken and formed

O
RGANIC

MOLECULES

(T
ERENCE

S
HEPPARD
)

Benzylacetone

in water

B
ENZYLACETONE

AT

W
ATER
/T
I
O
2

INTERFACE

Sasha
Lozovoi

and

Terence Sheppard

S
UMMARY


In developing a force field there are 3 main aspects to
consider. In order of relevance:


Designing the model


Choosing the properties to be reproduced


Choosing a methodology to fit those properties



Fitting procedures cannot work out their magic if the
model is not good enough



To become model
-
independent we need to introduce
electrons explicitly:
ab initio and DFT methods



These are
expensive
. Simplifications like
semi
-
empirical,

tight
-
binding and QM/MM

methods bridge the gap



Further simplifications lead to extremely efficient
bond
-
order potentials (BOP)