VASP Computing Bulk Modulus of Au

VASP: Computing Bulk Modulus of Au

Modified by Yanming Wang (Sep, 2015)

Input files

Here we give an example of how to use VASP to compute the bulk modulus of LDA-Au. We performed this calculation on MC-CC in serial model in the ~/Codes/VASP/runs/Au/LDA/perfect.21x21x21 directory. This directory contains the following files.

INCAR

PREC = High

ISTART = 0 
ICHARG = 2  
ISMEAR = 1
SIGMA = 0.1

EDIFF  = 1E-09  
NELM = 40      
ENMAX = 500     
ENCUT = 500    

ISIF = 2 
NSW = 100  
IBRION = 2  

KPOINTS

21x21x21
0        0 = automatic generation of k-points
Monkhorst
21 21 21
0 0 0

POSCAR

POSCAR for FCC Au (created manually)
4.068 
1.0  0    0
0    1.0  0
0    0    1.0
4
selective dynamics
direct  (real coordinates r)
0       0       0       T T T
0.5000  0.5000  0       T T T
0       0.5000  0.5000  T T T
0.5000  0       0.5000  T T T

POTCAR A POTCAR file that specifies the pseudopotential is required to be in the same directory. In this example, we will present the results from two calculations with US-LDA PP and PAW-LDA PP.

Run VASP

We can write one .pbs script to submit the VASP simulation request to the cluster.

Using MC2 cluster as an example, we provide the following vasp.pbs for the calculation of equilibrium lattice constant, cohesive energy and bulk modulus of Au. The goal of this script is to run vasp in iterations with different lattice constants.

#!/bin/bash
#PBS -N vasp.Au
#PBS -j oe
#PBS -l nodes=1:ppn=8,walltime=48:00:00
#PBS -V

### ---------------------------------------
### BEGINNING OF EXECUTION
### ---------------------------------------

echo The master node of this job is `hostname`
echo The working directory is `echo $PBS_O_WORKDIR`
echo This job runs on the following nodes:
echo `cat $PBS_NODEFILE`

ncpu=`cat $PBS_NODEFILE | wc -w`
echo "Number of processors = $ncpu "

### end of information preamble

cd $PBS_O_WORKDIR

echo $PWD

cd $PBS_O_WORKDIR

module load vasp

rm WAVECAR
for a in 4.064 4.065 4.066 4.067 4.068 4.069 4.070 4.071 4.072
do
cat > POSCAR << FIN
POSCAR for FCC Au (created manually)
   $a 
   1.0  0    0
   0    1.0  0
   0    0    1.0
   4
selective dynamics 
direct  (real coordinates r)
   0    0    0    T T T
   0.5  0.5  0    T T T
   0    0.5  0.5  T T T
   0.5  0    0.5  T T T
FIN

echo "a=$a  ncpu=$ncpu"

cmd="mpiexec -np $ncpu vasp"
$cmd >& vasp.log

E=`tail -1 OSZICAR`
echo $a $E | sed -s 's/F=//; s/E0=//; s/d E =//;' >> Elatt.B.dat

p=`grep pressure OUTCAR | cut -b 25-34`
echo $a $p >> platt.B.dat

done

To submit the job, make sure POTCAR, INCAR and KPOINTS are placed in the same folder of your PBS script, and then type the following command.

qsub vasp.pbs

Analyze data

After the simulation is finished, it will create data files Elatt.B.dat and platt.B.dat.

For simulation with US-LDA PP, the result files should be similar to Elatt.B.dat

4.064 1  -.17546329E+02  -.17546306E+02 -.175463E+02
4.065 1  -.17546473E+02  -.17546449E+02 -.175465E+02
4.066 1  -.17546574E+02  -.17546550E+02 -.175466E+02
4.067 1  -.17546633E+02  -.17546607E+02 -.175466E+02
4.068 1  -.17546649E+02  -.17546622E+02 -.175466E+02
4.069 1  -.17546622E+02  -.17546594E+02 -.175466E+02
4.070 1  -.17546553E+02  -.17546524E+02 -.175466E+02
4.071 1  -.17546441E+02  -.17546412E+02 -.175464E+02
4.072 1  -.17546288E+02  -.17546257E+02 -.175463E+02

platt.B.dat

4.064 5.34 5.34 5.34
4.065 3.95 3.95 3.95
4.066 2.56 2.56 2.56
4.067 1.18 1.18 1.18
4.068 -0.20 -0.20 -0.20
4.069 -1.57 -1.57 -1.57
4.070 -2.93 -2.93 -2.93
4.071 -4.29 -4.29 -4.29
4.072 -5.63 -5.64 -5.64

Launch Matlab and run the following functions fit_a0EB.m and fit_a0B.m,

[a0, Ecoh, B] = fit_a0EB('Elatt.B.dat',4,4);
fit_a0B ('platt.B.dat');

The first line fits the energy data to a quadratic curve and computes the equilibrium lattice constant, cohesive energy and bulk modulus. The second line fits the pressure data to a linear curve and computes the equilibrium lattice constant and bulk modulus. For simulation with US-LDA PP, the result is a0 = 4.068 $\mathrm {\AA}$, Ecoh = -4.39 eV, B = 186 GPa.