VASP Computing Elastic Constants of Si: Difference between revisions

Revision as of 05:56, 2 September 2015

VASP: Computing Elastic Constants of Si

Adriano Sanchez，Yanming Wang and Wei Cai

Input files

Here we give an example of how to use VASP to compute the elastic constants C11, C12 and C44 of Si. We performed this calculation on MC2 in PARALLEL mode in ~/Codes/VASP/runs/Si_bulk_C11 and ~/Codes/VASP/runs/Si_bulk_C44 two directories. Each 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

Here's an example for Si C11 and C12 calculation, which will be changed later in PBS script every iteration.

POSCAR for DC Si
5.402856        accurate equilimrium lattice constant
1.003  0.0 0.0  first  Bravais lattice vector
0.0 1.0 0.0     second Bravais lattice vector
0.0 0.0 1.0     third  Bravais lattice vector
8               number of atoms per species
selective dynamics
direct           direct or cart (only first letter is significant)
0.00000  0.00000  0.00000 T T T
0.50000  0.50000  0.00000 T T T
0.00000  0.50000  0.50000 T T T
0.50000  0.00000  0.50000 T T T
0.25000  0.25000  0.25000 T T T
0.75000  0.75000  0.25000 T T T
0.25000  0.75000  0.75000 T T T
0.75000  0.25000  0.75000 T T T

Here's an example for Si C44 calculation, which will be changed later in PBS script every iteration.

POSCAR for DC Si
5.402856         accurat equilibrium lattice constant
1.000  1.000 0.0 first  Bravais lattice vector
-1.0 1.0 0.0     second Bravais lattice vector
0.0 0.0 1.0      third  Bravais lattice vector
16               number of atoms per species
selective dynamics
direct           direct or cart (only first letter is significant)
0.0000    0.5000    0.0000  T T T
0.2500    0.7500    0.5000  T T T
0.0000    0.0000    0.0000  T T T
0.5000    0.0000    0.0000  T T T
0.2500    0.2500    0.5000  T T T
0.5000    0.5000    0.0000  T T T
0.7500    0.7500    0.5000  T T T
0.7500    0.2500    0.5000  T T T
0.2500    0.5000    0.2500  T T T
0.5000    0.7500    0.7500  T T T
0.2500    0.0000    0.2500  T T T
0.7500    0.0000    0.2500  T T T
0.5000    0.2500    0.7500  T T T
0.7500    0.5000    0.2500  T T T
1.0000    0.7500    0.7500  T T T
1.0000    0.2500    0.7500  T T T

POTCAR

In this example, we choose the ultrasoft-LDA PP of Si provided by VASP.

Run VASP

To compute the elastic constants we use the following script auto.B.serial to run vasp repeated with different values for the first C11 component of the first Bravais lattice vector, $a. Make sure you have previously calculated the equilibrium lattice constant for this supercell with 4 atoms, otherwise errors in the calculations could arise. The script creates automatically the POSCAR file as in the example above for every a especified in the script in the range from 0.997 to 1.003. For C44 script the value of "$a" is changed for the first and second component of the first Bravais lattice vector.

Script for C11

#!/bin/bash
rm WAVECAR 
for a in 0.997 0.998 0.999 1.000 1.001 1.002 1.003
do
cat > POSCAR << FIN
POSCAR for DC Si
5.402856        accurate equilibrium lattice constant
$a  0.0 0.0     first  Bravais lattice vector
0.0 1.0 0.0     second Bravais lattice vector
0.0 0.0 1.0     third  Bravais lattice vector
8               number of atoms per species
selective dynamics
direct           direct or cart (only first letter is significant)
0.0000  0.0000  0.0000 T T T
0.5000  0.5000  0.0000 T T T
0.0000  0.5000  0.5000 T T T
0.5000  0.0000  0.5000 T T T  
0.2500  0.2500  0.2500 T T T
0.7500  0.7500  0.2500 T T T
0.2500  0.7500  0.7500 T T T
0.7500  0.2500  0.7500 T T T
FIN
echo "a=$a"
./vasp

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

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

done

Script for C44

