求助：重复Gromacs官网上的伞形抽样教程得到的summary_distances.dat文件出现问题 - 计算化学公社

#!/bin/bash
# A awk script for creating batch script of Umbrella Sample, GMX Version should >= 5.1
#vvvvvvvvvvvvvvvvvvvvvvvvvvvCHANGE SECTIONvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv#
#From which step to process, support 1 and 2, 1 represents need tpr, xtc/trr and ndx file
#not need summary_distances.dat file; 2 represents only need summary_distances.dat file
# summary_distances.dat includs:
# two columns data
# 1st column: frame number, from gmx trjconv -sep
# 2nd column: distance [nm] or angle [degree], from gmx distance or gmx gangle
step=1 # starting step
# For example, tpr and xtc/trr file name, pull group name
# for calculating center of mass distance/dihedral should be changed by yourself
tpr="pull.tpr" # pull tpr file
xtc="pull.xtc" # pull xtc/trr file
groupA="Chain_A" # first group name for distance pull or dihedral pull(noly)
groupB="Chain_B" # second group name for distance pull
# gmx name, interval, ndx, pull type
# if step=2, only need input file [summary_distances.dat], also can see gromacs tutorials
gmx="gmx" # Maybe gmx, gmx_mpi etc
inputfile="summary_distances.dat" # temp file, do not change
ndx="index.ndx" # ndx file, including groupA, groupB
pull_type="distance" # distance pull, maybe "dihedral"
interval=0.2 # Unit nm for distance, degree for dihedral, you can change it
outdist="outdistance.dat" # distance/dihedral output file
outbatch="do_batch.bsh" # batch output file
#^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^#
if [[ $step -le 1 ]];then
#============================DO NOT CHANGE==============================#
# generate conf*.gro file in order to US simulation
echo -e "\n>>Generating a series of conf*.gro file..."
echo 0 | $gmx trjconv -s $tpr -f $xtc -o conf.gro -sep -quiet > /dev/null 2>&1
# generate distance/dihedral of each frame and change the first column into frame index(from 0 start)
if [[ $pull_type"x" == "distance""x" ]]; then
echo -e "\n>>Calculating center of mass distance between two groups..."
$gmx distance -s $tpr -f $xtc -n $ndx -select "com of group "$groupA" plus com of group "$groupB"" -oall temp.xvg -xvg none -quiet > /dev/null 2>&1
elif [[ $pull_type"x" == "dihedral""x" ]];then
echo -e "\n>>Calculating dihedral ..."
echo $groupA |$gmx gangle -s $tpr -f $xtc -n $ndx -g1 dihedral -oall temp.xvg -xvg none -quiet > /dev/null 2>&1
fi
awk '{printf("%4d %8f\n", NR-1, $2)}' temp.xvg > $inputfile
fi
if [[ $step -le 2 ]];then
# read input file and generate distance/dihedral file and batch file for simulation
echo -e "\n>>Creating batch script according to given interval..."
awk '
BEGIN {
inputfile = "'$inputfile'"; interval = "'$interval'"+0;
outdist = "'$outdist'"; outbatch = "'$outbatch'";
gmx = "'$gmx'"; ndx = "'$ndx'"; pull_type = "'$pull_type'"
print "\n>>Input file: ", inputfile
print ">>Input interval: ", interval, "(nm)"
#check file
Check_File(inputfile)
#read file and obtain data
ll = sampledist(inputfile, interval, sample_indexs, distances)
#output data file, including distance and batch script
if(pull_type == "distance")
printf("%10s %10s %10s\n", "frame", "dist", "d_dist") > outdist
else if(pull_type == "dihedral")
printf("%10s %10s %10s\n", "frame", "angle", "d_angle") > outdist
printf("#!/bin/bash\n\n") > outbatch
for(i = 0; i < ll; i++) {
frame = distances[sample_indexs[i], 0]
dist = distances[sample_indexs[i], 1]
if(i == 0) {delta_dist = "NA"}
else {
prev_dist = distances[sample_indexs[i-1], 1]
delta_dist = dist - prev_dist
}
out_dist(frame, dist, delta_dist, outdist)
out_batch(frame, outbatch)
}
Thanks()
}
#==========================SOME FUNCTION============================#
#output batch file, you can change it if you understand it means
function out_batch(fram, outfile) {
printf("#Do the %s frame simulation\n", fram) > outfile
printf("%s grompp -f npt_umbrella.mdp -c conf%s.gro -r conf%s.gro -p topol.top -n %s -o npt%s.tpr\n",
gmx, fram, fram, ndx, fram) > outfile
printf("%s mdrun -deffnm npt%s\n", gmx, fram) > outfile
printf("%s grompp -f md_umbrella.mdp -c npt%s.gro -r npt%s.gro -t npt%s.cpt -p topol.top -n %s -o umbrella%s.tpr\n",
gmx, fram, fram, fram, ndx, fram) > outfile
printf("%s mdrun -deffnm umbrella%s\n\n", gmx, fram) > outfile
}
#output distance file
function out_dist(fram, dist, delta_dist, outdist) {
printf("%10d %10.3f %10s\n", fram, dist, delta_dist) > outdist
}
#main function for finding appropriate windows
function sampledist(table, interval, sample_indexs, distances, _ARGVEND_, i, n) {
n = 0;
while(getline<table>0) {
distances[n, 0] = $1
distances[n, 1] = $2+0
n++
}
close(table)
current_index = 0
sample_indexs[0] = 0
iindex = 0
while(current_index < n) {
target_distance = distances[current_index, 1] + interval
tt = 0
for(i = current_index; i < n; i++) {
onward[tt++] = abs(target_distance - distances[i, 1])
}
mim_value = min(onward)
for(j in onward) {
if(onward[j] == mim_value) {
next_index = j + current_index
break
}
}
if(current_index == next_index) {break}
else {
sample_indexs[++iindex] = next_index
current_index = next_index
}
delete onward
}
return length(sample_indexs)
}
#abs function
function abs(a) { return a>0? a : -a }
#find minimum
function min(x, _ARGVEND_, value, i) {
value = 99999
for(i in x) {
value = value < x[i] ? value : x[i]
}
return value
}
#output information
function Thanks() {
print "\n=============RESULT=============="
print outdist " has been created!"
print outbatch " has been created!"
print "================================="
}
# check file
function Check_File(file)
{
if (getline < file < 0)
{
printf("\nError! file [%s] not found!", file)
exit 1
}
else
close(file)
}
'
fi

复制代码