计算化学公社

标题: vmd脚本求助 分析分子周围的粒子分布 [打印本页]
作者
Author: diaok    时间: 2015-12-4 16:37
标题: vmd脚本求助 分析分子周围的粒子分布
本帖最后由 diaok 于 2015-12-4 16:46 编辑

这个分析粒子分布的脚本是仿照sob老师的氢键分析写的
  1. proc rgwithin {sel sel2 fps1 fps2} {
  2. set sum 0.0
  3. set cut  6  
  4. #  0.1:0.1:6.1
  5. set piece  0.1
  6. set in  0
  7. for {set j 1} {$in<=$cut} {incr j} {
  8. set in  [expr $j*0.1]
  9. set selin [atomselect top  "$sel && within $in of $sel2"]
  10. set sum  0
  11. for {set i $fps1} {$i<=$fps2} {incr i} {
  12. $selin frame $i
  13. $selin update
  14. set tmp [$selin  num]
  15. set sum [expr $sum+$tmp]
  16. }
  17. #echo  $in
  18. lappend rg    [expr $sum/[expr $fps2-$fps1+1.0]]
  19. $selin  delete
  20. }
  21. return  $rg
  22. }
复制代码


但是这个计算速度特别慢，还会导致电脑特别卡
我的体系里面分析600帧就用了6分钟
而且不知道是什么原因进行第二次计算会越来越慢
脚本中明明已经用了$selin  delete

用什么方法计算一个原子到另一个大分子的最近距离比较快？
用gromacs或者vmd来实现都行


作者
Author: ruanyang    时间: 2015-12-4 18:16
tcl脚本计算本省就比较慢，你改成向量计算可能会快一点！
作者
Author: sobereva    时间: 2015-12-4 18:23
把轨迹转换成xyz格式，自行用Fortran或C代码来载入和分析会很快
注意运行过程中看看内存占用率，有没有越来越高的问题
作者
Author: diaok    时间: 2015-12-4 20:16
sobereva 发表于 2015-12-4 18:23
把轨迹转换成xyz格式，自行用Fortran或C代码来载入和分析会很快
注意运行过程中看看内存占用率，有没有越 ...

sob老师回复神速啊

第二次提取的时候主要是cpu占用升到100%了
vmd的内存占用100m一直没怎么变

很奇怪批处理为什么会越来越慢，6分钟一次还算可以接受
而且删除了分子重新载入都没效果，只能重启vmd


一次批处理的代码
  1. set outfile [open {W:\susedata\kbff\na8.txt} w]
  2. set  tmp  [rgwithin   "type NA"  "resname BGC"   0 499  ]
  3. puts $outfile $tmp
  4. set  tmp  [rgwithin   "type OO"  "resname BGC"   0 499   ]
  5. puts $outfile $tmp
  6. set  tmp  [rgwithin   "type HH"  "resname BGC"   0 499  ]
  7. puts $outfile $tmp
  8. set  tmp  [rgwithin   "type OW"  "resname BGC"   0 299     ]
  9. puts $outfile $tmp
  10. close $outfile
复制代码


Fortran和c语言不怎么熟悉
用matlab的textscan可行吗？轨迹500m，有30000个frame
作者
Author: diaok    时间: 2015-12-4 20:30
ruanyang 发表于 2015-12-4 18:16
tcl脚本计算本省就比较慢，你改成向量计算可能会快一点！

有什么好的方法可以在gromacs或者vmd里实现吗？
gromacs我只想到求每个原子与大分子的所有距离，求最小值
然后对所有帧所有原子做直方图
作者
Author: sobereva    时间: 2015-12-4 22:36
gromacs里自带了g_mindist工具，可以看看是否合用
作者
Author: diaok    时间: 2015-12-6 16:53
sobereva 发表于 2015-12-4 22:36
gromacs里自带了g_mindist工具，可以看看是否合用

非常感谢

批量处理可以提取出来
  1. mkdir oho
  2. file=oho/oho.ndx
  3. cat index.ndx >  $file
  4. for i in `seq 1 45`
  5. do
  6. echo [oho$i] >> $file
  7. ((t=$i+4672-2))
  8. echo $t >> $file
  9. done
  10. #批量生成所有离子的group

  12. for i in `seq 1 45`
  13. do
  14. gmx mindist -s md.tpr  -f md.xtc  -n oho/oho.ndx -od oho/o$i.xvg  << EOF
  15. oho$i
  16. 2   #用字母匹配会出错，原因未知
  17. EOF
  18. done
  19. # 批量提取每个粒子到大分子的最小距离

  21. 获得了所有距离
  22. 剩下的就交给matlab去循环读取了
复制代码
作者
Author: ruanyang    时间: 2015-12-6 18:17
强，脚本写的很牛逼
作者
Author: compume    时间: 2015-12-17 21:40
好厉害，我是刚入行菜鸟一枚，想问一下，直接用python之类的语言我应该怎么从trr中读取出每个分子或者原子的位置信息呢？
以为trr是二进制文件，不知道怎么看啊……
作者
Author: sobereva    时间: 2015-12-18 00:05
compume 发表于 2015-12-17 21:40
好厉害，我是刚入行菜鸟一枚，想问一下，直接用python之类的语言我应该怎么从trr中读取出每个分子或者原子 ...

自己用python读trr很难实现，建议你把trr用VMD转换为xyz轨迹，是二进制文件，这种轨迹格式非常简单，之后自己写程序分析处理很方便。当然体积比trr大不少。
作者
Author: compume    时间: 2015-12-18 10:01
sobereva 发表于 2015-12-18 00:05
自己用python读trr很难实现，建议你把trr用VMD转换为xyz轨迹，是二进制文件，这种轨迹格式非常简单，之后 ...

THX，可是我的trr文件有几个G……这样处理岂不会很慢
作者
Author: sobereva    时间: 2015-12-18 12:18
compume 发表于 2015-12-18 10:01
THX，可是我的trr文件有几个G……这样处理岂不会很慢

没法避免
作者
Author: compume    时间: 2015-12-18 14:53
sobereva 发表于 2015-12-18 12:18
没法避免

好嘞，多谢sob老师~
作者
Author: zhangj    时间: 2015-12-19 20:56
您好，请问sob老师氢键分析的脚本在哪里可以找到，非常感谢
作者
Author: tjuchan    时间: 2015-12-19 21:28
没有考虑周期边界是不是会少很多氢键？
作者
Author: sobereva    时间: 2015-12-19 22:03
zhangj 发表于 2015-12-19 20:56
您好，请问sob老师氢键分析的脚本在哪里可以找到，非常感谢

计算不同z位置水能形成氢键数的VMD Tcl脚本
http://sobereva.com/54
作者
Author: zhangj    时间: 2015-12-23 10:48
sobereva 发表于 2015-12-19 22:03
计算不同z位置水能形成氢键数的VMD Tcl脚本
http://sobereva.com/54

非常非常感谢sob老师！



欢迎光临 计算化学公社 (http://bbs.keinsci.com/) Powered by Discuz! X3.3