Bryce's VASP Quickstart Guide

Example calculation: Rocksalt FeO

The first step is to define the geometry of your system. Usually you will tend to think of your compound in terms of its conventional unit cell, but when running VASP, you'll want to minimize computational time by utilizing a potentially much less intuitive primitive unit cell that contains fewer atoms but still fully describes the structure. For an illustration of the conventional-primitive distinction for bcc and fcc, see:

http://books.google.com/books?id=5zVozqxOqLwC&pg=PA28&dq=bcc+primitive+cell&sig=eOJ0MfXfzLVrraieJA3JoiLHdbY

Three different codes are used in this tutorial. GULP is a molecular simulation code available at http://projects.ivec.org/gulp/. You will also need convasp, a software written by Dane Morgan and Stefono Curtarolo and a new version can be downloaded from Curtarolo Group's website. The other one is ezvasp, which is a part of the ATAT package by Axel van de Walle.

GULP

The GULP utility will help you convert from conventional to primitive. GULP requires an input file of the format (name it feo.in):

$ cat feo.in

FeO # system name

cell # unit cell parameters go on the next line

4.3 4.3 4.3 90 90 90 # lattice parameters a, b, c, alpha, beta, gamma

frac # to use fractional coordinates for atoms

O 0.5 0.5 0.5 # oxygen position for rocksalt

Fe 0 0 0 # iron position for rocksalt

space # space group goes on the next line

225 # specify the space group (see note)

output xyz feo # output style

The one item in this file that deserves additional explanation is the space group designation. FeO's rocksalt structure is described by what is ordinarily referred to as the Fm3m space group. However, GULP needs a numerical designation for this structure: 225, in this case. There are 230 space groups in all; their numbering scheme and corresponding traditional labels (i.e., Fm3m) can be viewed at:

http://en.wikipedia.org/wiki/List_of_the_230_crystallographic_3D_space_groups

Note: IF you have a non-standard space group, GULP can also do this for you - under space group enter the name of the traditional label, with spaces in the middle, but no spaces for the part (1,2,3,4,6)/(a,b,c)...etc. For example, instead of 225, you can type "F M 3 M". An example of the special case, I2/a, which is a non-standard form of C 2/c (No. 15), I can enter the full space group name like so: "I 1 2/a 1"

You can also transform the a nonstandard crystal structure to the standard crystal structure using this utility: http://www.cryst.ehu.es/cryst/setstru.html

We're now ready to have GULP find the primitive cell:

$ gulp < feo.in > feo.out

