VASP Computing Bulk Modulus of Si: Difference between revisions

From Micro and Nano Mechanics Group
Jump to navigation Jump to search
 
(16 intermediate revisions by the same user not shown)
Line 3: Line 3:
VASP: Computing Bulk Modulus of Si</STRONG></font></P>
VASP: Computing Bulk Modulus of Si</STRONG></font></P>
<DIV>
<DIV>
<P ALIGN="CENTER"><STRONG> </STRONG></P>
<P ALIGN="CENTER"><STRONG>Adriano Sanchez，Yanming Wang and Wei Cai</STRONG></P>
</DIV>
</DIV>




===Input files===
===Input files===
Here we give an example of how to use VASP to compute the equilibrium lattice constant, cohesive energy and bulk modulus of Si. We performed this calculation on <tt>MC2</tt> in parallel mode in the <tt>~/Codes/VASP/runs/Si/LDA/perfect.21x21x21</tt> directory. This directory should contain the following files.
Under construction!
Here we give an example of how to use VASP to compute the bulk modulus of LDA-Si. We performed this calculation on <tt>MC-CC</tt> in serial model in the <tt>~/Codes/VASP/runs/Si/LDA/perfect.21x21x21</tt> directory. This directory contains the following files.


'''INCAR'''
'''INCAR'''


<pre>
<pre>

ENCUT = 400
ISMEAR = 1
PREC = High

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

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

ISIF = 2
NSW = 100
IBRION = 2
</pre>
</pre>


Line 29: Line 41:
</pre>
</pre>


'''POSCAR example for Si'''
'''POSCAR'''


Here's one example for Si, which will be automatically rewritten every iteration.
<pre>
<pre>
POSCAR for DC Si (created manually)
POSCAR for DC Si (created manually)
5.431 universal scaling factor
???
0 0.5 0.5
0.0 0.5 0.5
0.5 0 0.5
0.5 0.0 0.5
0.5 0.5 0
0.5 0.5 0.0
2
1
Direct
Cartesian (real coordinates r)
-0.125 -0.125 -0.125
0 0 0
0.125 0.125 0.125
</pre>
</pre>


'''POTCAR'''
To do this calculation, you also need to put the LDA pseudopotential file as <tt>POTCAR</tt> in this directory.


It requires to put one pseudopotential (PP) file <tt>POTCAR</tt> in this directory. There are different choices of PP and these potential files are provided by VASP. In this example, we choose the ultrasoft-LDA PP. The following is an article explains the differences of different PPs.
===Run VASP===


https://www.vasp.at/vasp-workshop/slides/pseudoppdatabase.pdf
Now we are ready to run


===Run VASP===
vasp


We can write one .pbs script to submit the VASP simulation request to the cluster.
To compute the equilibrium lattice constant, cohesive energy and bulk modulus, we use the following script <tt>auto.B.serial</tt> to run <tt>vasp</tt> repeated with different lattice constants.

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


<pre>
<pre>
#!/bin/bash
#!/bin/bash
#PBS -N vasp.Si
rm WAVECAR
#PBS -j oe
for a in 5.375 5.38 5.385 5.39 5.395 5.40 5.405 5.41 5.415 5.42 5.425 5.430 5.435
#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 5.382 5.383 5.384 5.385 5.386 5.387 5.388 5.389 5.390 5.391 5.392 5.393 5.394 5.395 5.396 5.397
do
do
cat > POSCAR << FIN
cat > POSCAR << FIN
cubic diamond
POSCAR for FCC Si (created manually)
$a universal scaling factor
$a
0.5 0.5 0.0
0.0 0.5 0.5
0.0 0.5 0.5
0.5 0.0 0.5
0.5 0.0 0.5
0.5 0.5 0.0
2
2
Direct
cart
0.00 0.00 0.00
-0.125 -0.125 -0.125
0.25 0.25 0.25
0.125 0.125 0.125
FIN
FIN


echo "a=$a ncpu=$ncpu"


