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)
Comments 0
Log in to post a comment