Graphene with PHONON - Francesco Buonocore CMAST Website 3.0

Francesco Buonocore
ENEA - TERIN-ICT Division, HPC Lab
F. Buonocore
Go to content

Graphene with PHONON


Graphene phonons calculation using phonon (ph.x in the QUANTUM ESPRESSO package)

In the following tutorial it is shown how to calculate phonon of graphene with Pwscf & Phonon codes included in the QUANTUM ESPRESSO package installed on CRESCO.
The details on keywords are found in the QUANTUM ESPRESSO web-page and ph.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.
1) Build up the graphene structure with your preferred atomistic model editor.
2) Make a relax of graphene following the example of the relax calculation in the previous tutorials. Here we provide the relaxed atomic coordinates and unit cell parameters.
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).

3) Create your pwscf input file as follows:

      calculation = 'scf',
      restart_mode = 'from_scratch',
      pseudo_dir = './',
      outdir = './tmp',
      prefix = 'graphene',
            ibrav = 4,
                a = 2.4629,
                c = 10
              nat = 2,
             ntyp = 1,
      occupations = 'smearing',
         smearing = 'methfessel-paxton',
          degauss = 0.02,
          ecutwfc = 120,
nbnd = 8,
  mixing_mode = 'plain' ,
  conv_thr = 1.0d-10,
  mixing_beta = 0.7,
C 12.011 C_ONCV_PBE-1.0.upf
C  0.333333333  0.666666666  0.250000000
C  0.666666666  0.333333333  0.250000000
K_POINTS (automatic)
12 12 1 0 0 0
i) The convergence threshold of SCF calculation conv_thr must be very low (1.0D-10), smaller than for total energy calculation.
ii) GGA norm-conserving pseudopotentials have been used (C_ONCV_PBE-1.0.upf generated using ONCVPSP code by D. R. Hamann).
3) Launch the scf calculation using the input file above with EXE= pw .

4) Create your phonon input file as follows:

phonons of graphene on a grid      ! title line
prefix="graphene", ! must be the same than scf calculation
tr2_ph = 1d-20, ! the threshold
ldisp = .true., ! when true the calculation is on a grid
nq1=9, nq2=9, nq3=1, ! uniform q-point grid
fildyn='graphene.dyn', ! name of output file

5) Launch the phonon calculation using the input file above with EXE= ph .

6) Check the ouput files. In the same directory than input file you will find the files graphene.dyn#, with # ranging from 0 to 12. The files from 1 to 12 are the dynamical matrices on the grid of q-points. graphane.dyn0 is the list of inequivalent q-point (12 in this case).
7) Create the following input for q2r.x:

fildyn='graphene.dyn', ! the dynamical matrix from phonon
flfrc='graphene.fc', ! the interatomic force constants

8) Launch the conversion of the dynamical matrices in a supercell 9x9x1 of the reciprocal space from the phonon calculation (q) to the interatomic force constants in a supercell 9x9x1 of real space (r), using the input file above with EXE= q2r .

9) Create the input file for matdyn.x:

   flfrc='graphene.fc', ! the interatomic force constants
   flfrq='graphene.freq', ! output file with frequencies
   q_in_band_form=.true. ! band format of k-points
0       0       0 40 ! Gamma point
0       0.577350269     0 20 ! M point
-0.333333333    0.577350269     0 40 ! K point
0       0       0 20 ! Gamma point

i) Note: setting q_in_band_form=.true. the list of k-points must be given as follows:
# of high-symmetry k-points
k1(A) k2(A) k3(A) np(A-B)
k1(B) k2(B) k3(B) np(B-C)

where np(A-B) is the number of k-points along the path from A to B, and so on, A, B, ... being the high-symmetry k-points.

10) Calculate the phonon dispersion using the interatomic force constants by calculating phonons at a generic q point, using the input file above with the EXE= matdyn .

11) Generate a plot of the phonon bands. You can use plotband.x inserting graphene.disp.freq as input file. Frequency units are cm-1. Use plotband.x in /afs/<your_qe_version>.

You can plot the graphene.disp.freq.dat data file with xmgrace.
However, latest versions of QE gives the .gp file for direct plotting.

12) Create also the input file for matdyn.x:

dos=.true., ! option for dos phonon calculation
fldos='graphene.phdos', ! output file
nk1=50,nk2=50,nk3=1, grid of q-points for DOS calculation
deltaE=1.0d0 ! interval of energy

13) Calculate the vibrational density of states, using the input file above with EXE= matdyn .

14) Plot the density of states from graphene.phdos output file. Frequency units are cm-1.
Casaccia Research Center, Rome, Italy
TERIN-ICT Division, HPC Lab
Back to content
Cookies Policy