cmd="mpiexec -np $ncpu vasp"
$MPIRUN -np $nprocs -hostfile $PBS_NODEFILE $VASP5
$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
"sielastic.pbs" [dos] 49L, 814C echo "a=$a"
./vasp


E=`tail -1 OSZICAR`
E=`tail -1 OSZICAR`
Line 89: Line 126:
</pre>
</pre>


To submit the job, prepare the '''vasp.pbs''' file with the above content and put it in the same folder that contains the '''INCAR''', '''POTCAR''', '''KPOINTS''' files. From this folder, type

qsub vasp.pbs


===Analyze data===
===Analyze data===


After running it as <tt>./auto.B.serial</tt>, it will create data files <tt>Elatt.B.dat</tt> and <tt>platt.B.dat</tt>.
When the simulation is finished, it will create data files <tt>Elatt.B.dat</tt> and <tt>platt.B.dat</tt>.


When choosing US-LDA PP, your simulation data should be similar to
Launch <tt>octave</tt> and run the following functions [[media:fit_a0EB.m.txt | fit_a0EB.m]] and [[media:fit_a0B.m.txt | fit_a0B.m]],


fit_a0EB('Elatt.B.dat');
<tt>Elatt.B.dat</tt>
<pre>
fit_a0B ('platt.B.dat');
5.382 1 -.11953275E+02 -.11953296E+02 -.119533E+02
5.383 1 -.11953325E+02 -.11953346E+02 -.119533E+02
5.384 1 -.11953368E+02 -.11953389E+02 -.119534E+02
5.385 1 -.11953403E+02 -.11953424E+02 -.119534E+02
5.386 1 -.11953431E+02 -.11953452E+02 -.119534E+02
5.387 1 -.11953451E+02 -.11953473E+02 -.119535E+02
5.388 1 -.11953465E+02 -.11953486E+02 -.119535E+02
5.389 1 -.11953471E+02 -.11953492E+02 -.119535E+02
5.390 1 -.11953470E+02 -.11953491E+02 -.119535E+02
5.391 1 -.11953462E+02 -.11953483E+02 -.119535E+02
5.392 1 -.11953446E+02 -.11953468E+02 -.119534E+02
5.393 1 -.11953424E+02 -.11953445E+02 -.119534E+02
5.394 1 -.11953394E+02 -.11953415E+02 -.119534E+02
5.395 1 -.11953357E+02 -.11953378E+02 -.119534E+02
5.396 1 -.11953313E+02 -.11953334E+02 -.119533E+02
5.397 1 -.11953261E+02 -.11953283E+02 -.119533E+02
</pre>


<tt>platt.B.dat</tt>
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. In this example, the result is a0 = 5.389 angstrom, Ecoh = -5.977 eV, B = 95 GPa.
<pre>
5.382 3.89
5.383 3.35
5.384 2.81
5.385 2.28
5.386 1.74
5.387 1.20
5.388 0.67
5.389 0.14
5.390 -0.39
5.391 -0.92
5.392 -1.45
5.393 -1.98
5.394 -2.50
5.395 -3.03
5.396 -3.55
5.397 -4.07
</pre>


Launch <tt>Matlab</tt> and run the following functions [[media:fit_a0EB.m.txt | fit_a0EB.m]] and [[media:fit_a0B.m.txt | fit_a0B.m]],
===Parallel computation===

To run <tt>vasp</tt> in parallel, you need to submit [[media:vasp.pbs.txt | vasp.pbs]] as

qsub vasp.pbs


[a0, Ecoh, B] = fit_a0EB('Elatt.B.dat',2,8);
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. In this example with US-LDA PP, the results are a0 = 5.389 angstrom, Ecoh = -5.977 eV, B = 95 GPa.
You will need the following two files to do this calculation in parallel on SU-AHPCRC:
[[media:auto.B.par.Au.txt | auto.B.par]] and
[[media:B.pbs.Au.txt | B.pbs]].

