VASP How-tos

(updated )

Comments

Tutorials

  • Work through tutorials SharePoint folder from a 2003 VASP workshop.
  • Tutorials in SharePoint folder for first-principle simulation codes.
  • Papers in the group Zotero library, inside folders labeled classics.
  • Classic texts listed under Computational modeling and Catalysis.
  • Official YouTube channel of VASP.

Set up vasp-tools on the group server

  • SSH (PuTTY for Windows) to the group server and complete the steps below.
cd ~
ln -sf /share/apps/vasp-tools .
cp -a /home/common/sample-input/vasp/my-vasp-tools .
module add vasp-tools/linux

Edit ~/my-vasp-tools/user_specs to suit your system and use your Email address.

  • To find out the available options setrun supports, locate VASPROOT from the output of the following command:
module show vasp-tools

Typically, the result will be $HOME/vasp-tools, you can then examine the set of system and job types:

find $HOME/vasp-tools/runfiles/system/ -mindepth 2 -maxdepth 2 -type d  # for system types, or
find $HOME/vasp-tools/runfiles/jobtype/ -mindepth 2 -maxdepth 2 -type d  # for job types

You will see something like $HOME/vasp-tools/runfiles/system/gas/molecule and $HOME/vasp-tools/runfiles/jobtype/nrg/low, which means you can use -s gas/molecule and -j nrg/low for isolated molecular calculation and low-precision single-point energy calculation, respectively.

To create new systems or job types, or override existing ones, simply create the necessary folder hierarchy in ~/my-vasp-tools. For example, if we want to use a different INCAR template for gas/molecule, copy the existing files from $HOME/vasp-tools to use as the basis:

mkdir -p ~/my-vasp-tools/runfiles/system/gas
cp -a ~/vasp-tools/runfiles/system/gas/molecule ~/my-vasp-tools/runfiles/system/gas/

We then make the necessary changes to INCAR files in all sub directories. We can also add KPOINTS files to these sub directories and skip the -K option when running setrun.

Similarly, create/copy and edit files in ~/my-vasp-tools/runfiles/jobtype to customize new job types to use with the -j option of setrun.

  • Make sure to understand Step 4 under Shotgun basic workflow. Remember to copy the VASP executable to the target HPC machine. Where the VASP executable should be located is determined by the submission script. For example, /home/common/bin/sub-vasp-bridges.script has the following line PROG="$HOME/vasp-tools/hosts/mesabi/vaspexec/vasp5.4.1_ifort_impi_mkl_vtst_vaspsol", which means a binary executable named vasp5.4.1_ifort_impi_mkl_vtst_vaspsol should exist in your home directory on Bridges under vasp-tools/hosts/mesabi/vaspexec. We can rsync the executable as follows:
rsync -avzP /share/apps/vasp-tools/hosts/mesabi bridges:.
ssh bridges 'mkdir -p ~/vasp-tools/hosts; mv mesabi ~/vasp-tools/hosts/'
  • VASP executable: If VASP has not been compiled on a new HPC system, we can use a build for a different HPC machine assuming they have the same microarchitecture (see this list). Most VASP builds were compiled with Intel compilers + MKL. Remember to module add the latest Intel compiler suite.

Set up Open Studio

  • Log in with your UMass email at Open Studio .
  • Download the token script from the drop-down menu at the top right corner and copy it to your ~/my-vasp-tools/bin folder on the group server.
  • Give the script executable permission and create the ~/website folder:
chmod +x ~/my-vasp-tools/bin/post2web.sh
mkdir -p ~/website
  • After the setup, results_opt, results_neb, and results_dim scripts will upload VASP results to your Google Drive and return a link for web visualization.
  • You can customize the look and feel, make measurements, or even use Jmol console commands (right click on the visualization pane). You can save the customization by clicking the Export Online button.

Single-point energy calculations / geometry optimizations

This video gives a good introduction on geometry optimization, covering the theory, and discussing the different algorithms implemented in VASP. Below is the sequence of steps for running a VASP calculation to compute a single-point energy or to optimize a structure.

  • On your workstation
  1. Use Schrodinger or any other suitable software to prepare the structure to be optimized, export it as a .pdb file.
  2. Use rsync / scp / WinSCP (Windows only) to upload the file to group server.
  • On the group server
  1. Convert the .pdb file to VASP POSCAR format. It may be necessary to edit the .pdb file by deleting all the “CONECT” lines using a text editor. 
    cd /path/to/projects/ABC/moleculeX/reactant  # go to the directory where the .pdb file is uploaded.
