Graphane with PWSCF
Tutorials
Graphane electronic calculation using PWSCF (pw.x in the QUANTUM ESPRESSO package)
In the following tutorial it is shown how to calculate electronic structure of graphane with PWSCF code included in the QUANTUM ESPRESSO package installed on CRESCO.
The details on keywords are found in the QUANTUM ESPRESSO web-page and pw.x input file description.
You can check the best suited queue and scratch directory from the Cresco website. E.g., there the information for cresco4 and cresco6 clusters.
More information about how to use the QUANTUM ESPRESSO package on CRESCO can be found at this page.
To run a job on CRESCO6 the command line is:
bsub -n NP -q QUEUE -o lsf_file.out -e file.err qe 6.3.0 par EXE input.inp output.out
where NP is the number of processors (NP=48 is one node on CRESCO6), QUEUE=cresco6_48h24 (24h queue), cresco6_48h48 (48h queue, NP >= 48), cresco6_h144 (144h queue, NP can be less than 48), EXE = executable (e.g., pw, pp, bands, projwfc, ph).
Relax Calculation
1) Build up the graphane structure with your preferred atomistic model editor. Make note of the unit cell parameters and atomic coordinates.
2) Peform a vc-relax as done for graphene in the previous tutorial. Then prepare the relax input file graphane.pbe.relax.in as follows:
&CONTROL
calculation='relax'
title='graphane'
prefix='graphane'
verbosity='high'
restart_mode='from_scratch'
nstep=1000
iprint=1
tstress=.true.
tprnfor=.true.
disk_io='default'
pseudo_dir='/afs/enea.it/software/qu_esp/pseudo/'
outdir='./tmp'
wf_collect=.true.
tstress=.true.
tprnfor=.true.
forc_conv_thr=1.0d-4
etot_conv_thr=1.0d-5
/
&SYSTEM
ibrav = 12,
a = 2.538,
b = 2.538,
c= 20.000,
cosab=-0.500000,
nat = 4,
ntyp = 2,
ecutwfc = 60.0 ,
ecutrho = 600.0 ,
input_DFT = 'PBE' ,
/
&ELECTRONS
electron_maxstep = 100,
conv_thr = 1.0d-10 ,
mixing_mode = 'plain' ,
mixing_beta = 0.3d0 ,
/
&IONS
ion_dynamics='bfgs'
upscale=20.0
/
ATOMIC_SPECIES
C 12.011 C.pbe-van_ak.UPF
H 1.0079 H.pbe-van_ak.UPF
ATOMIC_POSITIONS angstrom
C -0.000000127 1.465548921 10.000000000
C 1.269000127 0.732774350 10.460000000
H -0.000000127 1.465548921 8.888000000
H 1.269000127 0.732774350 11.572000000
K_POINTS automatic
192 192 1 0 0 0
3) Launch the calculation using the input file above with the EXE=pw.
4) Check from graphane.pbe.relax.out that total force and stress are below the desired threshold.
They are found to be:
total stress (Ry/bohr**3) (kbar) P= 0.49
Total force = 0.000678
confirming the achived optimization of the cell and the atomic coordinates.
Self-Consistent Calculation
1) Extract the optimized atomic corrdinates from graphane.pbe.relax.out and prepare the input file graphane.pbe.scf.in for the self-consistent scf calculation as follows:
&CONTROL
calculation='scf'
title='graphane'
prefix='graphane'
verbosity='high'
restart_mode='from_scratch'
nstep=1000
iprint=1
tstress=.true.
tprnfor=.true.
disk_io='default'
pseudo_dir='/afs/enea.it/software/qu_esp/pseudo/'
outdir='./tmp'
wf_collect=.true.
tstress=.true.
tprnfor=.true.
/
&SYSTEM
ibrav = 12,
a = 2.538,
b = 2.538,
c= 20.000,
cosab=-0.500000,
nat = 4,
ntyp = 2,
ecutwfc = 60.0 ,
ecutrho = 600.0 ,
nbnd=10
input_DFT = 'PBE' ,
/
&ELECTRONS
electron_maxstep = 100,
conv_thr = 1.0d-10 ,
mixing_mode = 'plain' ,
mixing_beta = 0.3d0 ,
/
ATOMIC_SPECIES
C 12.011 C.pbe-van_ak.UPF
H 1.0079 H.pbe-van_ak.UPF
ATOMIC_POSITIONS angstrom
C -0.000000127 1.465495893 10.000196385
C 1.269000127 0.732827575 10.459793118
H -0.000000127 1.465530658 8.888518366
H 1.269000127 0.732792416 11.571492131
K_POINTS automatic
192 192 1 0 0 0
Notes:
I) nbnd is set to 10, so that 5 unoccupied bands are considered. By the way, we will get the lowest unoccupied state.
3) When the calculation is done, you can extract the total energy with the command:
tail -200 graphane.pbe.scf.out | grep !
You will get:
! total energy = -25.15033179 Ry
You can extract the PBE bandgap of graphane by typing:
tail -200 graphane.pbe.scf.out | grep high
You will find:
highest occupied, lowest unoccupied level (ev): -2.7721 0.6132
Therefore the PBE bandgap is 3.38 eV.
Density of States
1) Now we calculate the density of states. First we perform a non self-consistent calculation (nscf) based on the results of the self-consistent calculation in the previous step. In the following input file graphane.pbe.dos.nscf.in the number of bands (nbnd) and the k-point mesh are changed:
&CONTROL
calculation='nscf'
title='graphane'
prefix='graphane'
verbosity='high'
restart_mode='from_scratch'
nstep=1000
iprint=1
tstress=.true.
tprnfor=.true.
disk_io='default'
pseudo_dir='/afs/enea.it/software/qu_esp/pseudo/'
outdir='./tmp'
wf_collect=.true.
tstress=.true.
tprnfor=.true.
forc_conv_thr=1.0d-4
etot_conv_thr=1.0d-5
/
&SYSTEM
ibrav = 12,
a = 2.538,
b = 2.538,
c= 20.000,
cosab=-0.500000,
nat = 4,
ntyp = 2,
ecutwfc = 60.0 ,
ecutrho = 600.0 ,
nbnd = 48 ,
input_DFT = 'PBE' ,
occupations = 'tetrahedra'
/
&ELECTRONS
electron_maxstep = 100,
conv_thr = 1.0d-10 ,
mixing_mode = 'plain' ,
mixing_beta = 0.3d0 ,
/
ATOMIC_SPECIES
C 12.011 C.pbe-van_ak.UPF
H 1.0079 H.pbe-van_ak.UPF
ATOMIC_POSITIONS angstrom
C -0.000000127 1.465495893 10.000196385
C 1.269000127 0.732827575 10.459793118
H -0.000000127 1.465530658 8.888518366
H 1.269000127 0.732792416 11.571492131
K_POINTS automatic
384 384 1 0 0 0
Notes:
I) nbnd is set to 48, so a large number of occupied bands is taken into account.
II) A dense mesh of k-points is recommended.
III) occupations = 'tetrahedra' is best suited for DOS calculations.
2) Launch the nscf calculation using the input file above with EXE= pw .
3) To calculate the DOS here we use projwfc.x in the QUANTUM ESPRESSO package that will allow, in general, to calculate also the projected density of states. The input file graphane.pbe.proj.in is:
&projwfc
prefix='graphane'
outdir='./tmp'
ngauss=0
degauss=0.04
Emin=-20.0, Emax=20.0, DeltaE=0.1
filproj='graphane.pbe.proj'
filpdos='graphane.pbe.pdos'
/
4) Launch with EXE= projwfc .
Notes:
i) We can launch the job on a different number of processor (from 32 to 16) because we set the keyword wf_collect=.true. . Otherwise, the same number of processors would be mandatory. However, it can be useful to decrease the number of processors when the "too many processors" error is found.
The output files will be the total DOS graphane.pbe.pdos.pdos_tot and the projected DOS on the atomic orbitals graphane.pbe.pdos.pdos_atm#*(C)_wfc#*(*) and graphane.pbe.pdos.pdos_atm#*(H)_wfc#*(*).
You can plot the total DOS with xmgrace (the highest occupied level has been shifted to zero).
Band Structure
1) Now we calculate the band structure. We perform again a non self-consistent calculation (nscf) where k-point mesh is substituted with paths along the lines connecting the high symmetry points in the Brillouin zone. The input file graphene.pbe.bands.nscf.in follows:
&CONTROL
calculation='bands'
title='graphane'
prefix='graphane'
verbosity='high'
restart_mode='from_scratch'
nstep=1000
iprint=1
tstress=.true.
tprnfor=.true.
disk_io='default'
pseudo_dir='/afs/enea.it/software/qu_esp/pseudo/'
outdir='./tmp'
wf_collect=.true.
tstress=.true.
tprnfor=.true.
forc_conv_thr=1.0d-4
etot_conv_thr=1.0d-5
/
&SYSTEM
ibrav = 12,
a = 2.538,
b = 2.538,
c= 20.000,
cosab=-0.500000,
nat = 4,
ntyp = 2,
ecutwfc = 60.0 ,
ecutrho = 600.0 ,
nbnd = 48 ,
input_DFT = 'PBE' ,
/
&ELECTRONS
electron_maxstep = 100,
conv_thr = 1.0d-10 ,
mixing_mode = 'plain' ,
mixing_beta = 0.3d0 ,
/
ATOMIC_SPECIES
C 12.011 C.pbe-van_ak.UPF
H 1.0079 H.pbe-van_ak.UPF
ATOMIC_POSITIONS angstrom
C -0.000000127 1.465495893 10.000196385
C 1.269000127 0.732827575 10.459793118
H -0.000000127 1.465530658 8.888518366
H 1.269000127 0.732792416 11.571492131
K_POINTS crystal_b
4
-0.333333333 0.666666667 0.000000000 64 ! K
0.0000000000 0.000000000 0.000000000 64 ! G
0.000000000 0.500000000 0.000000000 64 ! M
-0.333333333 0.666666667 0.000000000 64 ! K
2) Launch the bands nscf calculation using the input file above with EXE= pw .
3) Now, the bands are to be organized for plotting with the bands.x program. Create the following input file graphane.pbe.bands.in:
&bands
outdir = './tmp'
prefix = 'graphane'
filband = 'graphane.bands'
/
4) Launch the bands calculation using the input file above with EXE= bands .
5) Plot the bands from the graphane.gnu output file with gnuplot or xmgrace (energy units are eV, the highest occupied level has been shifted to zero).