Latest revision as of 08:10, 4 September 2015

VASP: Computing Bulk Modulus of Si

Adriano Sanchez，Yanming Wang and Wei Cai


Input files

Here we give an example of how to use VASP to compute the equilibrium lattice constant, cohesive energy and bulk modulus of Si. We performed this calculation on MC2 in parallel mode in the ~/Codes/VASP/runs/Si/LDA/perfect.21x21x21 directory. This directory should contain 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

Here's one example for Si, which will be automatically rewritten every iteration. 

POSCAR for DC Si (created manually)
5.431           universal scaling factor 
0.0    0.5     0.5
0.5    0.0     0.5
0.5    0.5     0.0
2
Direct
-0.125 -0.125 -0.125
 0.125  0.125  0.125

POTCAR

It requires to put one pseudopotential (PP) file POTCAR in this directory. There are different choices of PP and these potential files are provided by VASP. In this example, we choose the ultrasoft-LDA PP. The following is an article explains the differences of different PPs. 

https://www.vasp.at/vasp-workshop/slides/pseudoppdatabase.pdf

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. The goal of this script is to run vasp in iterations with different lattice constants. 

#!/bin/bash
#PBS -N vasp.Si
#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 5.382 5.383 5.384 5.385 5.386 5.387 5.388 5.389 5.390 5.391 5.392 5.393 5.394 5.395 5.396 5.397
do
cat > POSCAR << FIN
cubic diamond
$a             universal scaling factor
0.0    0.5     0.5
0.5    0.0     0.5
0.5    0.5     0.0
2
Direct
 -0.125 -0.125 -0.125
  0.125  0.125  0.125
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, prepare the vasp.pbs file with the above content and put it in the same folder that contains the INCAR, POTCAR, KPOINTS files. From this folder, type 

qsub vasp.pbs

Analyze data

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

When choosing US-LDA PP, your simulation data should be similar to

Elatt.B.dat 

5.382 1  -.11953275E+02  -.11953296E+02 -.119533E+02
5.383 1  -.11953325E+02  -.11953346E+02 -.119533E+02
5.384 1  -.11953368E+02  -.11953389E+02 -.119534E+02
5.385 1  -.11953403E+02  -.11953424E+02 -.119534E+02
5.386 1  -.11953431E+02  -.11953452E+02 -.119534E+02
5.387 1  -.11953451E+02  -.11953473E+02 -.119535E+02
5.388 1  -.11953465E+02  -.11953486E+02 -.119535E+02
5.389 1  -.11953471E+02  -.11953492E+02 -.119535E+02
5.390 1  -.11953470E+02  -.11953491E+02 -.119535E+02
5.391 1  -.11953462E+02  -.11953483E+02 -.119535E+02
5.392 1  -.11953446E+02  -.11953468E+02 -.119534E+02
5.393 1  -.11953424E+02  -.11953445E+02 -.119534E+02
5.394 1  -.11953394E+02  -.11953415E+02 -.119534E+02
5.395 1  -.11953357E+02  -.11953378E+02 -.119534E+02
5.396 1  -.11953313E+02  -.11953334E+02 -.119533E+02
5.397 1  -.11953261E+02  -.11953283E+02 -.119533E+02

platt.B.dat 

5.382 3.89
5.383 3.35
5.384 2.81
5.385 2.28
5.386 1.74
5.387 1.20
5.388 0.67
5.389 0.14
5.390 -0.39
5.391 -0.92
5.392 -1.45
5.393 -1.98
5.394 -2.50
5.395 -3.03
5.396 -3.55
5.397 -4.07

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

[a0, Ecoh, B] = fit_a0EB('Elatt.B.dat',2,8);
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. In this example with US-LDA PP, the results are a0 = 5.389 angstrom, Ecoh = -5.977 eV, B = 95 GPa.

Retrieved from "http://micro.stanford.edu/mediawiki/index.php?title=VASP_Computing_Bulk_Modulus_of_Si&oldid=6282"