module add python/conda3.6  # in order to use the ASE python module
ase convert moleculeX-reactant.pdb POSCAR  # converts the .pdb file to POSCAR
  2. Set up other input files. 
    module add vasp-tools  # in order to use the setrun script (and others)
setrun -f PBE -F D3BJ -p PAW54 -K 1x1x1 -s gas/molecule -j opt/high

    setrun generates 3 directories: STEP1, STEP2, STEP3, for a 3-stage optimization. STEP3 uses high K-point but low precision. It’s mostly for metals and not relevant for isolated molecules or zeolite systems. Remove it by

    rm -r STEP3

    In the above command, -s gas/molecule specifies the type of system and -j opt/high specifies the job type.

  3. Prepare submission scripts. The best practice is to write a generic submission script and put it in the local asset directory (e.g., /path/to/projects/ABC/bin/). Then create symbolic links: 
    ln -sf ../../bin/sub-vasp.script .  # Run in current directory /path/to/projects/ABC/moleculeX/reactant

    NOTE: these symlinks must use relative paths instead of absolute paths, because the local asset directory will be synced to the remote target under your home directory, which may have a different absolute path.

  4. Submit the job using shotgun 
    shotgun submit -c 8 -m 5200 -n 1 -t 24 -r ghpcc06 -q long -j ../../bin/sub-vasp-ghpcc.script  # submitting the job

    NOTE: the ./ part in ./sub-vasp-ghpcc.script is required. See shotgun --help submit for a list of options that can be used.

NEB calculations

The nudged elastic band method uses a series of structures to approximate the states along a reaction pathway and optimize them in a way so that they approach the minimum energy pathway. Below is the sequence of steps for setting up a NEB calculation.

  1. On your workstation, use Schrodinger to gradually modify your reactant structure and morph it into the desired product structure. Save 4 — 8 intermediate states, depending on how complicated the reaction is. Since we will use linear interpolation of Cartesian coordinates to construct the initial guess for the reaction pathway, translation is therefore easy while rotation is difficult.

    NOTE: make sure you understand this point. To help with understanding this, construct two “reactions”: pure translation and pure rotation of a single methane molecule. Use only reactant and product to construct a 8-image pathway for both cases (see below for the steps). Visualize the movie of the two pathways created and pay attention to how bonds are compressed in the rotation “reaction”. We would want to avoid feeding obviously unphysical structures to DFT calculations.

  2. Optimize the reactant and product structures, following the steps under Geometry Optimizations.
  3. Export the intermediate states as .pdb files and upload them to group server.
  4. SSH to group server. On the group server, convert all the .pdb files to POSCARs: 
    cd /path/to/projects/ABC/moleculeX/NEB  # go to the directory where the intermediate images are uploaded.
module add python/conda3.6  # in order to use the ASE python module
ase convert moleculeX-intermediate1.pdb POSCAR1
ase convert moleculeX-intermediate2.pdb POSCAR2
# ... say you have saved five intermediate images ...
ase convert moleculeX-intermediate5.pdb POSCAR5

    Copy the optimized reactant and product structures as POSCAR0 and POSCAR6. Make sure to use the FINAL optimized structure from STEP2, not the stage 1 result in STEP1.

    cp -a /path/to/projects/ABC/moleculeX/reactant/STEP2/CONTCAR POSCAR0
cp -a /path/to/projects/ABC/moleculeX/product/STEP2/CONTCAR POSCAR6
  5. Set up the input files: 
    module add vasp-tools  # in order to use the setrun script (and others)
setrun -f PBE -F D3BJ -p PAW54 -K 1x1x1 -s zeolite/surface -j neb/high -i 8

    -i 8 specifies that this is a 8-image setup. For challenging reactions, you can increase the number of images to, say, 16, using -i 16

  6. Visualize the output movie file POSCAR.res that represents the initial guess of the reaction pathway, and make sure that it matches your expectation.
  7. Edit INCAR to replace IMAGES = $IMAGES with the actual number of images used; e.g., IMAGES = 8. You can also do it from the command line: 
    sed -i 's/$IMAGES/8/' INCAR
  8. Link the submission scripts: 
    ln -sf ../../bin/sub-vasp.script .  # Run in current directory /path/to/projects/ABC/moleculeX/NEB
  9. Submit the job using shotgun 
    shotgun submit -n 2 -t 24 -r unity-40Cx2 -j ../../bin/sub-vasp-ghpcc.script  # submitting the job

    Here we are requesting 2 full nodes, all 40 cores on each node, for a total of 80 (physical) cores. See the hyper-threading section here for physical vs. logical cores. With 80 cores, each image will be run on 80/8 = 10 cores. Always aim for an even number of cores per image by 1) targeting an HPC with nodes that have physical cores that can divide your images evenly or 2) adjusting the number of cores per image to fit the HPC’s architecture. Use 8-12 cores per image, usually no more than 12 cores. For example, a 16-image job could run on Anvil which has nodes with 128 cores, so the job will use 8 cores per image. Since NEB calculations typically use many cores, there is no good reason to request -c N with N less than the full core count per node.

  10. After the calculation converges, copy OUTCARs for reactant and product to the corresponding directories. Then, perform analysis and move the NEB trajectory movie from the ~/website folder. Make sure the copied movie file corresponds to the current calculation since the files may be in a sub-directory such as ~/website/X where X is an integer. One way to verify is to compare file time-stamps using ls -alhF. Another way is to compare file contents using vim -d ~/website/1/neb2.dat /path/to/calc/neb2.dat 
    cp -a /path/to/projects/ABC/moleculeX/reactant/STEP2/OUTCAR 00/