#!/bin/bash
rm WAVECAR 
for a in 0.997 0.998 0.999 1.000 1.001 1.002 1.003
do
cat > POSCAR << FIN
POSCAR for DC Si
5.402856      accurate equilibrium lattice constant
$a  $a  0.0   first  lattice vector
-1.0 1.0 0.0  second lattice vector
0.0 0.0 1.0   third  lattice vector
16            number of atoms per species
selective dynamics
direct           direct or cart (only first letter is significant)
0.0000    0.5000    0.0000  T T T
0.2500    0.7500    0.5000  T T T
0.0000    0.0000    0.0000  T T T
0.5000    0.0000    0.0000  T T T
0.2500    0.2500    0.5000  T T T
0.5000    0.5000    0.0000  T T T
0.7500    0.7500    0.5000  T T T
0.7500    0.2500    0.5000  T T T    
0.2500    0.5000    0.2500  T T T
0.5000    0.7500    0.7500  T T T
0.2500    0.0000    0.2500  T T T
0.7500    0.0000    0.2500  T T T
0.5000    0.2500    0.7500  T T T
0.7500    0.5000    0.2500  T T T
1.0000    0.7500    0.7500  T T T
1.0000    0.2500    0.7500   T T T
FIN
echo "a=$a"
./vasp

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

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

done

This script calculates C44 employing the shear method

Script for C44 shear method

#!/bin/bash
rm WAVECAR 
for a in -0.003 -0.002 -0.001 0.00 0.001 0.002 0.003
do
cat > POSCAR << FIN
POSCAR for FCC Si
5.402856      universal scaling factor
1.0  0.0 0.0  first  Bravais lattice vector
$a 1.0 0.0    second Bravais lattice vector
0.0 0.0 1.0   third  Bravais lattice vector
8             number of atoms per species
selective dynamics
direct           direct or cart (only first letter is significant)
0.00000  0.00000  0.00000 T T T
0.50000  0.50000  0.00000 T T T
0.00000  0.50000  0.50000 T T T
0.50000  0.00000  0.50000 T T T  
0.25000  0.25000  0.25000 T T T
0.75000  0.75000  0.25000 T T T
0.25000  0.75000  0.75000 T T T
0.75000  0.25000  0.75000 T T T
FIN

Analyze data

After running the two scripts above for C11 and C44, ./auto.B.serial it will create the data files: Elatt.C11.dat,platt.C11.dat, Elatt.C44.dat, and platt.C44.dat.

Launch octave and run the following functions fit_a0EB.m, fit_C11EB.m, fit_C44EB.m,

Run cal_VASP.m in Matlab,

fit_a0EB('Elatt.B.dat',4,4);
fit_C11EB('Elatt.C11.dat',a0,B);
fit_C44EB('Elatt.C44.dat',a0,C11,C12);

In MATLAB click "change folder" when you open the file and remove the semicolon to run it from the command line. The first line fits the energy data to a quadratic curve and computes the equilibrium lattice constant, cohesive energy and bulk modulus([a0, Ecoh, B]) respectively. Refer to previous calculation?!) The second lines computes the C11 and C12 elastic constants and the third line calculates C44. You should obtain C11 = 161 GPa, C12 = 179.2 GPa and C44 = 42 GPa.

Parallel computation

To run vasp in parallel, you need to submit vasp.pbs as

qsub vasp.pbs

You will need the following two files to do this calculation in parallel on SU-AHPCRC: auto.B.par and B.pbs.

@@ Line 43: / Line 43: @@
 '''POSCAR '''
-Here's an example for Si C11 and C12 calculation, later will be changed in PBS script every iteration.
+Here's an example for Si C11 and C12 calculation, which will be changed later in PBS script every iteration.
 <pre>
 POSCAR for DC Si
@@ Line 63: / Line 63: @@
 </pre>
-Here's an example for Si C44 calculation, later will be changed in PBS script every iteration.
+Here's an example for Si C44 calculation, which will be changed later in PBS script every iteration.
 <pre>

VASP Computing Elastic Constants of Si: Difference between revisions

Revision as of 05:56, 2 September 2015

Contents

Input files

Run VASP

Analyze data

Parallel computation

Navigation menu

VASP Computing Elastic Constants of Si: Difference between revisions

Revision as of 05:56, 2 September 2015

Input files

Run VASP

Analyze data

Parallel computation

Navigation menu

Search