| 查看: 1977 | 回复: 4 | |||
[交流]
【讨论】关于不同z值处的水能形成氢键数的脚本讨论 已有3人参与
|
|
《转载》 氢键判据用的是常用的35度3.5埃的几何判据，当然也可以直接在脚本里改。计算方法是，只要有水分子有一个原子在某一层里，则这个水分子就认为属于这一层的水。对于每一帧，计算属于每一层的水selin与其它物质selbig之间的氢键数，氢键包括了这一层中的水作为氢键受体和供体两种情况，其数目分别为代码中的变量a和b。并且加上这一层水内部之间的氢键数（变量c）的2倍。a+b+2c除以这一层的水数，作为这一帧这一层的每个水的平均氢键数。脚本中循环轨迹中的每一帧，最终得到这一层平均氢键数。nonum变量记录有多少帧在所设定的范围里没有水，这些帧不计算。#后面那行用于调试目的，要考察每帧结果就去掉开头的#。 首先运行下面的脚本，来加载实现这个功能的子程序 proc numhbavg {sel fps1 fps2} { set selin [atomselect top $sel] set selbig [atomselect top "same resid as exwithin 3.5 of $sel"] set k 0.0 set nonum 0 for {set i $fps1} {$i<=$fps2} {incr i} { $selin frame $i $selin update $selbig frame $i $selbig update if {[$selin num]!=0} { set a [llength [lindex [measure hbonds 3.5 35 $selbig $selin] 0]] #measure hbonds 3.5 35 $selbig $selin这个命令返回什么样的值？怎么又先用了lindex、llength命令？ set b [llength [lindex [measure hbonds 3.5 35 $selin $selbig] 0]] #这句跟上句的区别？ set c [llength [lindex [measure hbonds 3.5 35 $selin] 0]] #这句为何只有一个$selin？一般都是像上2句那样有2个$selin $selbig或$selbig $selin？ set k [expr $k+($a+$b+2*$c)*3.0/[$selin num]] #puts "fps:$i $a+$b+[expr 2*$c] num_water:[expr [$selin num]/3.0] avg:[expr $k+($a+$b+2*$c)*3.0/[$selin num]]" } else {incr nonum} } if {[expr $fps2-$fps1+1]==$nonum} {return "no result"} return [expr $k/[expr $fps2-$fps1+1-$nonum]] } 然后下面的循环会调用这个子程序来输出每一层的平均氢键数，这里假设要计算z=4.0~5.6埃的数据，间隔为0.1埃，且限定20 set k [expr $i*0.1] set now [numhbavg "same resid as resname SOL and x<40 and x>20 and y<40 and y>20 and z<[expr $k+0.1] and z>=$k" 100 150] puts [format "%4.2f %4.2f %5.3f" $k [expr $k+0.1] $now] } 我这里随便算一个主要由水构成的普通的体系，水形成的氢键数大概在3.1左右，标准放宽到4.0埃，40度，则可形成氢键数约为3.6。在冰中由于结构十分有序，可形成4个氢键，在液态情况下分子的动能必然造成氢键的破坏，所以结果是很合理的。 帧数范围越大、xy平面越大计算越慢，这个脚本计算速度比较慢，不要一下将范围设得太大。 输出结果如下，前两列代表统计的z值范围，第三列是水的平均氢键数 4.00 4.10 3.015 4.10 4.20 3.161 4.20 4.30 3.201 4.30 4.40 3.159 4.40 4.50 3.237 4.50 4.60 3.130 4.60 4.70 3.201 4.70 4.80 3.338 4.80 4.90 3.201 4.90 5.00 3.041 5.00 5.10 3.122 5.10 5.20 3.182 5.20 5.30 3.160 5.30 5.40 3.309 5.40 5.50 3.189 5.50 5.60 3.327 请各位指点一下，我把自己看不懂的写在上面了，不过可以对程序解释的越详细越好，更便于本人的正确理解！！期待各位高手指点！尤其是bay__gulf ！ 另：氢键几何判据中角度到底指的哪个？vmd中的measure hbonds 的命令中的角度指的是：“the angle formed by the donor, hydrogen, and acceptor must be less than angle from 180 degrees. “ 而我看到的文献Phys. Chem. Chem. Phys. , 2004, 6, 829–835 中氢键几何判据：“The bond angle a between the O–O direction and the molecular O–H direction of the donor, where H is the hydrogen which forms the bond, is lower than a threshold angle aC ”这两个好像不一致啊？ [ Last edited by jianchaoyv on 2010-6-29 at 10:54 ] |
回复此楼» 猜你喜欢
重要消息，中午系统在维护
已经有12人回复
-
哪位老哥知道今年的国自然具体哪一天放榜？
已经有20人回复
-
言之凿凿，8月21号(明天)放榜
已经有14人回复
-
时间戳变了，能看出什么问题？
已经有17人回复
-
只有每年这种时候来逛逛小木虫
已经有21人回复
-
欢迎发来filecode的Mz6后的代码验证其规律
已经有125人回复
-
今天基金会出结果吗？20260819
已经有15人回复
-
filecode，4个jtjc了
已经有13人回复
-
为什么余额宝的年化利率越来越低？主要原因有哪些？
已经有5人回复
-
今天放榜没戏了吧
已经有5人回复
» 本主题相关商家推荐: (我也要在这里推广)
» 本主题相关价值贴推荐，对您同样有帮助:
- 旋蒸溶液，析出的固体在瓶壁上形成规则图案
已经有3人回复
- 从outcar中提取原子坐标的脚本
已经有14人回复
- 虫虫们，分步沉淀法制备复合物的相关文献能分享下么？我不是很懂其中形成机理的
已经有4人回复
- 平面内不同区域不同材料的设置
已经有6人回复
- 水上加氢精修后，氢键的键长总有变化，如何解决？
已经有9人回复
- Discover dyanimics 模拟醇的水溶液时的氢键能为什么总是零？
已经有7人回复
- 水中水分子与水分子之间的氢键键长
已经有5人回复
- 如何加游离水上的氢键？请达人讨论一下吧。有实例，有图有真相。
已经有4人回复
- 【求助】酰胺基与水形成氢键
已经有4人回复
- 【讨论】什么影响氢键的形成？
已经有11人回复
- 【请教】三个结晶水上的氢找不到，可以写氢键么？
已经有14人回复
coolrainbow
木虫 (著名写手)
未来国家冻凉
- 应助: 1 (幼儿园)
- 贵宾: 0.1
- 金币: 3979.5
- 散金: 10
- 红花: 50
- 帖子: 1137
- 在线: 272.3小时
- 虫号: 239560
- 注册: 2006-04-07
- 性别: GG
- 专业: 理论和计算化学
2楼2010-06-29 11:03:21
bay__gulf
金虫 (著名写手)
刘苏州
- 模拟EPI: 8
- 应助: 9 (幼儿园)
- 贵宾: 4.85
- 金币: 2332.8
- 红花: 1
- 帖子: 1344
- 在线: 271小时
- 虫号: 592012
- 注册: 2008-09-03
- 专业: 理论和计算化学
- 管辖: 分子模拟
3楼2010-06-29 11:20:23
|
程序来源于http://hi.baidu.com/sobereva/blo ... 595ef7fc037fbf.html 麻烦最上面程序中还有问题请教！请给解释一下 另外，bay__gulf 那个vmd图画的与vmd的ug中描述好像不符吧？？ [ Last edited by jianchaoyv on 2010-6-29 at 14:43 ] |
4楼2010-06-29 14:15:05
bay__gulf
金虫 (著名写手)
刘苏州
- 模拟EPI: 8
- 应助: 9 (幼儿园)
- 贵宾: 4.85
- 金币: 2332.8
- 红花: 1
- 帖子: 1344
- 在线: 271小时
- 虫号: 592012
- 注册: 2008-09-03
- 专业: 理论和计算化学
- 管辖: 分子模拟
★ ★ ★ ★ ★
小木虫(金币+0.5):给个红包，谢谢回帖交流
ghcacj(金币+4):谢谢 2010-06-29 15:07:14
小木虫(金币+0.5):给个红包，谢谢回帖交流
ghcacj(金币+4):谢谢 2010-06-29 15:07:14
|
hbonds cutoff angle selection1 [selection2]: Find all hydrogen bonds in the given selection(s), using simple geometric criteria. Donor and acceptor must be within the cutoff distance, and the angle formed by the donor, hydrogen, and acceptor must be less than angle from 180 degrees. Only non-hydrogen atoms are considered in either selection. If both selection1 and selection2 are given, the selection1 is considered the donor and selection2 is considered the acceptor. If only one selection is given, all non-hydrogen atoms in the selection are considered as both donors and acceptors. The two selections must be from the same molecule. The function returns three lists; each element in each list corresponds to one hydrogen bond. The first list contains the indices of the donors, the second contains the indices of the acceptors, and the third contains the index of the hydrogen atom in the hydrogen bond. Known Issue: The output of hbonds cannot be considered 100% accurate if the donor and acceptor selection share a common set of atoms. ============= $selbig $selin这个命令返回什么样的值？怎么又先用了lindex、llength命令？ 请在上面找答案 #这句跟上句的区别？ 逐字比对 这句为何只有一个$selin？一般都是像上2句那样有2个$selin $selbig或$selbig $selin？ 请在上面找答案 氢键几何判据中角度到底指的哪个？ 请在上面找答案 vmd 的ug 较长, 全部打印不现实 但8，10，11三章一定要打印下来仔仔细细看 |
5楼2010-06-29 15:05:34