(If you're interested, you can explicitly view the differences between the primitive and conventional cells in the contents of feo.out; for further discussion of how to use GULP, see Stefano Curtarolo's thesis: http://burgaz.mit.edu/PUBLICATIONS/THESIS/stefano.pdf).

We have two output files of interest: feo.out and feo.xyz. feo.xyz will form the basis for our VASP input geometry file, or POSCAR. Right now, it looks like this:

2

SCF Done 0.00000000

O 2.150000000 2.150000000 2.150000000

Fe 0.000000000 0.000000000 0.000000000

These are the Cartesian (not fractional!) atom positions. The 2 above indicates the total number of atoms. We also need to be able to tell VASP the lattice vectors it should use. GULP stored that information in the feo.out file. Specifically, search for the part of the file that looks like this:

Cartesian lattice vectors (Angstroms) :

0.000000 2.150000 2.150000

2.150000 0.000000 2.150000

2.150000 2.150000 0.000000

These are your lattice vectors. Paste them (just the three lines of numbers, not the first line) into the top of your feo.xyz file. Then rearrange the .xyz file such that it has this format (note that your file cannot actually have the comments; they have been added for clarity):

FeO #1

1 #2

0.000000 2.150000 2.150000

2.150000 0.000000 2.150000

2.150000 2.150000 0.000000 #3

1 1 #4

Cartesian #5

2.150000000 2.150000000 2.150000000 O

0.000000000 0.000000000 0.000000000 Fe #6

#1: label

#2: scaling factor; leave it as 1 unless you know what you're doing

#3: three lattice vectors

#4: number of O atoms, number of Fe atoms

#5: specify coordinate style; we're using Cartesian

#6: atom positions

An important thing to note here is that the names of the atoms now come after their Cartesian positions.

One line of code that can do this is:

$ sed -i 's/\([A-Za-z]*\) \([0-9. -]*\)/\2 \1/g' NAME.xyz

This moves all the atom positions in the file NAME.xyz to the end of the line and you should verify that the format is correct. Pay attention to the line #4 where the number of atoms is specified. It might be messed up !

Anyway, the above is the format required for another utility, convasp.

convasp

convasp will convert your atom positions back into fractional (yeah, this is a bit circular) and produce a file that looks just like a VASP POSCAR:

$ convasp -direct < feo.xyz > feo.vaspxyz

This command produces the feo.vaspxyz file, which looks like this:

$ cat feo.vaspxyz

FeO

1.0000000000

0.0000000000 2.1500000000 2.1500000000

2.1500000000 0.0000000000 2.1500000000

2.1500000000 2.1500000000 0.0000000000

1 1

D

0.5000000000 0.5000000000 0.5000000000 O

0.0000000000 0.0000000000 0.0000000000 Fe

Give yourself a pat on the back! This is all the geometrical information VASP needs. But what should VASP do with these atoms? You need to give it some more help. Copy your feo.vaspxyz file into a new file, called vasp.in. Add some additional lines, and delete the “1 1” line near the end, so that vasp.in looks like this:

[INCAR]

SYSTEM = FeO

ISTART = 0

ISMEAR = 1

SIGMA = 0.1

ISIF = 3

ISPIN = 2

PREC = HIGH

IBRION = 2

NSW = 30

ENCUT = 400

NPAR = 2

KSCHEME= Monkhorst-Pack

KPPRA = 10000

DOGGA

[POSCAR]

FeO

1.0000000000

0.0000000000 2.1500000000 2.1500000000

2.1500000000 0.0000000000 2.1500000000

2.1500000000 2.1500000000 0.0000000000

D

0.5000000000 0.5000000000 0.5000000000 O

0.0000000000 0.0000000000 0.0000000000 Fe

The VASP INCAR file contains configuration items that tell VASP how it should perform calculations on your atoms. For a good description of the parameters, see http://hidra.iqfr.csic.es/vasp/node92.html. Right now, don't worry about changing them. Finally, we're ready to generate actual INCAR and POSCAR files for VASP to read. Another tool, called ezvasp, can do this based on the contents of your vasp.in:

$ ezvasp -n vasp.in

The -n flag tells ezvasp not to launch VASP immediately.

There are some elements for which normal potentials do not exist. For VASP, it is required that you replace the name of the element in the vasp.in file with the accurate potential for that element. These elements are listed below, along with their replacement:

Ba/Ba_sv, Cs/Cs_sv, Dy/Dy_3, Er/Er_3, Eu/Eu_3, Gd/Gd_3, Ho/Ho_3, Lu/Lu_3, Nd/Nd_3,

Nb/Nb_sv, K/K_sv, Pr/Pr_3, Pm/Pm_3, Rb/Rb_sv, Sm/Sm_3, Sr/Sr_sv, Tb/Tb_3, Y/Y_sv

Also, occasionally, you may want to use 'hard' potentials for more accurate results, instead of the 'normal' potentials. Generally this can be done by adding a '_h' after the element name, for example, you can use a hard potential for hydrogen using H_h. See the potpaw_GGA directory for a full list of the potentials in VASP, at /usr/local/vasp/potpaw_GGA/. Note that radon and radium do not have GGA potentials in VASP.

After all of this work, we ought to check to make sure our generated POSCAR file still represents the system correctly. A useful way to do this is to verify its space group symmetry using convasp, together with a tool called platon and a script called platonSG:

$ convasp -platon < POSCAR | platonSG

Fm-3m #225

Fm-3m #225

Hooray! Looks like everything is just as it should be.

PBS job submission

Our last step before launching VASP is to prepare a small file that will submit our job to the minicluster's scheduling queue. Create a file called FeO.q and structure it like this:

#PBS -l ncpus=3

#PBS -l walltime=1:00:00

cd $PBS_O_WORKDIR

mpirun -np 3 /usr/local/bin/vasp

The only two items you'll need to change here are ncpus (along with -np below; the two must agree) and walltime. ncpus represents the number of processors you're requesting; on the minicluster, you can ask for up to 7. Walltime is how long the cluster should work on the calculation. Keep in mind that as you request more resources, your priority in the queue will drop! At last, it's time to run! Submit your job to the minicluster like this:

$ qsub FeO_test.q

747.encina.northwestern.edu

747 is your job's number. To watch VASP's progress,

$ qcat 747

encina

Working in /home/bmeredig/feo-dongwon/benching

stdout

running on 3 nodes

distr: one band on 1 nodes, 16 groups

vasp.4.6.34 5Dec07 complex

POSCAR found : 2 types and 2 ions

LDA part: xc-table for Ceperly-Alder, standard interpolation

POSCAR, INCAR and KPOINTS ok, starting setup

WARNING: wrap around errors must be expected

FFT: planning ... 1

reading WAVECAR

entering main loop

N E dE d eps ncg rms rms(c)

DAV: 1 0.599819007010E+02 0.59982E+02 -0.65888E+03 9786 0.115E+03

DAV: 2 -0.117810265339E+02 -0.71763E+02 -0.65751E+02 9564 0.228E+02

DAV: 3 -0.159528643008E+02 -0.41718E+01 -0.40187E+01 10671 0.639E+01

The bottom of the output represents energy relaxation steps—VASP is working!

Effects of Kpoints

Luckily, ezvasp will automatically determine your kpoints file for a cell given some value for KPPRA, or Kpoints per reciprocal atom. There is a value of KPPRA for which your energies will be converged to within 1 meV. For most metals, this value is around 5000-6000 KPPRA, although it is within the Wolverton mentality to be conservative and safe, and so I personally tend to use around 7500 KPPRA. For some materials, such as oxides, convergence can occur at much lower KPPRA. Run tests to observe what these are.

Note that doubling your Kpoint mesh can increase your computation time by an order of magnitude. Please be considerate of the group resources when determining ideal Kpoints.

Encut

Generally, one can determine Encut by multiplying the greatest ENMAX of the atoms in your POTCAR by 1.25. You can find this ENMAX value by using the line "grep ENMAX POTCAR" in the directory of your POTCAR. Remember that if you are modeling a system, you must use the exact same INCAR, and therefore the same ENCUT, for all your calculations in the system. E.g. if you have one element whose ENCUT dominates the system, such as hard potential Hydrogen (875 eV), but some part of your calculation doesn't include it, your ENCUT for that phase/structure still must be done at Encut = 875.

Automation

If you are doing this process many times (or even just a few!), Wenhao and Danny have written a script that automates this process from beginning to end, using a .cif file as input, and can be found here: *cif to VASP. If you are interested, please contact Wenhao. There isn't exactly a final, one-case-fits-all version - it's constantly being edited and thus making it hard to post on the wiki.

Also, if you are doing a ground state crystal structure database search, the script can also replace atoms for you, for example, if you are searching for the ground state of Ca(AlH4)2, you can replace any structure with AB2X8 format with the calcium atoms in the A position, Aluminum atoms in the B position, and H atoms in the X position. Again, contact Wenhao for details.

Other things to be written into this guide

mesh style, magnetism, ediff, ediffg, LDA vs GGA