|||
【静电修改】
直接去改top文件即可。
【vdw修改】
由于vdw是nonbond相互作用,不能通过直接修改top文件来达到。
需要在include"**/forcefield.itp"和[moleculetypes]之间加入重新定义的原子类型。
比如修改某原子类型与Na离子的vdw相互作用(离子通道蛋白)。若不是所有的该原子类型全部都改,那么需要单独定义一个新的原子类型,进而需要重新定义与该原子相关的ffnonbonded.itp和ffbonded.itp中的相互作用参数。
一个例子脚本如下:
#need file ffbonded.itp kcsA_vdw.top, obtain add.top lbond=6337 lpair=12298 langle=27727 ldihedral1=38499 ldihedral2=54241 lmap=55111 echo "processing the kcsA_vdw.top file" awk -v lbond=$lbond -v lpair=$lpair '{if(NR>lbond+1 && NR< lpair-2) print}' kcsA_vdw.top > bond.top awk -v langle=$langle -v ldihedral1=$ldihedral1 '{if(NR>langle+1 && NR< ldihedral1-2) print}' kcsA_vdw.top >angle.top awk -v ldihedral1=$ldihedral1 -v ldihedral2=$ldihedral2 '{if(NR>ldihedral1+1 && NR <ldihedral2-2) print }' kcsA_vdw.top > dihedral1.top awk -v ldihedral2=$ldihedral2 -v lmap=$lmap '{if(NR>ldihedral2+1 && NR <lmap-2) print}' kcsA_vdw.top >dihedral2.top echo "extract the change lines" for cx in 899 900 901 902 906 907 908 909 2368 2369 2370 2371 2375 2376 2377 2378 3837 3838 3839 3840 3844 3845 3846 3847 5306 5307 5308 5309 5313 5314 5315 5316 do grep $cx bond.top >> bond_change.top done for cx in 899 900 901 902 906 907 908 909 2368 2369 2370 2371 2375 2376 2377 2378 3837 3838 3839 3840 3844 3845 3846 3847 5306 5307 5308 5309 5313 5314 5315 5316 do grep $cx angle.top >> angle_change.top done for cx in 899 900 901 902 906 907 908 909 2368 2369 2370 2371 2375 2376 2377 2378 3837 3838 3839 3840 3844 3845 3846 3847 5306 5307 5308 5309 5313 5314 5315 5316 do grep $cx dihedral1.top >> dihedral1_change.top done for cx in 899 900 901 902 906 907 908 909 2368 2369 2370 2371 2375 2376 2377 2378 3837 3838 3839 3840 3844 3845 3846 3847 5306 5307 5308 5309 5313 5314 5315 5316 do grep $cx dihedral2.top >> dihedral2_change.top done mv bond_change.top bond.top mv angle_change.top angle.top mv dihedral1_change.top dihedral2.top mv dihedral2_change.top dihedral2.top echo "obtain the index and name of atom, then substitue the numbers with name" lstart=71 lend=6336 awk -v lstart=$lstart -v lend=$lend '{if(NR>=lstart && NR<=lend) {if($1!=";") print $1,$2}}' kcsA_vdw.top > refer.top awk 'BEGIN{printf("sed \47")} {printf("s/ %s / %s /g;",$1,$2)}END{printf("\47 bond.top >temp\nmv temp bond.top")}' refer.top > refer.sh echo "do the substitution of bond,may need a little time" bash refer.sh awk 'BEGIN{printf("sed \47")} {printf("s/ %s / %s /g;",$1,$2)}END{printf("\47 angle.top >temp\nmv temp angle.top")}' refer.top > refer.sh echo "do the substitution of angle,may need a little time" bash refer.sh awk 'BEGIN{printf("sed \47")} {printf("s/ %s / %s /g;",$1,$2)}END{printf("\47 dihedral1.top >temp\nmv temp dihedral1.top")}' refer.top > refer.sh echo "do the substitution of dihedral1,may need a little time" bash refer.sh awk 'BEGIN{printf("sed \47")} {printf("s/ %s / %s /g;",$1,$2)}END{printf("\47 dihedral2.top >temp\nmv temp dihedral2.top")}' refer.top > refer.sh echo "do the substitution of dihedral2,may need a little time" bash refer.sh rm refer.top rm refer.sh echo "finish the substituion" awk '/CX|HX/{printf("%s %s ",$1,$2);gsub(/CX/,"CA");gsub(/HX/,"HP");printf("%s %s\n",$1,$2)}' bond.top >temp sort -n temp |uniq > bond.top awk '/CX|HX/{printf("%s %s %s ",$1,$2,$3);gsub(/CX/,"CA");gsub(/HX/,"HP");printf("%s %s %s\n",$1,$2,$3)}' angle.top >temp sort -n temp |uniq > angle.top awk '/CX|HX/{printf("%s %s %s %s ",$1,$2,$3,$4);gsub(/CX/,"CA");gsub(/HX/,"HP");printf("%s %s %s %s\n",$1,$2,$3,$4)}' dihedral1.top >temp sort -n temp |uniq > dihedral1.top awk '/CX|HX/{printf("%s %s %s %s ",$1,$2,$3,$4);gsub(/CX/,"CA");gsub(/HX/,"HP");printf("%s %s %s %s\n",$1,$2,$3,$4)}' dihedral2.top >temp sort -n temp |uniq > dihedral2.top echo "extracting the lines from ffbonded.itp" lbondt=1 lconstr=187 langlet=227 ldihed1=701 ldihed2=1303 awk -v lbondt=$lbondt -v lconstr=$lconstr '{if(NR>lbondt+1 && NR<lconstr) print $1,$2,$3,$4,$5}' ffbonded.itp > bond.itp awk -v lbondt=$lbondt -v lconstr=$lconstr '{if(NR>lbondt+1 && NR<lconstr) print $2,$1,$3,$4,$5}' ffbonded.itp >> bond.itp sort -n bond.itp | uniq > temp mv temp bond.itp awk -v langlet=$langlet -v ldihed1=$ldihed1 '{if(NR>langlet && NR<ldihed1) print $1,$2,$3,$4,$5,$6,$7,$8}' ffbonded.itp > angle.itp awk -v langlet=$langlet -v ldihed1=$ldihed1 '{if(NR>langlet && NR<ldihed1) print $3,$2,$1,$4,$5,$6,$7,$8}' ffbonded.itp >> angle.itp sort -n angle.itp | uniq > temp mv temp angle.itp awk -v ldihed1=$ldihed1 -v ldihed2=$ldihed2 '{if(NR>ldihed1 && NR<ldihed2) print $1,$2,$3,$4,$5,$6,$7,$8}' ffbonded.itp > dihedral1.itp awk -v ldihed1=$ldihed1 -v ldihed2=$ldihed2 '{if(NR>ldihed1 && NR<ldihed2) print $4,$3,$2,$1,$5,$6,$7,$8}' ffbonded.itp >> dihedral1.itp sort -n dihedral1.itp | uniq > temp mv temp dihedral1.itp awk -v ldihed2=$ldihed2 '{if(NR>ldihed2+1) print $1,$2,$3,$4,$5,$6,$7}' ffbonded.itp > dihedral2.itp awk -v ldihed2=$ldihed2 '{if(NR>ldihed2+1) print $4,$3,$2,$1,$5,$6,$7}' ffbonded.itp >> dihedral2.itp sort -n dihedral2.itp | uniq > temp mv temp dihedral2.itp echo "processing bondtypes" echo "[ bondtypes ]" > add.top wc -l bond.top > lt.txt lt=`awk '{print $1}' lt.txt` i=0 until [ $i -eq $lt ]; do pair=`awk -v i=$i '{if(NR==i+1) printf("%s %s",$3,$4)}' bond.top` grep "$pair" bond.itp > temp awk -v i=$i '{if(NR==i+1) print $1,$2}' bond.top > temp1 paste temp temp1 > to_add awk '{print $6,$7,$3,$4,$5}' to_add >> add.top i=$(($i+1)) done rm bond.top rm bond.itp echo "processing angletypes" echo "[ angletypes ]" >> add.top wc -l angle.top > lt.txt lt=`awk '{print $1}' lt.txt` i=0 until [ $i -eq $lt ]; do pair=`awk -v i=$i '{if(NR==i+1) printf("%s %s %s",$4,$5,$6)}' angle.top` grep "$pair" angle.itp > temp pair1=`cat temp` if [ -n "$pair1" ]; then awk -v i=$i '{if(NR==i+1) print $1,$2,$3}' angle.top > temp1 paste temp temp1 > to_add awk '{print $9,$10,$11,$4,$5,$6,$7,$8}' to_add >> add.top fi i=$(($i+1)) done rm angle.top rm angle.itp echo "processing dihedraltypes" echo "[ dihedraltypes ]" >> add.top wc -l dihedral1.top > lt.txt lt=`awk '{print $1}' lt.txt` i=0 until [ $i -eq $lt ]; do pair=`awk -v i=$i '{if(NR==i+1) printf("%s %s %s %s",$5,$6,$7,$8)}' dihedral1.top` grep "$pair" dihedral1.itp > temp pair1=`cat temp` if [ -n "$pair1" ]; then awk -v i=$i '{if(NR==i+1) print $1,$2,$3,$4}' dihedral1.top > temp1 paste temp temp1 > to_add awk '{print $9,$10,$11,$12,$5,$6,$7,$8}' to_add >> add.top fi pair=`awk -v i=$i '{if(NR==i+1) printf("X %s %s X",$6,$7)}' dihedral1.top` grep "$pair" dihedral1.itp > temp pair1=`cat temp` if [ -n "$pair1" ]; then awk -v i=$i '{if(NR==i+1) printf("X %s %s X\n",$2,$3)}' dihedral1.top > temp1 paste temp temp1 > to_add awk '{print $9,$10,$11,$12,$5,$6,$7,$8}' to_add >> add.top fi i=$(($i+1)) done rm dihedral1.top rm dihedral1.itp echo "[ dihedraltypes ]" >> add.top wc -l dihedral2.top > lt.txt lt=`awk '{print $1}' lt.txt` i=0 until [ $i -eq $lt ]; do pair=`awk -v i=$i '{if(NR==i+1) printf("%s %s %s %s",$5,$6,$7,$8)}' dihedral2.top` grep "$pair" dihedral2.itp > temp pair1=`cat temp` if [ -n "$pair1" ]; then awk -v i=$i '{if(NR==i+1) print $1,$2,$3,$4}' dihedral2.top > temp1 paste temp temp1 > to_add awk '{print $8,$9,$10,$11,$5,$6,$7}' to_add >> add.top fi pair=`awk -v i=$i '{if(NR==i+1) printf("X %s %s X",$6,$7)}' dihedral2.top` grep "$pair" dihedral2.itp > temp pair1=`cat temp` if [ -n "$pair1" ]; then awk -v i=$i '{if(NR==i+1) printf("X %s %s X\n",$2,$3)}' dihedral2.top > temp1 paste temp temp1 > to_add awk '{print $8,$9,$10,$11,$5,$6,$7}' to_add >> add.top fi i=$(($i+1)) done rm dihedral2.top rm dihedral2.itp rm temp1 rm lt.txt
这里面包含了awk,grep,sed以及bash的if, until, for语句,可以用来参考。
for循环的使用,完全可以替换掉until:
http://gaodr.blog.163.com/blog/static/10461500820093793933887/
建议使用类C的语法:
for ((i=0; i<5; i++)) do ... done
Archiver|手机版|科学网 ( 京ICP备07017567号-12 )
GMT+8, 2024-11-22 23:54
Powered by ScienceNet.cn
Copyright © 2007- 中国科学报社