Skip to content

Structural optimization and band structure calculations using PWmat

1 Introduction

PWmat is a density functional theory (DFT) simulation software that runs on graphics processing units (GPUs) Item [1,2], which uses the plane wave pseudopotential method to support electronic structure, structural relaxation, ab initio molecular dynamics Multiple basic computational functions such as transition state search, real-time time-dependent density functional, and non adiabatic molecular dynamics. Simply put, there are three steps to using PWmat: 1. Prepare the input file; 2. Run PWmat; 3. Post processing and visualization.

Users need to prepare three types of input files: 1. Pseudo potential file: Only supports mode conservative pseudo potential and ultra soft pseudo potential. Users can download the pseudo potential package from the PWmat official website. If the structure to be calculated contains several elements, it is necessary to copy the pseudopotential files of these elements to the calculation directory; 2. Structure file: It must be in the format specified by PWmat, and users can use the structure conversion tool on the terminal Convert common crystal structure formats such as. cif and. xsf to PWmat format; 3. Parameter file: must be named etot.input, which specifies the names of the pseudopotential file and structure file, and Parameters such as truncation energy and k-point can be set.

2 Lattice optimization and atomic relaxation calculation

Different relaxation requirements are mainly achieved by modifying the parameters of atom.config and etot.input files. Example: * (Fixed) lattice optimization and atomic relaxation * (Partial) atomic relaxation * Optimization of electrified system

Step 1. Main input files:

atom.config： 

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
4
LATTICE
2.46596583        0.00000000          0.00000000
-1.23298291       2.13558905          0.00000000
0.00000000        0.00000000          6.41151115
POSITION
6    0.00000000    0.00000000    0.75000000 1 1 1
6    0.66667000    0.33333000    0.75000000 1 1 1
6    0.00000000    0.00000000    0.25000000 1 1 1
6    0.33333000    0.66667000    0.25000000 1 1 1      #Atoms can relax along three directions: a, b, and c. Yes (1)/No (0)

Note: Two dimensional materials can be optimized by fixing the lattice, which can be set in atom.config

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
4
LATTICE
2.46596583        0.00000000          0.00000000
-1.23298291       2.13558905          0.00000000
0.00000000        0.00000000          6.41151115
POSITION
6    0.00000000    0.00000000    0.75000000 1 1 1
6    0.66667000    0.33333000    0.75000000 1 1 1
6    0.00000000    0.00000000    0.25000000 1 1 1
6    0.33333000    0.66667000    0.25000000 1 1 1       #Atoms can relax along three directions: a, b, and c. Yes (1)/No (0)
STRESS_MASK
1 1 0
1 1 0
0 0 0   #  Note:Lattice optimization with fixed c-axis direction, but atoms can still relax along the c-axis

etot.input

1
2
3
4
5
6
7
8
JOB = RELAX                    #What are you going to do, such as scf/delax/nonssf
IN.PSP1=C.SG15.PBE.UPF         #pseudo potential
IN.ATOM= atom.config           #Read structural files
RELAX DETAIL= 1 100 0.01 1 0.02
VDW = DFT-D2                   #Van der Waals correction for two-dimensional materials
Ecut = 70 
Ecut2=280 
MP_N123=12 12 4 0 0 0          #K-point setting

Partial atomic relaxation：

atom.config： 

1
2
3
4
5
6
7
8
2
LATTICE
2.46596583        0.00000000          0.00000000
-1.23298291       2.13558905          0.00000000
0.00000000        0.00000000          6.41151115
POSITION
14    0.00000000    0.00000000    0.75000000 0 0 0
14    0.66667000    0.33333000    0.75000000 1 1 1
etot.input

1
2
3
4
5
6
7
8
JOB = RELAX                  
IN.PSP1=C.SG15.PBE.UPF        
IN.ATOM= atom.config          
RELAX DETAIL= 1 100 0.05      #Do not set lattice optimization parameters, default not to optimize lattice
VDW = DFT-D2                 
Ecut = 70 
Ecut2=280 
MP_N123=12 12 4 0 0 0

RELAX_DETAIL=IMTH NSTEP FORCE_TOL ISTRESS TOL_STRESS 1. IMTH represents optimization methods

  • IMTH=1(default),conjugated gradient;
  • IMTH=2,BFGS method;
  • IMTH=3,steepest decent;
  • IMTH=4,preconditioned conjugate gradient;
  • IMTH=5,Limited-memory BFGs method;
  • IMTH=6,FIRE: Fast Inertial Relaxation Engine.

  • NSTEP represents the maximum number of relaxation steps

  • FORCE_TOL represents the convergence criterion for relaxation (eV/Å)

  • ISTRESS indicates whether lattice optimization is performed

  • TOL STRESS represents the standard for lattice optimization (eV/Natom)

