Transcript Modeling from the nanoscale to the macroscale

Section I: ab initio electronic structure
methods
Part II: introduction to the
computational methodology
Instructor: Marco Buongiorno Nardelli
516C Cox - tel. 513-0514 - e-mail: [email protected]
Office hours: T,Th, 2:00-3:00PM
http://ermes.physics.ncsu.edu
Suggested readings:
• R. M. Martin, Electronic Structure, Cambridge, 2004, http://www.electronicstructure.org
• W.E. Pickett, Pseudopotential methods in condensed matter applications, Comp. Phys.
Rep. 9, 115 (1989)
• J.I. Gersten, F.W. Smith, The Physics and Chemistry of Materials, Wiley, 2001
• http://www.pwscf.org
Modeling from the nanoscale to the macroscale - Fall 2004
Electronic ground state
•
Stable structure of solids are classified on the basis of their electronic ground
state, which determines the minimum energy equilibrium structure, and thus the
characteristics of the bonding between the nuclei
•
•
•
•
Closed-shell systems: rare gases and
molecular crystals. They remain atom-like
and tend to form close-packed solids
Ionic systems: compound formed by
elements of different electronegativity.
Charge transfer between the elements
thus stabilizes structures via the strong
Coulomb (electrical) interaction between
ions
Covalent bonding: involves a complete
change of the electronic states of the
atoms with pair of electrons forming
directional bonds
Metals: itinerant conduction electrons
spread among the ion cores. Electron
“gas” as electronic glue of the system
Modeling from the nanoscale to the macroscale - Fall 2004
The many-body problem
•
•
•
How do we solve for the electronic ground state? Solve a many-body problem:
the study of the effects of interaction between bodies, and the behavior of a
many-body system
The collection of nuclei and electrons in a piece of a material is a formidable
many-body problem, because of the intricate motion of the particles in the manybody system:
Electronic structure methods deal with solving this formidable problem starting
from the fundamental equation for a system of electrons ({ri}) and nuclei ({RI})
2
h
Hö  
2me
2
Z
e
1
e2
2
I
 i   | r  R |  2  | r  r |
i
i,I
i j
i
I
i
j
h2 2 1
e2

I  
2 I J | R I  R J |
I 2MI
Modeling from the nanoscale to the macroscale - Fall 2004
Electronic structure methods
•
In independent electron approximations, the electronic structure problem
involves the solution of a Schroedinger-like equation for each of the electrons in
the system
 h2 2
 s
s
ö
Heff  (r)   
 Veff (r)   i (r)   is is (r)
 2me

s
i
•
•
•
In this formalism, the ground state energy is found populating the lowest
eigenstates according to the Pauli exclusion principle
Central equation in electronic structure theory. Depending on the level of
approximation we find this equation all over:
• Semi-empirical methods (empirical pseudopotentials, tight-binding)
• Density Functional Theory
• Hartree-Fock and beyond
Mathematically speaking, we need to solve a generalized eigenvalue problem
using efficient numerical algorithms
Modeling from the nanoscale to the macroscale - Fall 2004
Towards Density Functional Theory
•
The fundamental tenet of Density Functional Theory is that the complicated
many-body electronic wavefunction  can be substituted by a much simpler
quantity, that is the electronic density
n(r) 
•
•
ö 
 n(r)
 

3
3
3
 d r2d r3 L d rN  ({r})
2
s
3
3
3
d
r
d
r
L
d
 1 2 rN ({r})
2
This means that a scalar function of position, n(r), determines all the information
in the many-body wavefunction for the ground state and in principle, for all
excited states
n(r) is a simple non-negative function subject to the particle conservation sum
rule
3
n(r)d
r N

