Get-irc-path.sh

From Wiki
Jump to navigation Jump to search

go back to Main Page, Computational Resources, Chemistry & More, Computational Codes, VASP, IDM or Scripts for VASP

DESCRIPTION[edit]

Use the script "get-irc-path.sh" (vide infra) to process an output of VASP IRC calculation (more excatly: the OSZICARs in subdirectories dir1 and dir2) in order to obtain a full energetic trajectory of the IRC reaction pathway:

get-irc-path.sh dir1 dir2

or without arguments simply

get-irc-path.sh

if the subdirectories with forward and reverse IRC runs are called "1" and "2".

The output is written into "trajectory.dat" file. The 1st column is an arbitrary unit (# of the ionic step), the 2nd column lists the corresponding energy.

CODE[edit]

#!/bin/bash

# Script for extracting the IRC (Intrinsic Reaction Coordinate) trajectory
# from 2 subdirectories and 
# (C) David Karhanek 2010-09-29, ICIQ Tarragona

# Unless the dir names are the defaults "1" and "2", submit them to the script as 2 arguments
if [ $# != 0 ] ; then if [ $# == 2 ] ; then dir1=$1 ; dir2=$2 ; fi ; else dir1=1 ; dir2=2 ; fi

# Parse energies and their last values
grep E0 ${dir1}/OSZICAR | awk '{printf "%14.7f\n",$5}' > 1.dat
grep E0 ${dir2}/OSZICAR | awk '{printf "%14.7f\n",$5}' > 2.dat
E01=`sed '$!d' 1.dat` ; E02=`sed '$!d' 2.dat`

# Decide which directory relates to products, put the energies together
if [ $(echo "$E02 < $E01" | bc -q) == 1 ] ; then ini=1 ; fin=2 ; else ini=2 ; fin=1 ; fi
tac $ini.dat | sed '$d' > 3.dat ; cat $fin.dat >> 3.dat
awk '{printf "%03u%14.7f\n",NR,$1}' 3.dat > trajectory.dat

# Perform simple statistics on the initial (IS), transition (TS) and final (FS) states
# (the TS energy check : the values should not differ by more than 0.01 eV)
ETSINI=`sed q $ini.dat` ; ETSFIN=`sed q $fin.dat` ; EIS=`sed '$!d' $ini.dat` ; EFS=`sed '$!d' $fin.dat`
if [ $(echo "a=($ETSFIN-($ETSINI))/1;if(0>a)a*=-1;a<0.01" | bc -ql) == 1 ] ; then LTSM="OK" ; else LTSM="NOT_OK_!!!" ; fi
ETS=$(echo "($ETSFIN+$ETSINI)/2" | bc -ql)
EA=$(echo "$ETS-$EIS" | bc -ql) ; ER=$(echo "$EFS-$EIS" | bc -ql)
of=statistics.txt ; if [ -f $of ] ; then mv $of old-$of ; fi
printf "TS energies in both runs matching? .. %s\n" $LTSM >> $of
printf "E(IS) = %14.7f eV  Initial State    (subdir %s)\n" $EIS $ini >> $of
printf "E(TS) = %14.7f eV  Transition State\n" $(echo "($ETSFIN+$ETSINI)/2" | bc -ql) >> $of
printf "E(IS) = %14.7f eV  Final State      (subdir %s)\n" $EFS $fin >> $of
printf "E(A)  = %14.7f eV  Energy Barrier\n" $EA >> $of
printf "E(R)  = %14.7f eV  Reaction Heat\n" $ER >> $of

# List output on-screen, clean-up
printf "\n== IRC traj. ==\n" ; cat trajectory.dat
printf "\n== Statistics: ==\n" ; cat $of
rm -rf ?.dat > /dev/null

COMMON ERRORS[edit]

If you find some bugs, report them to David.

--Dkarhanek 20:29, 29 September 2010 (CEST)