|
|
本帖最后由 lyj714 于 2021-3-29 19:45 编辑
如果你用的是最原始的那个perl脚本,那是因为以前的这个脚本有些问题。如果是用的现在原英文教程中的那个bash脚本就没啥问题,不过速度慢。本人提供一个我自用的脚本,可以一步得到US的批处理脚本
- #!/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
复制代码 |
|