Title: Quantum-ESPRESSO: The SCF Loop and Some Relevant Input Parameters
1Quantum-ESPRESSOThe SCF Loop and Some Relevant
Input Parameters
- Sandro Scandolo
- ICTP
- (most slides courtesy of Shobhana Narasimhan)
2www.quantum-espresso.org
3Quantum-ESPRESSO the project
- Quantum ESPRESSO stands for opEn Source Package
for Research in Electronic Structure, Simulation,
and Optimization. It is freely available to
researchers around the world under the terms of
the GNU General Public License. - Quantum ESPRESSO is an initiative of the
DEMOCRITOS National Simulation Center (Trieste)
and SISSA (Trieste), in collaboration with the
CINECA National Supercomputing Center in Bologna,
the Ecole Polytechnique Fédérale de Lausanne,
Princeton University, the Massachusetts Institute
of Technology, and Oxford University. Courses on
modern electronic-structure theory with hands-on
tutorials on the Quantum ESPRESSO codes are
offered on a regular basis in developed as well
as developing countries, in collaboration with
the Abdus Salam International Centre for
Theoretical Physics in Trieste.
4What can QE do (I)
- Ground-state calculations
- Self-consistent total energies, forces,
stresses - Electronic minimization with iterative
diagonalization techniques, damped-dynamics,
conjugate-gradients - Kohn-Sham orbitals
- Gamma-point and k-point sampling, and a
variety of broadening schemes (Fermi-Dirac,
Gaussian, Methfessel-Paxton, and
Marzari-Vanderbilt) - Separable norm-conserving and ultrasoft
(Vanderbilt) pseudo-potentials, PAW (Projector
Augmented Waves) - Several exchange-correlation functionals
from LDA to generalized-gradient corrections
(PW91, PBE, B88-P86, BLYP) to meta-GGA, exact
exchange (HF) and hybrid functionals (PBE0,
B3LYP, HSE) - Hubbard U (LDAU)
- Berry's phase polarization
- Spin-orbit coupling and noncollinear
magnetism - Maximally-localized Wannier functions using
WANNIER90 code. - Structural Optimization
- GDIIS with quasi-Newton BFGS
preconditioning - Damped dynamics
- Ionic conjugate-gradients minimization
- Projected velocity Verlet
5What can QE do (II)
- Born-Oppenheimer Molecular Dynamics
- o Microcanonical (Verlet) dynamics
- o Isothermal (canonical) dynamics -
Anderson, Berendsen thermostats - o Isoenthalpic, variable cell dynamics
(Parrinello-Rahman) - o Constrained dynamics
- o Ensemble-DFT dynamics (for
metals/fractional occupations) - Response properties (density-functional
perturbation theory) - Phonon frequencies and eigenvectors at any
wavevector - Full phonon dispersions inter-atomic force
constants in real space - Translational and rotational acoustic sum
rules - Effective charges and dielectric tensors
- Electron-phonon interactions
- Third-order anharmonic phonon lifetimes
- Infrared and (non-resonant) Raman
cross-sections - EPR and NMR chemical shifts using the GIPAW
method
6QE platforms
- Platforms
-
- Runs on almost every conceivable current
architecture from large parallel machines (IBM SP
and BlueGene, Cray XT, Altix, Nec SX) to
workstations (HP, IBM, SUN, Intel, AMD) and
single PCs running Linux, Windows, Mac OS-X,
including clusters of 32-bit or 64-bit Intel or
AMD processors with various connectivity (gigabit
ethernet, myrinet, infiniband...). Fully exploits
math libraries such as MKL for Intel CPUs, ACML
for AMD CPUs, ESSL for IBM machines.
7How to thank the authors
8Useful information about input variables
9INPUT_PW.html
10 The Kohn-Sham problem
- Want to solve the Kohn-Sham equations
- Note that self-consistent solution necessary, as
H depends on solution - Convention
H
11Self-consistent Iterative Solution
How to solve the Kohn-Sham eqns. for a set of
fixed nuclear (ionic) positions.
Yes
12Plane Waves Periodic Systems
- For a periodic system
- The plane waves that appear in this expansion can
be represented as a grid in k-space
where G reciprocal lattice vector
ky
- Only true for periodic systems that grid is
discrete. - In principle, still need infinite number of plane
waves.
kx
13Truncating the Plane Wave Expansion
- In practice, the contribution from higher Fourier
components (large kG) is small. - So truncate the expansion at some value of kG.
- Traditional to express this cut-off in energy
units
ky
Input parameter ecutwfc
14Step 0 Defining the (periodic) system
Namelist SYSTEM
15How to Specify the System
- All periodic systems can be specified by a
Bravais Lattice and an atomic basis.
16How to Specify the Bravais Lattice / Unit Cell
Input parameter ibrav
- - Gives the type of Bravais lattice (SC, BCC,
Hex, etc.)
Input parameters celldm(i)
- Give the lengths directions, if
necessary of the lattice vectors a1, a2, a3
a2
a1
- Note that one can choose a non-primitive unit
cell - (e.g., 4 atom SC cell for FCC structure).
17Atoms Within Unit Cell How many, where?
Input parameter nat
- Number of atoms in the unit cell
Input parameter ntyp
- Number of types of atoms
FIELD ATOMIC_POSITIONS
- - Initial positions of atoms (may vary when
relax done). - Can choose to give in units of lattice vectors
(crystal) - or in Cartesian units (alat or bohr or
angstrom)
18Step 1 Obtaining Vnuc
Yes
19Nuclear Potential
- Electrons experience a Coulomb potential due to
the nuclei. - This has a known and simple form
- But this leads to computational problems!
20Problem for Plane-Wave Basis
- Core wavefunctions sharply peaked near
nucleus.
Valence wavefunctions lots of wiggles near
nucleus.
High Fourier components present i.e., need large
Ecut ?
21Solutions for Plane-Wave Basis
- Core wavefunctions sharply peaked near
nucleus.
Valence wavefunctions lots of wiggles near
nucleus.
High Fourier components present i.e., need large
Ecut ?
Dont solve for the core electrons!
Remove wiggles from valence electrons.
22Pseudopotentials for Quantum Espresso - 1
- Go to http//www.quantum-espresso.org Click on
PSEUDO
23Pseudopotentials for Quantum Espresso - 2
- Click on element for which pseudopotential wanted.
24Pseudopotentials for Quantum Espresso - 3
- Pseudopotentials name gives information about
- type of exchange-correlation functional
- type of pseudopotential
- e.g.
25Element Vion info for Quantum Espresso
ATOMIC_SPECIES Ba 137.327 Ba.pbe-nsp-van.UPF Ti
47.867 Ti.pbe-sp-van_ak.UPF O 15.999
O.pbe-van_ak.UPF
- NOTE
- Should have same exchange-correlation functional
for all pseudopentials. - ecutwfc, ecutrho depend on type of
pseudopotentials used (should test for system
property of interest).
26Step 2 Initial Guess for n(r)
Yes
27Initial Choice of n(r)
- Various possible choices, e.g.,
- Superpositions of atomic densities.
- Converged n(r) from a closely related calculation
- (e.g., one where ionic positions slightly
different). - Approximate n(r) , e.g., from solving problem in
a smaller/different basis. - Random numbers.
Input parameter startingwfc
28Initial Choice of n(r)
29Step 3 VH VXC
Yes
30Exchange-Correlation Potential
- VXC ? dEXC/dn contains all the many-body
information. - Known numerically, from Quantum Monte Carlo
various analytical approximations for
homogeneous electron gas. - Local Density Approximation
- Excn ? n(r) VxcHOMn(r) dr
- -surprisingly successful!
- Generalized Gradient Approximation(s) Include
terms involving gradients of n(r)
Replace
pz
(in name of pseudopotential)
pw91, pbe
(in name of pseudopotential)
31Step 4 Potential Hamiltonian
Yes
32Kohn-Sham equations in plane wave basis
- Eigenvalue equation is now
- Matrix elements are
- Ionic potential given by
33Step 5 Diagonalization
Expensive! ?
Yes
34Exact Diagonalization is Expensive!
- Have to diagonalize (find eigenvalues
eigenfunctions of) H kG,kG - Typically, NPW gt 100 x number of atoms in unit
cell. - Expensive to store H matrix NPW2 elements to be
stored. - Expensive (CPU time) to diagonalize matrix
exactly, - NPW3 operations required.
- Note, NPW gtgt Nb number of bands required
Ne/2 or a little more (for metals). - So ok to determine just lowest few eigenvalues.
35Iterative Diagonalization
- Can recast diagonalization as a minimization
problem. - Then use well-established techniques for
iterative minimization by searching in the space
of solutions, - e.g., Conjugate Gradient.
- Another popular iterative diagonalization
technique is the Davidson algorithm.
Input parameter diagonalization
-which algorithm used for iterative
diagonalization
Input parameter nbnd
-how many eigenvalues computed
36Step 6 New Charge Density
Yes
37Brillouin Zone Sums
0
- Many quantities (e.g., n, Etot) involve sums over
k. - In principle, need infinite number of ks.
- In practice, sum over a finite number BZ
Sampling. - Number needed depends on band structure.
- Typically need more ks for metals.
- Need to test convergence wrt k-point sampling.
eF
e
k
38Types of k-point meshes
0
- Special Points Chadi Cohen
- Points designed to give quick convergence
for particular crystal structures. - Monkhorst-Pack
- Equally spaced mesh in reciprocal space.
- May be centred on origin non-shifted or
not shifted
K_POINTS tpiba automatic crystal gamma
1st BZ
If automatic, use M-P mesh
nk1, nk2, nk3, k1, k2, k3
b1
shift
nk1nk24
39Step 7 Test for Convergence
Yes
40How To Decide If Converged?
- Check for self-consistency. Could compare
- New and old wavefunctions / charge densities.
- New and old total energies.
- Compare with energy estimated using
Harris-Foulkes functional.
Input parameter conv_thr
-Typically OK to use 1.e-08
Input parameter electron_maxstep
-Maximum number of scf steps performed
41Step 8 Mixing
Can take a long time to reach self-consistency!
?
Yes
42Mixing - 1
- Iterations n of self-consistent cycle
- Successive approximations to density
- nin(n) ? nout(n) ? nin(n1).
- nout(n) fed directly as nin(n1) ?? No, usually
doesnt converge. - Need to mix, take some combination of input and
output densities (may include information from
several previous iterations). - Goal is to achieve self consistency (nout nin )
in as few iterations as possible.
43Mixing - 2
- Simplest prescription linear mixing
- nin(n1) b nout(n) (1-b) nin(n).
- There exist more sophisticated prescriptions
(Broyden mixing, modified Broyden mixing of
various kinds)
Input parameter mixing_mode
Input parameter mixing_beta
- Typically use value between 0.1 0.7
- (depends on type of system)
44Mixing - 3
45Output Quantities
Yes
46Output Quantities
- (Converged) Diagonalization ? ei, yi
- Strictly speaking, only n(r) Etot are
correct. - yy ? n
- To get Etot Sum over eigenvalues, correct for
double-counting of Hartree XC terms, add
ion-ion interactions. - Very useful quantity!
- Can use to get structures, heats of formation,
adsorption energies, diffusion barriers,
activation energies, elastic moduli, vibrational
frequencies,
47Geometry Optimization Using Etot
- Simplest case only have to vary one degree of
freedom - - e.g., structure of diatomic molecule
- - e.g., lattice constant of a cubic (SC, BCC,
FCC) crystal - Can just look for minimum in binding curve
48 Forces
- Need for geometry optimization and molecular
dynamics. - Could get as finite differences of total energy -
too expensive! - Use force (Hellmann-Feynman) theorem
- - Want to calculate the force on ion I
- - Get three terms
- When is an eigenstate,
- -Substitute this...
-
49 Forces (contd.)
0
- The force is now given by
- Note that we can now calculate the force from a
calculation at ONE configuration alone huge
savings in time. - If the basis depends upon ionic positions (not
true for plane waves), would have extra terms
Pulay forces. - should be exact eigenstate, i.e., scf
well-converged!
0
Input parameter tprnfor
50An Outer Loop Ionic Relaxation
Inner SCF loop for electronic iterations
Move ions
Outer loop for ionic iterations
Forces 0?
Structure Optimized!
51Geometry Optimization With Forces
0
- Especially useful for optimizing internal degrees
of freedom, surface relaxation, etc. - Choice of algorithms for ionic relaxation, e.g.,
steepest descent, BFGS.
calculation relax
NAMELIST IONS
Input parameter ion_dynamics
52Structure of the input file
53Input file namelists
54Input file namelists
55Input file input_cards
56Input file input_cards
57Input file a simple example