![[PukiWiki]](image/c60.png "[PukiWiki]")
[started by Pourya 2014.12.3] **Introduction [#intro] This is a Quantum Espresso tutorial in Saito Lab. We also learn from other tutorial provided by Saito(Susumu)-sensei's group at Titech, Tokyo: This is a Quantum Espresso tutorial in Saito Lab. We can also learn from other tutorial provided by Saito(Susumu)-sensei's group at Titech, Tokyo: http://www.stat.phys.titech.ac.jp/SATL_qe_tutorial/ (in Japanese!) (Very rough) Translation is available here: #ref(QE-Tutorial.pdf) **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)