cp -a /path/to/projects/ABC/moleculeX/product/STEP2/OUTCAR 09/
module add vasp-tools  # in order to use the results_neb script (and others)
results_neb
mv ~/website/CONTCAR.res .

Dimer calculations

Dimer calculations use three images to represent the transition state: a center image for the transition-state structure, and two others for the direction of the reaction pathways. The following is the sequence of steps for setting up a dimer calculation.

  1. The initial structure for dimer optimizations needs to be reasonably close to a saddle point. It is usually obtained from a NEB calculation. Assuming NEB image 05 is closest to the transition state of interest, we copy it as the dimer center image. 
    cp -a /path/to/NEB/05/CONTCAR POSCAR
  2. Next we need to provide a MODECAR file to represent the reaction direction. Using chemical intuition, pick either images 04 -> 05 or 05 -> 06 from the NEB calculation, and generate MODECAR. 
    cp -a /path/to/NEB/04/CONTCAR POSCAR02  # or cp -a /path/to/NEB/06/CONTCAR POSCAR02
module add vasp-tools  # in order to use the modemake, dimmovie, and setrun scripts
modemake POSCAR02 POSCAR  # or modemake POSCAR POSCAR02, the order does not really matter
dimmovie

    dimmovie generates a movie file, which you can visualize to verify the structure and mode.

  3. Set up the input files: 
    setrun -f PBE -F D3BJ -p PAW54 -K 1x1x1 -s zeolite/surface -j dim/high

    setrun generates 3 directories: STEP1, STEP2, STEP3, for a 3-stage optimization. STEP3 uses high K-point but low precision. It’s mostly for metals and not relevant for isolated molecules or zeolite systems. Remove it by

    rm -r STEP3
  4. Link the submission scripts: 
    ln -sf ../../bin/sub-vasp.script .  # Run in current directory /path/to/projects/ABC/moleculeX/dimer
  5. Submit the job using shotgun 
    shotgun submit -c 8 -m 5200 -n 1 -t 24 -r ghpcc06 -q long -j ../../bin/sub-vasp-ghpcc.script  # submitting the job
  6. After the calculation converges, perform the analysis in the STEP2 folder. The trajectory movie files CONTCAR.res and Dimer_Mode.res are saved in ~/website and can be moved to the job folder, especially before running another results_dim analysis since both files will be overwritten. 
    module add vasp-tools
cd /path/to/dim/STEP2
results_dim
mv ~/website/{CONTCAR.res,Dimer_Mode.res} .

Metadynamics

Metadynamics is used to estimate the free energy and other state functions of a system, where ergodicity is hindered by the form of the system’s energy landscape.

VASP Settings

  1. Standard MD-related tags.
  2. To continue the Metadynamics run. 
    cp HILLSPOT PENALTYPOT

Plumed Analysis

  1. Free energy calculated using bias potential (sum hills)
  2. Free energy calculated by reweighting

Center-of-mass drifting

  1. Our simulation box is an isolated system that experiences no external force and has no initial velocity, so its center of mass (COM) should not move. If it does, it is usually a symptom of bad numerics, i.e., energy and forces are not computed precisely enough.
  2. Say if the COM moves with a velocity of v and each atom moves relative to the COM with velocity v_i, then v_i should be the velocity used to compute the system kinetic energy and temperature. However, with almost all MD codes, (v + v_i) will be used instead. So if the system is thermostatted at T (using the inflated velocity), the real temperature will be lower than T.
  3. The ideal solution to this problem may be quite expensive. Oftentimes it can be permissible to remove the center of mass velocity.

Comments

Comments powered by Disqus