Tutorial2 script

for i in "pure_methanol" "pure_water" "mixture_1" "mixture_2" "mixture_3" "mixture_4" "mixture_5" "mixture_6" "mixture_7" "mixture_8" "mixture_9" "mixture_10"; do

if [ $i == "pure_methanol" ]; then
mols="MET    400"
tc_grps="MET"
energygrps="MET"
tau_t=0.1; ref_t=298.
elif [ $i == "pure_water" ]; then
mols="SOL    800"
tc_grps="SOL"
energygrps="SOL"
tau_t=0.1; ref_t=298.
elif [ $i == "mixture"* ]; then
mols="MET    200\nSOL    400"
tc_grps="MET SOL"
energygrps="MET SOL"
tau_t="0.1 0.1"; ref_t="298. 298"
fi

cat>topol.top<<eof
#include "oplsaa.ff/forcefield.itp"
#include "oplsaa.ff/spce.itp"
#include "methanol.itp"

[ system ]
Water+methanol

[ molecules ]
; name  number
eof
echo -e $mols >> topol.top

awk -v tc_grps=$tc_grps -v energygrps=$energygrps -v tau_t=$tau_t -v ref_t=$ref_t \
'{if($1 == "Tcoupl"){print "Tcoupl = v-rescale"}
  else if($1 == "Pcoupl"){print "Pcoupl = Berendsen"}
  else if($1 == "Pcoupl"){print"gen_vel = yes"}
  else if($1 == "nsteps"){print "nsteps = 250000"}
  else if($1 == "nstxout"){print "nstxout = 1000"}
  else if($1 == "nstvout"){print "nstvout = 1000"}
  else if($1 == "tc-grps"){print "tc-grps =", tc_grps}
  else if($1 == "energygrps"){print "energygrps =", energygrps}
  else if($1 == "tau_t"){print "tau_t =", tau_t}
  else if($1 == "ref_t"){print "ref_t =", ref_t}
  else {print $0}}' md.mdp > temp; cp temp md.mdp

gmx grompp -f em.mdp -c ${i}.gro -p topol.top -o em.tpr
gmx mdrun -v -deffnm em

gmx grompp -f md.mdp -c em.gro -p topol.top -o md.tpr
gmx mdrun -v -deffnm md

echo "Box-X" | gmx energy -f md.edr -o density.xvg; grep -v "\@" density.xvg | grep -v "\#" > box
gnuplot -e "set fit quiet; fit[250:500] a 'box' via a; set print 'l'; print a"
l=`awk '{print $1}' l`

awk -v l=$l '{if(NR==2403){printf "   %.5f   %.5f   %.5f", l, l, l} else {print $0}}' md.gro > temp; mv temp md.gro

awk '{if($1 == "Tcoupl"){print "Tcoupl = nose-hoover"}
 else if($1 == "Pcoupl"){print "Pcoupl = no"}
 else if($1 == "Pcoupl"){print"gen_vel = no"}
 else {print $0}}' md.mdp > temp; cp temp md.mdp

cp md.gro NPT_${i}.gro
gmx grompp -f md.mdp -c NPT_${i}.gro -p topol.top -o md.tpr
gmx mdrun -v -deffnm md

cp md.gro NVT_${i}.gro

awk '{if($1 == "nsteps"){print "nsteps = 10000"}
 else if($1 == "nstxout"){print "nstxout = 2"}
 else if($1 == "nstvout"){print "nstvout = 2"}
 else {print $0}}' md.mdp > temp; cp temp md.mdp

gmx grompp -f md.mdp -c NVT_${i}.gro -p topol.top -o md.tpr
gmx mdrun -v -deffnm md

echo "0" | gmx trjconv -f md.trr -s md.tpr -o traj_${i}.gro -dt 0.004 -ndec 5 -vel

done