- The added line is THIS COLOR.
- The deleted line is THIS COLOR.
[started by 2014.12.3]
**Prerequisites [#pre]
- We will install and run the program in our own home directory.
- Therefore it's better to prepare a directory specific for quantum espresso calculation, e.g. "espresso"
mkdir ~/espresso
- Download quantum espresso package from (for version 5.0.2):
http://qe-forge.org/gf/download/frsrelease/116/403/espresso-5.0.2.tar.gz
- If you use quantum espresso(v.5.0.2), please download the patch (e.g. as "espresso-5.0.2-5.0.3.diff"):
http://qe-forge.org/gf/project/q-e/frs/?action=FrsReleaseView&release_id=128
**How to Install [#install]
- We assume that espresso-5.0.2.tar.gz is already in the "espresso" directory
- Login to FLEX and then TUBE system, e.g.
ssh tube50
- Run the following command under the "espresso" directory (to make sure, don't forget execute "cd ~/espresso" firstly):
tar -xzvf espresso-5.0.2.tar.gz
mv espresso-5.0.2 00main
cd 00main
patch -p1 < espresso-5.0.2-5.0.3.diff
You will be asked for input file to be patched, neglect it ans press enter some times, then
./configure --with-scalapack=NO
make all
(note that we don't use scalapack in TUBE, therefore we set the option as "NO", otherwise the compilation will be failed)
**Structure of the input data [#input]
- The structure of the input file of quantum espresso is
&CONTROL
/
&SYSTEM
/
&ELECTRONS
/
[ &IONS
...
/ ]
[ &CELL
...
/ ]
ATOMIC_SPECIES
X Mass_X PseudoPot_X
Y Mass_Y PseudoPot_Y
Z Mass_Z PseudoPot_Z
ATOMIC_POSITIONS { alat | bohr | crystal | angstrom }
X 0.0 0.0 0.0 {if_pos(1) if_pos(2) if_pos(3)}
Y 0.5 0.0 0.0
Z O.0 0.2 0.2
K_POINTS { tpiba | automatic | crystal | gamma | tpiba_b | crystal_b | tpiba_c | crystal_c }
if (gamma)
nothing to read
if (automatic)
nk1, nk2, nk3, k1, k2, k3
if (not automatic)
nks
xk_x, xk_y, xk_z, wk
[ CELL_PARAMETERS { alat | bohr | angstrom }
v1(1) v1(2) v1(3)
v2(1) v2(2) v2(3)
v3(1) v3(2) v3(3) ]
[ OCCUPATIONS
f_inp1(1) f_inp1(2) f_inp1(3) ... f_inp1(10)
f_inp1(11) f_inp1(12) ... f_inp1(nbnd)
[ f_inp2(1) f_inp2(2) f_inp2(3) ... f_inp2(10)
f_inp2(11) f_inp2(12) ... f_inp2(nbnd) ] ]
[ CONSTRAINTS
nconstr { constr_tol }
constr_type(.) constr(1,.) constr(2,.) [ constr(3,.) constr(4,.) ] { constr_target(.) } ]
- These notations mean: Input data format: { } = optional, [ ] = it depends, | = or.
** Compile file in quantum espresso: [#compile]
- pw.x:
calculates electronic structure, structural optimization, molecular dynamics, barriers with NEB.
- ph.x:
calculates phonon frequencies and displacement patterns, dielectric tensors, effective charges (uses data produced by pw.x).
- pp.x:
extracts the specified data from files produced by pw.x, prepare
data for plotting by writing them into formats that can be read by
several plotting programs.
- bands.x:
extracts and reorders eigenvalues from files produced by pw.x for band structure plotting.
- projwfc.x:
calculates projections of wavefunction over atomic orbitals,performs L?owdin population analysis and calculates projected density of states. These can be summed using auxiliary code sumpdos.x
- plotrho.x:
produces PostScript 2-d contour plots.
- plotband.x:
reads the output of bands.x, produces band structure PostScript plots.
- dos.x:
calculates electronic Density of States (DOS).
** Where is the pseudopotentials ?: [#pseudo]
- http://www.quantum-espresso.org/pseudopotentials/
** Where can we find the Brillouin zone and High symmetry points? [#p52c2f19]
- Bilbao Crystallographic Server:
http://www.cryst.ehu.es/cryst/get_kvec.html
** Some Conversions: [#q42f9c19]
- These are some helpful conversions for calculation in quantum espresso
1 bohr = 0.529177249 angstrom
1 Rydberg (R ) = 13.6056981 eV
1 eV =1.60217733 x 10-19 Joules
**Graphene [#Graphene]
- A single layer graphene is shown as following figure:
#ref(graphene.png)
- The electron structure of graphene is : $1s^2$ $2s^2$ $2p^4$ and its radial wave function for each orbital is :
#ref(c-2.png)
- We use the following input file as an example in this calculation:
&CONTROL
calculation = 'scf',
restart_mode = 'from_scratch',
pseudo_dir = 'pseudo/',
outdir = 'tmp/',
prefix = 'graphene',
/
&SYSTEM
ibrav = 4,
a = 2.4623,
c = 10
nat = 2,
ntyp = 1,
occupations = 'smearing',
smearing = 'methfessel-paxton',
degauss = 0.02,
ecutwfc = 60,
ecutrho = 720,
nbnd = 8,
/
&ELECTRONS
conv_thr = 1.0d-10,
mixing_mode = 'plain',
mixing_beta = 0.7,
diagonalization = 'cg',
/
ATOMIC_SPECIES
C 12.0107 C.pbe-rrkjus.UPF
ATOMIC_POSITIONS {crystal}
C 0.333333333 0.666666666 0.500000000
C 0.666666666 0.333333333 0.500000000
K_POINTS {automatic}
42 42 1 0 0 0
** The meaning of each command [#command]
- &control:
declares the control block.
- calculation = 'scf':
tells PWSCF that this will be a self-consistent field calculation.
- restart_mode = 'from_scratch':
declares that we will be generation a new structure.
- pseudo_dir = 'pseudo/':
defines the location of the directory where you store the pseudo-potentials
- outdir='tmp/':
defines the location of the temporary files. This should always be a local scratch disk so that large I/O operations do not occur across the network.
- prefix='graphene':
declares the file name prefix to be used for temporary files.
- / :
denotes the end of a block.
- &SYSTEM:
declares the system block.
- ibrav:
gives the crystal system. ibrav=4 : refers to the hexagonal lattice.
- nat:
number of atoms (each individual unique atom). Note that graphene for graphene nat=2 because of quantum berry phase.
- ntyp:
number of types of atoms
- ecutwfc:
Energy cutoff for pseudo-potentials. This one is important; we change this to find stable structure.
- Occupations, Smearing, deguass:
these keywords are particular details for the Brillioun zone integration for metals. Since there is a discontinuity of the occupation number for the bands around the Fermi energy, total energy with respect to the number of k-points converges very slowly. Adding electronic temperature (degauss) smooth out the abrupt change of the occupation number and as a result total energy converges with fewer number of k-points.
- nbnd:
number of electronic states (bands) to be calculated.Note that in spin-polarized calculations the number of k-point, not the number of bands per k-point, is doubled
- Atomic species declaration.
- After the keyword ATOMIC_SPECIES, for each ntyp enter:
atomic symbol atomic weight pseudo-potential
- Atomic positions:
after the keyword ATOMIC_POSITIONS, for each nat enter
atomic symbol x y z
where x,y,z are given as fractional coordinates of the conventional cell.
- k-point selection:
after the keywork K_POINTS, 'automatic' tells PWSCF to automatically generate a k-point grid.
The format of the next line is
nkx nky nkz offx offy offz
where nk* is the number of intervals in a direction and off* is the offset of the origin of the grid.
- The more details can be found on this website: http://www.quantum-espresso.org/wp-content/uploads/Doc/INPUT_PW.html
**How to Run [#run]
- Create a directory for running the program inside the espresso directory, for example:
cd ~/espresso
mkdir graphene
- Under "graphene" directory, put the necessary input files and pseudopotential files needed for calculation. Therefore, at least in the graphene directory we should have
++ an executable file "pw.x" which is the main program of quantum espresso
++ a file "input.in"
++ a directory "pseudo" which consists of pseudopotential files
++ a temporary directory "tmp" for calculation.
- The calculation can be performed by the following command inside graphene directory:
cd ~/espresso/graphene
nohup ./pw.x < input.in > output.out &
(we use "nohup" command for running the calculation as a background process)
- To check the calculation is running, we can see the iteration process by following command:
grep iteration output.out
#ref(slg.bands.png)
- Some important output about this calculation:
bravais-lattice index = 4
lattice parameter (alat) = 4.6090 a.u.
unit-cell volume = 347.6569 (a.u.)^3
number of atoms/cell = 2
number of atomic types = 1
number of electrons = 8.00
number of Kohn-Sham states= 8
kinetic-energy cutoff = 40.0000 Ry
charge density cutoff = 480.0000 Ry
Exchange-correlation = SLA PZ NOGX NOGC ( 1 1 0 0 0)
EXX-fraction = 0.00
Largest allocated arrays est. size (Mb) dimensions
Kohn-Sham Wavefunctions 0.03 Mb ( 264, 8)
NL pseudopotentials 0.06 Mb ( 264, 16)
Each V/rho on FFT grid 0.45 Mb ( 29808)
Each G-vector array 0.08 Mb ( 10283)
G-vector shells 0.02 Mb ( 2599)
Largest temporary arrays est. size (Mb) dimensions
Each subspace H/S matrix 0.00 Mb ( 8, 8)
Each <psi_i|beta_j> matrix 0.00 Mb ( 16, 8)
** Example silicon [#si]
- The crystal shape of Silicon is plotted as followig:
#ref(si.png)
- let's make a simple example to reproduce silicon electronic structure. All calculations and results can be found in:
/home/students/pourya/espresso/silicon
- The structure of si is shown in si.eps which was ploted by xcrysden software.I will explain how to plot at this program later.
&control
prefix='silicon',
pseudo_dir='/home/pourya/Desktop/espresso/pseudo/'
outdir = '/home/pourya/tmp/',
/
&system
ibrav= 2, celldm(1) =10.2, nat= 2, ntyp= 1,
ecutwfc = 12.0,
/
&electrons
/
ATOMIC_SPECIES
Si 28.086 Si.pz-vbc.UPF
ATOMIC_POSITIONS
Si 0.00 0.00 0.00
Si 0.25 0.25 0.25
K_POINTS
2
0.25 0.25 0.75 3.0
0.25 0.25 0.25 1.0
- In order to compile this code, you need type following command:
/pw.x < si.scf.in > si.scf.out
- Here the cutoff energy is consider 12.0 R; but, the total energy versus different cutoff energy is plotted in figure si_cut.pdf in order to show how total energy convergence by different cuttoff value.
#ref(si_cut.png)
- In order to find band at high-symmetry point in silicon please compile the si.nscf.in. The symmetry point is add manually.
/ pw.x < si.nscf.in> si.nscf.out
- Now, the band structure of silicon will be calculated by following command. Please note that the k-point is add manually. It can be calculated in the BZ along high-symmetry points direction as defined by in k-point-path.
/ pw.x < si.bands.in> si.bands.out
- Collect the band results for plotting:
/ bands.x < bands.in > bands.out
- Therefore, there will be 2 files bands.out and bands.dat. for plotting:
/ plotband.x
input file > bands.dat
Range: -5.6680 16.4950eV Emin, Emax > -6.0 10.0
output file (xmgr) > si.bands.xmgr
output file (ps) > si.bands.ps
Efermi > 6.337
deltaE, reference E (for tics) 1.0, 6.337
- hence, you can see the band structure of silicon
/evince si.bands.ps
#ref(si.bands.png)
- Charge density plot for silicon plane:
/pp.x < si.pp_rho.in
/plotrho.x
input file > si.rho.dat
output file > si.rho.ps
Logarithmic scale (y/n)? > n
/evince si.rho.ps
- Finally, the general information about electron structure of silicon is produced.
- Some important outputs about this calculation:
bravais-lattice index = 2
lattice parameter (alat) = 10.2000 a.u.
unit-cell volume = 265.3020 (a.u.)^3
number of atoms/cell = 2
number of atomic types = 1
number of electrons = 8.00
number of Kohn-Sham states= 8
kinetic-energy cutoff = 12.0000 Ry
charge density cutoff = 48.0000 Ry
Exchange-correlation = SLA PZ NOGX NOGC ( 1 1 0 0 0)
EXX-fraction = 0.00
Largest allocated arrays est. size (Mb) dimensions
Kohn-Sham Wavefunctions 0.02 Mb ( 200, 8)
NL pseudopotentials 0.02 Mb ( 200, 8)
Each V/rho on FFT grid 0.05 Mb ( 3375)
Each G-vector array 0.01 Mb ( 1459)
G-vector shells 0.00 Mb ( 43)
Largest temporary arrays est. size (Mb) dimensions
Auxiliary wavefunctions 0.10 Mb ( 200, 32)
Each subspace H/S matrix 0.02 Mb ( 32, 32)
Each <psi_i|beta_j> matrix 0.00 Mb ( 8, 8)
** Charge Density of Silicon [#si_ph]
#ref(si.charge0012.png)
#ref(si.charge001.png)
- The charge density of silicon can be calculated by QE as following commands the source file can be found in :
- pw.x < si.scf.in > si.scf.out
- pp.x < si.pp.in > si.pp.out
** Phonon [#si_ph]
- One LA phonon mode and one LO phonon mode is shown at $\Gamma$ point. The vectors indicate of movement direction.
#ref(2.png)
- LO Mode :
#ref(4.png)
![[PukiWiki]](image/c60.png "[PukiWiki]")