where N is the total number of electrons in the system
Modeling from the nanoscale to the macroscale - Fall 2004
Kohn and Sham ansatz
•
•
•
•
•
H-K theory is in principle exact (there are no approximations, only two elegant
theorems) but impractical for any useful purposes
Kohn-Sham ansatz: replace a problem with another, that is the original manybody problem with an auxiliary independent-particle model
Ansatz: K-S assume that the ground state density of the original interacting
system is equal to that of some chosen non-interacting system that is exactly
soluble, with all the difficult part (exchange and correlation) included in some
approximate functional of the density.
Key steps:
• Definition of the non-interacting auxiliary system
• The auxiliary Hamiltonian contains the usual kinetic energy term and a local
effective potential acting on the electrons
Actual calculations are performed on this auxiliary Hamiltonian
HKS (r)  
1 2
 VKS (r)
2
through the solution of the corresponding Schroedinger equation for N
independent electrons
Modeling from the nanoscale to the macroscale - Fall 2004
Kohn and Sham ansatz
Non-interacting
auxiliary
particles in an
effective potential
Interacting electrons
+ real potential
•
The density of this auxiliary system is then:
n(r)    |  is (r) |2
s i 1,N
•
The kinetic energy is the one for the independent particle system:
1
1
s
2
s
Ts      i (r)   i (r)    |  is (r) |2
2 s i 1,N
2 s i 1,N
•
We define the classic electronic Coulomb energy (Hartree energy) as usual:
EHartree [n] 
1
n(r)n(r ')
3
3
d
rd
r
'
2 
|r r '|
Modeling from the nanoscale to the macroscale - Fall 2004
Kohn and Sham equations
•
Finally, we can rewrite the full H-K functional as
EKS [n]  Ts [n]   d 3rVext (r)n(r)  EHartree [n]  EII  E xc [n]
•
All many body effects of exchange and correlation are included in Exc
E xc [n]  FHK [n]  (Ts [n]  EHartree [n]) 
Tö  Ts [n]  Vöint  EHartree [n]
•
•
So far the theory is still exact, provided we can find an “exact” expression for the
exchange and correlation term
The minimization of this functional under the particle conservation constraint
leads to a set of Schroedinger-like equations
HKS is   is is
•
With an explicit effective potential
 EHartree  Exc
VKS (r)  Vext (r) 

 Vext (r) VHartree (r) Vxc (r)
 n(r)  n(r)
Modeling from the nanoscale to the macroscale - Fall 2004
Kohn and Sham equations
•
•
The great advantage of recasting the H-K functional in the K-S form is that
separating the independent particle kinetic energy and the long range Hartree
terms, the remaining exchange and correlation functional can be reasonably
approximated as a local or nearly local functionals of the electron density
Local Density Approximation (LDA): Exc[n] is a sum of contribution from each
point in space depending only upon the density at each point independent on
other points
LDA
3
E xc [n]   d rn(r) xc (n(r))
•
•
•
where  xc (n) is the exchange and correlation energy per electron.
 xc (n) is a universal functional of the density, so must be the same as for a
homogeneous electron gas of given density n
The theory of the homogeneous electron gas is well established and there are
exact expression (analytical or numerical) for both exchange and correlation
terms
0.458
Exchange as  x (n)  
where rs is defined as the average distance between
rs
4
1
electrons at a given density n :
•
rs3 
3
n
Correlation from exact Monte Carlo calculations (Ceperley, Alder, 1980)
Modeling from the nanoscale to the macroscale - Fall 2004
Kohn and Sham equations
•
•
•
•
•
Finally, the set of K-S equations with LDA for
exchange and correlation give us a
formidable theoretical tool to study ground
state properties of electronic systems
Set of self-consistent equations that have to
be solved simultaneously until convergence is
achieved
Note: K-S eigenvalues and energies are
interpreted as true electronic wavefunction
and electronic energies (electronic states in
molecules or bands in solids)
Note: K-S theory is a ground-state theory and
as such is supposed to work well for ground
state properties or small perturbations upon
them
Extremely successful in predicting materials
properties - golden standard in research and
industry
Modeling from the nanoscale to the macroscale - Fall 2004
Overview of solid-state concepts
•
•
•
•
•
In DFT calculations, one solves iteratively the Kohn-Sham equations for the
electron density in the external potential of the ions
We need to frame the K-S equations in a crystalline (or molecular or atomic)
environment and be able to formulate the problem in a computationally tractable
manner (a crystalline set-up is the most general for our purposes and for most
available computer codes)
A crystal is an ordered state of matter in which the position of the nuclei are
repeated periodically in space
Completely specified by the shape of one repeat unit (primitive unit cell) and
type and position of nuclei in that unit (basis)
The unit cell is repeated infinitely in space through a set of rules that describe
the repetition (translations) - the Bravais lattice
Crystal
Points or translation vectors
Modeling from the nanoscale to the macroscale - Fall 2004
Basis
Overview of solid-state concepts
•
Simple Bravais lattice (cubic cell): primitive cell with one atom basis:
Note: A primitive cell is defined
by 3 unit vectors directed
along the edges of the periodic
“box”
a3
a2
•
a1
The mimimal primitive cell that can generate the whole periodic lattice is called
the Wigner-Seitz cell, and depends on the particular symmetry of the crystal
(lattice + basis)
a2
a1
Wigner-Seitz Cell
•
Translations are defined as vectors T multiple of the primitive cell vectors:
T=la1+ma2+na3
Modeling from the nanoscale to the macroscale - Fall 2004
Overview of solid-state concepts
•
Any function defined in a crystal has to satisfy the periodicity of the system:
f (r  T(l,m,n))  f (r)
This applies to all the quantities that we have introduced so far: Hamiltonians,
wavefunctions, electron densities, etc.
It is well known that periodic functions can be represented by Fourier transform
in terms of Fourier components at specific wavevectors q. We define a Fourier
transform:
1
f (q) 
drf (r) exp(iq  r)

Crystal 
•
•
Crystal
•
f (q) 
For periodic functions the above expression can be written as:
1
Crystal
 
l.m.n 
Crystal
drf (r) exp(iq  (r  T(l,m,n)) 
1
Ncell
1
 exp(iq  T(l,m,n)) 
l.m.n

drf (r)exp(iq r)
cell cell
•
In a periodic crystal, Fourier components are restricted to those that are periodic
within a large volume of crystal made of Ncell:
exp(iq T(l,m,n))  1 or q T(l,m,n)  2 integer for all translations T
•
The set of Fourier components that satisfy the above condition form the
“reciprocal crystal lattice”. The reciprocal lattice is itself a Bravais lattice whose
axis bi are defined through the relation bi  a j  2 ij
Modeling from the nanoscale to the macroscale - Fall 2004
Overview of solid-state concepts
•
Since the reciprocal lattice is a Bravais lattice, it has a primitive cell. The minimal
primitive cell of the reciprocal lattice is called “Brillouin zone” and its principal
axis are given by
b1  2
a2  a3
a1 (a2  a3 )
and cyclic permutations
a2
b2
a1
b2 a 2
Wigner-Seitz Cell
•
b1
a1
b1
Brillouin Zone
The reciprocal space is comprised of a lattice of points defined by the reciprocal
lattice vectors G=lb1+mb2+nb3so that the Fourier transform of the periodic
function can be written as
1
f (G) 
drf (r) exp(iG  r)
cell 
cell
Modeling from the nanoscale to the macroscale - Fall 2004
Overview of solid-state concepts
•
•
•
•
•
Electrons in solids are subject to the external potential of the nuclei that has the
periodicity of the crystal
Independent electrons (such as the ones of DFT in the K-S approach) in a
periodic potentials have a very important property as a general consequence of
the periodicity of the potential
Bloch’s theorem: the eigenstates of the crystal Hamiltonian (the electrons’
wavefunctions) can be always chosen as the product of a plane wave times a
function with the periodicity of the Bravais lattice
 k (r)  eikruk (r) where uk (r)  uk (r  T)
with k within the primitive cell of the reciprocal lattice
Note: for the crystal we have r and T, for the reciprocal space we have k and G
All possible eigenstates of the Hamiltonian are specified by k within the primitive
cell of the reciprocal lattice. For each k there is a descrete set of eigenvalues
that form what are the energy bands of the crystal. Usually they are calculated
within the Brillouin zone and along specific directions in reciprocal space
For a complete discussion on reciprocal space vectors and periodic boundary
conditions see the Ashcroft and Mermin book or any other book on Solid-state
theory
Modeling from the nanoscale to the macroscale - Fall 2004
Example
•
Brillouin zone and band structure for a cubic fcc solid: Au
z
X
L
G
y
L
K
S
X
W
U
D
X
K
W
X
U
X
From: http://cst-www.nrl.navy.mil/ElectronicStructureDatabase/
Modeling from the nanoscale to the macroscale - Fall 2004
An electronic structure code: PWscf
•
•
PWscf is a state-of-the-art software package that we will use as an introduction
to the actual procedure of modeling the electronic structure of a solid
http://www.pwscf.org
Modeling from the nanoscale to the macroscale - Fall 2004
An electronic structure code: PWscf
Modeling from the nanoscale to the macroscale - Fall 2004
An electronic structure code: PWscf
Best first step: READ THE USER MANUAL!
Modeling from the nanoscale to the macroscale - Fall 2004
An electronic structure code: PWscf
•
Three basic preliminary steps. From the User’s manual:
•
Installation The PWscf package can be downloaded from the
http://www.pwscf.org site. Presently, only source files are provided. Some
precompiled executables (binary files) are provided only for the GUI. Providing
binaries for PWscf would require too much effort and would work only for a small
number of machines anyway. Uncompress and unpack the code in an empty
directory of your choice that will become the root directory of the distribution. On
Linux machines, you may use:
tar -xvzf pw.2.0.tgz.
On other Unix machines: gzip -dc pw.2.0.tgz | tar -xvf -
•
Automatic configuration: An experimental automatic configuration, using the
GNU "configure" utility, is available (thanks to Gerardo Ballabio, CINECA). From
the root directory, type:
./configure
•
•
The script will examine your hardware and software, generate dependencies
needed by the Makefile's, produce suitable configuration files make.sys and
make.rules. Presently it is expected to work for Linux PCs, IBM sp machines,
SGI Origin, some HP-Compaq Alpha machines. For more details, read the
INSTALL file.
Modeling from the nanoscale to the macroscale - Fall 2004
An electronic structure code: PWscf
•
•
•
•
•
•
•
Compilation There are a few adjustable parameters in
Modules/parameters.f90. The present values will work for most cases. All
other variables are dynamically allocated: you do not need to recompile your
code for a different system. You can compile the following codes:
make pw produces PW/pw.x and PW/memory.x. pw.x calculates electronic
structure, structural optimization, molecular dynamics, barriers with NEB.
memory.x is an auxiliary program that checks the input of pw.x for correctness
and yields a rough (under-)estimate of the required memory.
make ph produces PH/ph.x. ph.x calculates phonon frequencies and
displacement patterns, dielectric tensors, effective charges (uses data produced
by pw.x).
make d3 produces D3/d3.x d3.x calculates anharmonic phonon lifetimes
(third-order derivatives of the energy), using data produced by pw.x and ph.x.
make gamma produces Gamma/phcg.x. phcg.x is a version of ph.x that
calculates phonons at = 0 using conjugate-gradient minimization of the density
functional expanded to second-order. Only the ( = 0 ) point is used for Brillouin
zone integration. It is faster and takes less memory than ph.x, but does not
(yet) support Ultrasoft pseudopotentials.
make pp produces a variety of post-processing codes
make tools produces utility programs, mostly for phonon calculations
Modeling from the nanoscale to the macroscale - Fall 2004
An electronic structure code: PWscf
•
Running Pwscf for an electronic and ionic structure calculation:
•
Main information: the input file
Modeling from the nanoscale to the macroscale - Fall 2004
PWscf input file
•
•
The input data is organized as several namelists, followed by other fields
introduced by keywords. The namelists are
• &CONTROL: general variables controlling the run
• &SYSTEM: structural information on the system under investigation
• &ELECTRONS: electronic variables: self-consistency, smearing
• &IONS (optional): ionic variables: relaxation, dynamics
• &CELL (optional): variable-cell dynamics
• &PHONON (optional): information needed to produce data for phonon
calculations
Optional namelist may be omitted if the calculation to be performed does not
require them. This depends on the value of variable calculation in namelist
&CONTROL. Most variables in namelists have default values. Only the following
variables in &SYSTEM MUST be specified:
• ibrav (integer): bravais-lattice index
• celldm (real, dimension 6): crystallographic constants
• nat (integer): number of atoms in the unit cell
• ntyp (integer): number of types of atoms in the unit cell
• ecutwfc (real): kinetic energy cutoff (Ry) for wavefunctions.
Description of all the input cards can be found in the file INPUT_PW
Modeling from the nanoscale to the macroscale - Fall 2004
PWscf input file
•
Input file for an electronic structure
calculation for a bulk Si crystal
•
&control:
• scf = self-consistent solution of
the Kohn-Sham equations
(ground state energy, electron
densities) from an arbitrary initial
density (‘from_scratch’)
• prefix = prefix for file names
• tstress,tprnfor = compute
also forces and stresses in the
given geometry
• Specify working directories
(optional)
Modeling from the nanoscale to the macroscale - Fall 2004
Input parameters: &system
•
The &system namelist is the one where
we specify the geometrical parameters of
our simulation cell:
• ibrav = specifies the Bravais lattice
type
• celldm = specifies the
crystallographic constants (the
dimension of the simulation cell)
Modeling from the nanoscale to the macroscale - Fall 2004
Input parameters: &system
•
•
nat and ntyp specify the number of atoms in the simulation cell and how many
different atomic species are in the system we want to study
Together with the information contained under the keywords ATOMIC SPECIES
and ATOMIC POSITIONS they determine the basis in the primitive cell
•
ATOMIC SPECIES: list of all the atomic species in the system with information
on the atomic mass and a link to an external file that contains the information
necessary to describe that particular atom - pseudopotential file
•
ATOMIC POSITIONS: list of all the different atoms in the simulation cell with the
specification of their atomic coordinates
Modeling from the nanoscale to the macroscale - Fall 2004
Input parameters: the pseudopotential
Modeling from the nanoscale to the macroscale - Fall 2004
Numerical solution: plane waves
 h2 2

ö
Heff  i (r)   
  Veff (r)  i (r)   i  i (r)
 2me

•
•
•
•
•
•
Kohn-Sham equations are differential equations that have to be solved
numerically
To be tractable in a computer, the problem needs to be discretized via the
introduction of a suitable representation of all the quantities involved
Various discretization approeches. Most common are Plane Waves (PW) and
real space grids.
In periodic solids, plane waves of the form eiqr are most appropriate since they
reflect the periodicity of the crystal and periodic functions can be expanded in
the complete set of Fourier components through orthonormal PWs
1
 (r)   ci,q  e iqr   ci,q  q

q
q
In Fourier space, the K-S equations become

q' Höeff q ci,q   i  q' q ci,q   i ci,q'
q
q
We need to compute the matrix elements of the effective Hamiltonian between
plane waves
Modeling from the nanoscale to the macroscale - Fall 2004
Numerical solution: plane waves
•
Kinetic energy becomes simply a sum over q
1 2
1 2
 q  q  qq'
2
2
q
The effective potential is periodic and can be expressed as a sum of Fourier
components in terms of reciprocal lattice vectors
1
Veff (r)  Veff (Gm )exp(iGm r) where Veff (Gm ) 
drVeff (r)exp(iG r)


m
cell 

•
q' 
cell
•
•
Thus, the matrix elements of the potential are non-zero only if q and q’ differ by
a reciprocal lattice vector, or alternatively, q = k+Gm and q’ = k+Gm’
The Kohn-Sham equations can be then written as matrix equations
H
m,m'
(k)ci,m' (k)   i (k)ci,m (k)
m'
2
Hm,m ' (k)  k  Gm  m,m ' Veff (Gm  Gm ' )
•
where:
•
We have effectively transformed a differential problem into one that we can
solve using linear algebra algorithms!
Modeling from the nanoscale to the macroscale - Fall 2004
Input parameters: ecutwfc
•
In this representation both the potentials and the Bloch functions, solution of the
K-S problem, are expanded on a set of plane waves
 i,k (r)  eikr  ci,m (k)exp(iGm r)
•
•
•
•
•
m
In principle, the plane waves basis set is infinite, since I have an infinite number
of reciprocal lattice vectors (or, in other words, Fourier components).
In practice, we need to limit ourselves to a finite basis set for the practical
solution of the linear equations: approximation!
Remember: in Fourier space, that is in reciprocal space in a crystal, small q
components describe long-range features (wave-length), while large q
components, describe short-range features. Very sharp oscillations, for instance,
need to be described by a large number of plane waves with large G vectors
Increasing the dimension of the basis (number of plane waves, so larger
magnitudes of Gm) allows for a better description of short-range features in
either the potential or the density.
ecutwfc is the parameter that controls the number of PW’s in the basis, so it
affects directly the accuracy of the calculation: convergence parameter
|Gmax|2  ecutwfc
Modeling from the nanoscale to the macroscale - Fall 2004
Input parameters: &electrons
•
•
Kohn-Sham equations are always selfconsistent equations: the effective K-S
potential depends on the electron density that
is the solution of the K-S equations
In reciprocal space the procedure becomes:
H
m,m'
(k)ci,m' (k)   i (k)ci,m (k)
m'
Hm,m ' (k)  k  Gm
2
where
 m,m ' Veff [nk,i (Gm  Gm ' )]
and
nk,i (G)   ci (k)ci (k)
m
•
Iterative solution of self-consistent equations often is a slow process if particular tricks are
not used: mixing schemes
Modeling from the nanoscale to the macroscale - Fall 2004
Mixing schemes
•
•
•
Key problem: updating the potential and/or the
density between successive iteration steps
(loop in the solution of the K-S equations)
A direct approach, where we start from an
arbitrary density and use the solution of the KS loop as input of the next will not work instabilities in the solution of the minimization
problem - the numerical procedure does not
converges into the minimum
Simplest approach to resolve the issue is linear mixing: estimate an improved
density input nini+1 at step I+1 as the linear combinantion of input and output
densities nini and nouti at step i
i
nini1  nout
 (1  )nini
• mixing_mode = particular algorithm to mix
input and output density/potentials
• mixing_beta = numerical value of 
• conv_thr = convergence threshold for selfconsistency, estimated energy error < conv_thr
Modeling from the nanoscale to the macroscale - Fall 2004
Mixing schemes
•
•
 is the parameter that controls the rate of convergence
• Large : output density from previous iteration weighs most - fast
convergence (large density updates from one step to the other). Works well
for strongly bound, rigid systems (insulators or semiconductors, regions
around the ionic cores)
• Small : output density weighs lest - slower convergence (smaller density
updates from one step to the other) but more stable. Works well for “soft”
systems such as metals, alloys, surfaces or open systems (lots of empty
space).
Various mixing scheme beyond linear mixing have been developed and are
available in PWscf (Anderson, Broyden, etc.)
Modeling from the nanoscale to the macroscale - Fall 2004
Symmetry and special points
•
•
•
•
•
Crystals are very symmetric systems whose geometrical properties are best
described by the specification of their space group
Space group of a crystal is the group composed by the whole set of translations
and point symmetry operations (rotations, inversions, reflections and
combinations of the above) that leave the system invariant (including a particular
set of operations that combines a rotation with a non-integer translation, or glide,
of a fraction of a crystal translation vector, called non-symmorphic operations)
Any function that has the full symmetry of the crystal is invariant upon any
operation of the space group: Sng(r) = g(Snr)
In particular, since the Hamiltonian is invariant upon any symmetry operation,
any operation Sn leads to a new equation with r  Snr and k  Snk
The new solution of the transformed equation is still an eigenfunction of the
Hamiltonian with the same eigenvalue:
• A “high-symmetry” k point is defined by the identity relation: Snk  k. Helpful
in the classification of electronic states.
• One can define the Irreducible Brillouin Zone (IBZ), which is the smallest
fraction of the Brillouin Zone that is sufficient to determine all the information
on the electronic structure of the crystal. All the properties for k outside the
IBZ are obtainable via symmetry operations
Modeling from the nanoscale to the macroscale - Fall 2004
Symmetry and special points
•
•
This concept becomes particularly significant when we have to compute
properties that require an integration over the Brillouin Zone, such as energy or
electron density.
In general an integral over the BZ becomes a sum over a discrete set of states
corresponding to different k points
fi 
•
1
Nk
 f (k)
i
k
Two important issues:
• To have an accurate numerical integration the discrete set of k values has
to be dense enough to have a sufficient number of points in regions where
the integrand varies rapidly. Crucial difference between metals and
insulators. Since insulators have only filled bands, integrals can be
computed using a few well-chosen points in the BZ: special points
• Symmetry must be used to reduce the calculations to ones that involve only
k-points comprised in the IBZ
Modeling from the nanoscale to the macroscale - Fall 2004
Input parameters: K-POINTS
•
•
Set of special points in the BZ can be chosen for an efficient integration of
smooth periodic functions
Most general method has been proposed by Monkhorst and Pack (implemented
in PWscf). Uniform sets of points are chosen in reciprocal space according to
the formula:
k n ,n ,n 
1
2
3
2ni  Ni  1
 2N bi
i 1,3
i
where bi are the reciprocal space
basis vectors and Ni are the
parameters that determine the
mesh size (ni = 1,2,…,Ni)
•
•
They form a uniform mesh in reciprocal space
Can be specified in an automatic way by the program:
nk1,2,3 = N1,2,3 and k1,2,3 = 0 or 1
0 = centered grid, 1 = shifted grid
Modeling from the nanoscale to the macroscale - Fall 2004
Input parameters: K-POINTS
•
•
A sum over a uniform set of special k-points integrates exactly a periodic
function that has Fourier components that extend only to NiTi in each direction
This logic can be easily understood in one dimension:
the value of the following integral
2
I1 
 sin(k)dk  0
0
•
is given by the value of the integrand function f2(k)= sin(k) at the mid-point
k = , f2(k= )= sin()=0
If one has two integrand functions, f2(k)=A1sin(k) +A2sin(2k), the integral is given
by the sum over two points
2
3
)
0
2
2
For functions with the symmetry of the full crystal symmetry group, the density of
the special k-point mesh gives a direct measure of the accuracy with which we
compute the integrals: convergence parameter
I1 
•
A1 sin(k)  A2 sin(2k)dk  0  f2 (k1 

)  f2 (k2 
Modeling from the nanoscale to the macroscale - Fall 2004
Input parameters: K-POINTS
•
•
•
Integrals over the full Brillouin Zone (BZ) can be replaced by integrals over the
first Brillouin Zone (IBZ). In this way we need to perform summations on a
subset of the full k-point mesh in the BZ.
NOTE: All the k points of the full BZ can be obtained from the k points in the IBZ
via the symmetry operations of the crystal space group.
If we define the weight factor wk as the total number of distinguishable k points
related by symmetry to a k point in the IBZ divided by the total number of points
Nk.
With this proviso, any sum over the BZ can be written as a sum over the IBZ
that includes the appropriate weight for each point:
1
fi 
Nk
IBZ
 w f (k)
k i
k
Modeling from the nanoscale to the macroscale - Fall 2004
A few notes on “convergence”
•
•
•
The accuracy of an electronic structure calculation depends on various
approximation factors, both “physical” and numerical:
Physical approximations:
• description of nuclei/ions: all-electron vs. pseudopotential, type of
pseudopotential used
• choice of exchange and correlation functional (LDA vs. GGA. vs. hybrid
functionals or mixed schemes (this is usually done in the pseudopotential
input file)
Numerical approximations:
• The accuracy of the basis set:
• In a plane wave basis, the extension of the reciprocal space mesh
(Fourier basis)
• In a real space calculation, the density of the grid in real space
• The accuracy of the integrals in reciprocal space: density of the special kpoint mesh. Particularly important for metallic systems, where large sets of
k-points have to be used for a consistent description of the system
• The accuracy of the self-consistent solution of the K-S equations: mixing
schemes and convergence threshold
Modeling from the nanoscale to the macroscale - Fall 2004
Results
•
Calculation for a crystal of Si
• Diamond structure
• Equilibrium lattice parameter
• FCC cell with two atoms in the basis - 48 point symmetry operation with
non-symmorphic translations
a3
a2
a1
Modeling from the nanoscale to the macroscale - Fall 2004
Results
pw.x < si.scf.in > out
•
Preliminaries
•
First part: geometry
• Warning on the existence of non-symmorphic point group operations
Modeling from the nanoscale to the macroscale - Fall 2004
Results
•
More on geometry and computational parameters:
lattice parameter and volume of the cell, atoms and species, kinetic energy cutoff (ecutwfc), convergence threshold and mixing parameters, exchange and
correlation functional, real and reciprocal space basis vectors
Modeling from the nanoscale to the macroscale - Fall 2004
Results
•
Pseudopotential parameters
•
Symmetry operations and atomic coordinates
•
k-mesh for BZ integration
Modeling from the nanoscale to the macroscale - Fall 2004
Results
•
Iterate!
Modeling from the nanoscale to the macroscale - Fall 2004
Results
•
•
•
•
•
Convergence!
Electronic energies (eigenvalues of
the K-S equation) for the occupied
bands of the solid
Total energy = ground state energy
= minimum of the K-S functional
Individual contributions to the
ground state energy (given in an
alternative way, that groups
together various contributions and
uses the sum of the eigenvalues of
the K-S equations):
• band energy
• 1-electron contribution
• Hartree energy
• x-c energy
• Ewald (ionic) energy
Implicitly, we clearly know the
ground state electron density
Modeling from the nanoscale to the macroscale - Fall 2004
Results
•
•
•
From the knowledge of the electron density, we can compute all the ground
state properties of interest
Forces (in this case are zero because of symmetry)
Stress (zero, or negligible, since the calculation has been done at the
equilibrium lattice parameter)
Modeling from the nanoscale to the macroscale - Fall 2004
Results
•
Collect run-time statistics and finish the run
Modeling from the nanoscale to the macroscale - Fall 2004
Project
•
1.
2.
3.
4.
Using the input file given on your home directory on pharos.physics.ncsu.edu do
the following exercises:
Study the variation of the total ground state energy as a function of the
convergence parameters:
• ecutwfc = 2,5,8,12,18,24,32 Ry with the given k-point mesh. Discuss the
results and find the converged value for Etot.
• k-point mesh: using the converged value for ecutwfc from above, use the
option {automatic} and check convergence for various grid sizes both
centered and not:
1 1 1 0 0 0, 1 1 1 1 1 1, 2 2 2 0 0 0, 2 2 2 1 1 1, etc.
Discuss the results
Displace one of the atoms in the cell by a small amount along the diagonal of
the cell and repeat the study of point 1 above looking at the values of the forces
and stresses. Discuss the results
Using the ideal geometry, do a study of total ground state energy versus lattice
parameter (volume). The minimum of the curve will give you the theoretical
lattice parameter for Si corresponding to the pseudopotential used (GGA-PBE:
Si.pbe-rrkj.UPF)
Repeat the study of point 3 above using the LDA pseudopotential (Si.pzvbc.UPF). Discuss the results.
Modeling from the nanoscale to the macroscale - Fall 2004