220 lines
3.9 KiB
Bash
Executable File
220 lines
3.9 KiB
Bash
Executable File
# A script to calculate RESP charge based on ORCA and Multiwfn
|
|
# Written by Tian Lu (sobereva@sina.com), last update: 2022-Mar-8
|
|
# Examples:
|
|
# Calculating neutral singlet molecule in water: RESP_ORCA.sh H2O.xyz
|
|
# Calculating anionic singlet molecule in water: RESP_ORCA.sh ani.pdb -1 1
|
|
# Calculating neutral triplet molecule in ethanol: RESP_ORCA.sh nico.xyz 0 3 ethanol
|
|
# Calculating neutral singlet molecule in ETHYL ETHANOATE: RESP_ORCA.sh nico.xyz 0 1 "ETHYL ETHANOATE"
|
|
# Calculating cation doublet molecule in vacuum: RESP_ORCA.sh maki.xyz 0 2 gas
|
|
|
|
#!/bin/bash
|
|
|
|
#Set actual paths of ORCA and orca_2mkl utility here
|
|
ORCA="/root/orca_6_0_0_shared_openmpi416_avx2/orca"
|
|
orca_2mkl="/root/orca_6_0_0_shared_openmpi416_avx2/orca_2mkl"
|
|
|
|
#Set number of CPU cores used in calculation here
|
|
nprocs=12
|
|
maxcore=1000
|
|
|
|
keyword_opt="! B97-3c opt"
|
|
keyword_SP="! B3LYP/G D3 def2-TZVP def2/J RIJCOSX"
|
|
|
|
export inname=$1
|
|
filename=${inname%.*}
|
|
suffix=${inname##*.}
|
|
|
|
### Parse arguments
|
|
if [ $2 ];then
|
|
echo "Net charge = $2"
|
|
chg=$2
|
|
else
|
|
echo "Net charge was not defined. Default to 0"
|
|
chg=0
|
|
fi
|
|
|
|
if [ $3 ];then
|
|
echo "Spin multiplicity = $3"
|
|
multi=$3
|
|
else
|
|
echo "Spin multiplicity was not defined. Default to 1"
|
|
multi=1
|
|
fi
|
|
|
|
if [ "$4" ];then
|
|
if [[ $4 == gas ]] ; then
|
|
echo Calculations will be done in vacuum
|
|
solvent=""
|
|
else
|
|
echo Solvent is $4
|
|
solvent=$4
|
|
fi
|
|
else
|
|
echo "Solvent name was not defined. Default to water"
|
|
solvent="Water"
|
|
fi
|
|
|
|
### Convert current file to tmp.xyz
|
|
Multiwfn_noGUI $1 > /dev/null << EOF
|
|
100
|
|
2
|
|
2
|
|
tmp.xyz
|
|
0
|
|
q
|
|
EOF
|
|
|
|
### Create input file for optimization (opt.inp)
|
|
cat << EOF > opt.inp
|
|
$keyword_opt
|
|
%maxcore $maxcore
|
|
%pal nprocs $nprocs end
|
|
%geom Convergence loose end
|
|
EOF
|
|
if [[ $solvent != "" ]] ; then
|
|
cat << EOF >> opt.inp
|
|
%cpcm
|
|
smd true
|
|
SMDsolvent "$solvent"
|
|
end
|
|
EOF
|
|
fi
|
|
echo "* xyz $chg $multi" >> opt.inp
|
|
awk '{if (NR>2) print }' tmp.xyz >> opt.inp
|
|
echo "*" >> opt.inp
|
|
rm tmp.xyz
|
|
|
|
### Run optimization
|
|
echo Running optimization task via ORCA...
|
|
$ORCA opt.inp > opt.out
|
|
|
|
if grep -Fq "ORCA TERMINATED NORMALLY" opt.out
|
|
then
|
|
echo Done!
|
|
mv opt.xyz tmp.xyz
|
|
rm opt.* opt_*
|
|
else
|
|
echo The optimization task has failed! Please check content of opt.out to find reason
|
|
echo The script is terminated
|
|
mv opt.out tmp.out
|
|
rm opt.* opt_*
|
|
mv tmp.out opt.out
|
|
exit 1
|
|
fi
|
|
|
|
### Create input file for single point (SP.inp)
|
|
cat << EOF > SP.inp
|
|
$keyword_SP
|
|
%maxcore $maxcore
|
|
%pal nprocs $nprocs end
|
|
EOF
|
|
if [[ $solvent != "" ]] ; then
|
|
cat << EOF >> SP.inp
|
|
%cpcm
|
|
smd true
|
|
SMDsolvent "$solvent"
|
|
end
|
|
EOF
|
|
fi
|
|
echo "* xyz $chg $multi" >> SP.inp
|
|
awk '{if (NR>2) print }' tmp.xyz >> SP.inp
|
|
echo "*" >> SP.inp
|
|
rm tmp.xyz
|
|
|
|
### Run single point
|
|
echo Running single point task via ORCA...
|
|
$ORCA SP.inp > SP.out
|
|
|
|
if grep -Fq "ORCA TERMINATED NORMALLY" SP.out
|
|
then
|
|
echo Done!
|
|
else
|
|
echo The single point task has failed! Please check content of SP.out to find reason
|
|
echo The script is terminated
|
|
mv SP.out tmp.out
|
|
rm SP.* SP_*
|
|
mv tmp.out SP.out
|
|
exit 1
|
|
fi
|
|
|
|
### Convert to .molden file
|
|
echo Running orca_2mkl...
|
|
$orca_2mkl SP -molden > /dev/null
|
|
|
|
#Field containing number of valence electrons of def2 pseudopotential basis set
|
|
cat << EOF > Nval.txt
|
|
[Nval]
|
|
Rb 9
|
|
Sr 10
|
|
Y 11
|
|
Zr 12
|
|
Nb 13
|
|
Mo 14
|
|
Tc 15
|
|
Ru 16
|
|
Rh 17
|
|
Pd 18
|
|
Ag 19
|
|
Cd 20
|
|
In 21
|
|
Sn 22
|
|
Sb 23
|
|
Te 24
|
|
I 25
|
|
Xe 26
|
|
Cs 9
|
|
Ba 10
|
|
La 11
|
|
Ce 30
|
|
Pr 31
|
|
Nd 32
|
|
Pm 33
|
|
Sm 34
|
|
Eu 35
|
|
Gd 36
|
|
Tb 37
|
|
Dy 38
|
|
Ho 39
|
|
Er 40
|
|
Tm 41
|
|
Yb 42
|
|
Lu 43
|
|
Hf 12
|
|
Ta 13
|
|
W 14
|
|
Re 15
|
|
Os 16
|
|
Ir 17
|
|
Pt 18
|
|
Au 19
|
|
Hg 20
|
|
Tl 21
|
|
Pb 22
|
|
Bi 23
|
|
Po 24
|
|
At 25
|
|
Rn 26
|
|
EOF
|
|
|
|
cat Nval.txt SP.molden.input > SP.molden
|
|
|
|
### Calculate RESP charges
|
|
echo Running Multiwfn_noGUI...
|
|
Multiwfn_noGUI SP.molden -ispecial 1 > /dev/null << EOF
|
|
7
|
|
18
|
|
1
|
|
y
|
|
0
|
|
0
|
|
q
|
|
EOF
|
|
|
|
### Finalize
|
|
chgname=${1//$suffix/chg}
|
|
mv SP.chg $chgname
|
|
rm SP.* SP_* Nval.txt
|
|
|
|
echo Finished! The optimized atomic coordinates with RESP charges \(the last column\) have been exported to $chgname in current folder
|
|
echo Please properly cite Multiwfn_noGUI in your publication according to \"How to cite Multiwfn.pdf\" in Multiwfn package
|