Step 2. Main output files:

  • final.config: Relaxed structure
  • RELAXSTEPS: RELAX per step
  • MOVEMENT: Trajectory

final.config:

alt text

RELAXSTEPS:

alt text

  1. It: The number of steps for ion optimization
  2. NEW: represents a change in the search direction of the step force; CORR: indicates that the search direction of the step force is consistent with the direction of the NEW step on it
  3. E: The total energy of the ion step system
  4. AV_F: The average force exerted on each atom in this step
  5. M_F: The maximum force exerted on all atoms in this step
  6. DE, dRho: The absolute value of the difference between the last step of the ion step that is self consistent and the previous step
  7. SCF: The number of autonomous steps included in the ion step
  8. dL: Step size of ion step movement
  9. p*F: Projection value of atomic force along the search direction

Note: If * END appears in the last step, it means that normal operation has ended

Extract column data from RELAXSSTEP and view the convergence information of the relaxation process (separated by spaces)

Method 1: Draw using the gnuplot command (enter gnuplot first to enter command-line mode) Enter gnuplot to enter command-line mode 

1
2
plot 'RELAXSSTEPS'u 2:5 w lp # View energy changes
plot 'RELAXSSTEPS'u 2:9 w lp # View maximum force variation
Method 2: Extract information from RELAXATEPS using the gawk command: 
1
2
gawk '{print $2, $5}' RELAXSSTEPS>E.dat # Energy change information saved in E.dat file
gawk '{print $2, $9}' RELAXSSTEPS>M_F.dat # Maximum force change information saved in M_F.dat file

MOVEMENT:

You can convert MOVEMENT to an. xsf format file to view the relaxed trajectory animation Command line: 

1
convert_from_config.x MOVEMENT

3 band structure calculations

  • Perform self consistent computation first, followed by non autonomous computation (job=nonsref)
  • When using HSE functional calculations, K-point parallel settings cannot be used

Step 1. Self consistent calculation (job=scf):

etot.input:

1
2
3
4
5
6
1 4
JOB =SCF
IN.PSP1= Si.SG15.PBE.UPF
IN.ATOM = atom.config
ecut= 50
MP_N123 = 9 9 9 0 0 0

After the calculation is completed, keep OUT.FERMI and copy the potential file "OUT.VR" to "IN.VR“

Step 2. __ Prepare the high symmetry K-point file IN.KPT required for non autonomous computation__:

2.1 Generate k points based on the path of high symmetry points, create a gen.kpt file, write the coordinates of high symmetry points and the number of inserted points

gen.kpt: 

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
K-path generated by Q-Studio, density:0.01 ...       #Identification line
19                                                   #The number of k-points inserted between high symmetry points X and R 
0.000000  0.000000  0.000000  G                      # the fractional coordinates of high symmetry points in reciprocal space; label
0.500000  0.000000  0.500000  X
10
0.500000  0.000000  0.500000  X
0.500000  0.250000  0.750000  W
7
0.500000  0.250000  0.750000  W
0.375000  0.375000  0.750000  K
20
0.375000  0.375000  0.750000  K
0.000000  0.000000  0.000000  G
...

Note: The gen.kpt file can be generated using Q-Studio. Steps:
(1) Import structural file (convert to primitive cell)
(2) Click on settings, symmetry, Brillouin zone path in sequence
(3) Generate path according to prompts, enter K-point spacing, click "accept" to generate gen.kpt file

2.2 Generate IN.KPT file after running the command, command line: 

1
split_kp.x gen.kpt
·IN.KPT: 
1
2
3
4
5
6
136  #Number of K points
2   0   #Fractional coordinates; weight
0.000000   0.000000   0.000000   1.00000
0.025000   0.000000   0.025000   1.00000
0.050000   0.000000   0.050000   1.00000
0.075000   0.000000   0.075000   1.00000

Prepare parameter files for non autonomous calculation:

etot.input: 

1
2
3
4
5
6
7
1 4
JOB =NONSCF
IN.PSP1= Si.SG15.PBE.UPF
IN.ATOM = atom.config
ecut= 50
IN.VR =T    #Read IN.VR file
IN.KPT = T   #Read IN.KPT file

Result post processing

Run plot band structure.x to generate bandstructure 1.txt

Command line:
plot band structure.x

alt text

After pwmat completes the nonscf calculation, if there is an out.fermi file in the calculation directory, the txt and pictures obtained after running the script will reset the Fermi level to zero. In addition, if there is a report file, you can also use the gap read script to quickly calculate the band gap