24小时热门版块排行榜    

北京石油化工学院2026年研究生招生接收调剂公告
查看: 2074  |  回复: 15

jianchaoyv

金虫 (小有名气)

[交流] 【求助】向bay__gulf求助!! 已有4人参与

目前在做碳纳米管中水的相变,请教一下怎样写.tcl脚本来计算radial density of oxygen 和 orientation distributions of water molecules within CNTs?
虽然这个问题对你可能比较简单,但对我来说却很棘手,望赐教!
同时恭喜bay__gulf发JPCA!
也欢迎其他高手指教!
回复此楼
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
回帖支持 ( 显示支持度最高的前 50 名 )

bay__gulf

金虫 (著名写手)

刘苏州

★ ★ ★ ★ ★ ★ ★ ★ ★ ★
jianchaoyv(金币+10):谢谢指点 2010-06-22 18:29:11
ghcacj(金币+10):谢谢 2010-06-22 19:55:31
radial density, 也就是径向分布, 很常用, vmd 也提供了一个工具可以直接计算
就是measure gofr, 参考 http://www.ks.uiuc.edu/Research/vmd/current/ug/node133.html, 以及一个gui 界面, 参考 http://www.ks.uiuc.edu/Research/vmd/plugins/gofrgui/

至于取向分布, 用户较少只能自己写了,
vmd 的强大之处不在于内置了多少工具, 更不是做多漂亮的图形
而是提供了一个原子(集团)为数据结构的编程环境, 可以按照自己的需要实现任意的功能
下面是当时用的代码, 现在看来很幼稚了但拿它入门还是不错的

# count the dipolar-z angle and Legendre 2nd order polynomial in a pore,
# make sure the psf/dcd be the top molecular

set outf [open theta.dat w]

set numFrame [molinfo top get numframes]

set lgd 0
for { set fr 1 } { $fr < $numFrame } {incr fr } {
   set wt [atomselect top "water and name OH2 and z < 7.6 and z > -7.6 " frame $fr]
   set n [$wt num]
   set theta 0

   foreach resid [$wt get resid] {

      set wto  [atomselect top "resid $resid and name OH2" frame $fr]
      set wth1 [atomselect top "resid $resid and name H1" frame $fr]
      set wth2 [atomselect top "resid $resid and name H2" frame $fr]

      set coordo  [measure center  $wto]
      set coordh1 [measure center $wth1]
      set coordh2 [measure center $wth2]

      set coordh  [vecadd $coordh1 $coordh2]
      set coordh  [vecscale 0.5 $coordh]
      set vecd    [vecsub $coordh $coordo]

      set thetai  [expr acos ([expr [lindex $vecd 2] / [veclength $vecd]])]
      set thetai  [expr $thetai * 180.0 / 3.1416]
      set theta   [expr $theta + $thetai]

      $wto  delete
      $wth1 delete
      $wth2 delete
   }

   set theta [expr $theta/$n ]
   set costheta [expr cos ([expr $theta*3.1416/180.0])]
   set lgd [expr $lgd + [expr $costheta * $costheta]]
   puts $fr
   puts $outf [format "%d %f" $fr $theta]

   $wt delete
}
   set lgd [expr $lgd / $numFrame * 1.5 - 0.5]
   puts "3*<(cosθ)^2> - 1 = $lgd"
close $outf
2楼2010-06-22 16:44:02
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
普通回帖

jianchaoyv

金虫 (小有名气)

能否具体讲一下两种操作方法的步骤或给一个具体例子?你给的连接看了,但还是不知道怎么去操作。麻烦你帮忙写一个.tcl脚本,谢谢!

[ Last edited by jianchaoyv on 2010-6-22 at 20:24 ]
3楼2010-06-22 19:41:03
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

jianchaoyv

金虫 (小有名气)

请问上面提供的脚本中:
  set coordh  [vecadd $coordh1 $coordh2]
      set coordh  [vecscale 0.5 $coordh]
      set vecd    [vecsub $coordh $coordo]

      set thetai  [expr acos ([expr [lindex $vecd 2] / [veclength $vecd]])]
还有最后set lgd [expr $lgd / $numFrame * 1.5 - 0.5]
这几句看不懂,能不能解释一下?谢谢!

[ Last edited by jianchaoyv on 2010-6-23 at 09:56 ]
4楼2010-06-23 09:50:22
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

bay__gulf

金虫 (著名写手)

刘苏州

★ ★ ★ ★ ★
ghcacj(金币+5):谢谢 2010-06-24 10:06:00
set coordh  [vecadd $coordh1 $coordh2] ;#两个H的坐标之和
set coordh  [vecscale 0.5 $coordh] ;#两个H的平均坐标, 即中点
set vecd    [vecsub $coordh $coordo] ;#O原子到H中点的矢量, 即水的偶记(方向)
set thetai  [expr acos ([expr [lindex $vecd 2] / [veclength $vecd]])]
# 水的偶极u跟z轴(0 0 1)的夹角, cos(theta) = u(z)/|u|

还有最后set lgd [expr $lgd / $numFrame * 1.5 - 0.5]
勒让德项 (3-1)/2,
参考 http://www.mdbbs.org/viewthread.php?tid=11530
5楼2010-06-24 09:13:13
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

jianchaoyv

金虫 (小有名气)

请问vecadd、vecscale、vecsub、veclength这些命令在哪里可以学习一下?根本不知道它们的用途?请给一个关于这类命令的全面的学习资料!!!十分感谢!

另外,这句set thetai  [expr acos ([expr [lindex $vecd 2] / [veclength $vecd]])]的意思 #水的偶极u跟z轴(0 0 1)的夹角, cos(theta) = u(z)/|u|看不懂啊?
我的理解是:lindex $vecd 2 #是O原子的位置矢量
            veclength $vecd #O原子到H中点的矢量|u|
    所以,set thetai  [expr acos ([expr [lindex $vecd 2] / [veclength   $vecd]])] 并不是真正的#水的偶极u跟z轴(0 0 1)的夹角, cos(theta) = u(z)/|u|
正在研究你给的这个脚本中!

[ Last edited by jianchaoyv on 2010-6-25 at 16:12 ]
6楼2010-06-25 15:28:58
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

bay__gulf

金虫 (著名写手)

刘苏州

★ ★ ★
jianchaoyv(金币+3): 2010-06-26 08:49:39
lei0736(金币+3):谢谢 2010-06-26 09:17:42
请问vecadd、vecscale、vecsub、veclength这些命令在哪里可以学习一下?根本不知道它们的用途?请给一个关于这类命令的全面的学习资料!!!十分感谢!
====
vmd ug 的第10章讲的是Vectors and Matrices
这个ug 仅仅是vmd 特有的和常用的tcl 命令, 更深入了解需要Tcl 的资料

lindex $vecd 2 #是O原子的位置矢量
===
vecd 是偶极矢量,
lindex $vecd 2 是偶极的z 分量
veclength $vecd 是偶极的长度

很荣幸自己的工作被别人分析和引用
7楼2010-06-25 16:34:06
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

coolrainbow

木虫 (著名写手)

未来国家冻凉

★ ★
jianchaoyv(金币+2): 2010-06-26 08:49:58
lei0736(金币+2):谢谢 2010-06-26 09:18:03
引用回帖:
Originally posted by jianchaoyv at 2010-06-25 15:28:58:
请问vecadd、vecscale、vecsub、veclength这些命令在哪里可以学习一下?根本不知道它们的用途?请给一个关于这类命令的全面的学习资料!!!十分感谢!

另外,这句set thetai  [expr acos ([expr [lindex $vec ...

赞!bay_gulf 可是这里水平极高而又热心助人的高手之一啊

建议你可以用一天的时间学学TCL语言,再看看vmd中特有的命令,再学习几个namd中已有的tcl脚本作为范例,基本上你就可以自己写脚本了

纯tcl推荐学习:Practical Programming in Tcl and Tk, Fourth Edition
技术博客:http://hi.baidu.com/coolrainbow/blog
8楼2010-06-25 21:08:01
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

jianchaoyv

金虫 (小有名气)

bay__gulf确实是牛人,也很耐心!
我把你那段程序的理解写下来:
程序开始到此已理解
set theta [expr $theta/$n ]  #对每帧所有水分子偶极距求平均,另外,这个theta角是跟Z轴正半轴的夹角吧?而看到的文献中大都是算的与Z轴负半轴的夹角,感觉没什么差别,这样理解多吗?
   set costheta [expr cos ([expr $theta*3.1416/180.0])]#求cos(theta)
   set lgd [expr $lgd + [expr $costheta * $costheta]]#用于求不同帧的cos(theta)^2求和
   puts $fr #输出第几帧
   puts $outf [format " %f"  $theta] 输出该帧下水分子偶极距的平均值
   $wt delete
}
   set lgd [expr $lgd / $numFrame * 1.5 - 0.5] #对所有帧的2阶勒让德函数求平均
   puts "3*<(cosθ)^2> - 1 = $lgd" #输出P2的平均值
close $outf

对上面的解释若有错误请指正,还想请教一下算P2用来表征什么的?不知道算P2的目的。
若set theta [expr $theta/$n ]  #对每帧所有水分子偶极距求平均,理解正确的话,也就默认了每帧里的所有水分子的偶极距方向都一样了,感觉同一帧不同位置的水分子偶极距还是有可能不同的,有这种可能吗?

[ Last edited by jianchaoyv on 2010-6-26 at 10:09 ]
9楼2010-06-26 09:48:01
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

bay__gulf

金虫 (著名写手)

刘苏州

★ ★ ★ ★
小木虫(金币+0.5):给个红包,谢谢回帖交流
lei0736(金币+3):谢谢 2010-06-26 11:43:25
引用回帖:
Originally posted by jianchaoyv at 2010-06-26 09:48:01:
bay__gulf确实是牛人,也很耐心!
我把你那段程序的理解写下来:
程序开始到此已理解
set theta [expr $theta/$n ]  #对每帧所有水分子偶极距求平均,另外,这个theta角是跟Z轴正半轴的夹角吧?而看到的文献中 ...

理解很正确
引用回帖:
也就默认了每帧里的所有水分子的偶极距方向都一样了,感觉同一帧不同位置的水分子偶极距还是有可能不同的,有这种可能吗?

为使程序好理解, 修改了点数据
water and name OH2 and z < 7.6 and z > -7.6
这是整个管子中的水
改变z 的上下限, 可以求得一小段区域中, 水的平均偶极方向
多做几个位置, 可以求得不同位置的数值

[ Last edited by bay__gulf on 2010-6-26 at 10:37 ]
10楼2010-06-26 10:36:53
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
相关版块跳转 我要订阅楼主 jianchaoyv 的主题更新
普通表情 高级回复 (可上传附件)
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[考研] 336材料与化工085600求调剂 +8 水星记infp 2026-04-05 8/400 2026-04-06 06:09 by houyaoxu
[考研] 322求调剂 +3 嗯哼哼恒 2026-04-05 3/150 2026-04-05 19:52 by nepu_uu
[考研] 考研调剂生寻找导师 +3 顾瞻考研啊 2026-04-05 3/150 2026-04-05 18:18 by 啵啵啵0119
[考研] 调剂 +5 好好读书。 2026-04-01 5/250 2026-04-05 17:54 by liucky
[考研] 296求调剂 +3 汪!?! 2026-04-05 5/250 2026-04-05 17:38 by 蓝云思雨
[考研] 材料专硕 调剂 +14 CXN123456 2026-04-03 14/700 2026-04-05 17:18 by Hdyxbekcb
[考研] 272求调剂 +4 电气李 2026-04-05 4/200 2026-04-05 10:41 by lbsjt
[考研] 生物工程求调剂 +6 喜欢还是不甘心 2026-04-05 6/300 2026-04-05 10:28 by 唐沐儿
[考研] 282电子信息0854专硕调剂 +4 202451007219 2026-04-02 6/300 2026-04-04 21:55 by laoshidan
[考研] 359求调剂 +7 hhhhaaaa$ 2026-04-04 7/350 2026-04-04 18:49 by imissbao
[考研] 22408求调剂 354分 可跨专业 +3 hannnnnnn 2026-04-04 3/150 2026-04-04 14:35 by 土木硕士招生
[考研] 334求调剂 +8 曾仰之 2026-04-03 8/400 2026-04-04 11:16 by w_xuqing
[考研] 一志愿华南师范大学-22408计算机-292分-求华南师范大学调剂 +4 爱读书的小鳄鱼 2026-04-02 4/200 2026-04-02 18:35 by 求调剂zz
[考研] 08开头看过来!!! +4 wwwwffffff 2026-03-31 6/300 2026-04-02 11:42 by 均值回归
[考研] 085410 一志愿211 22408分数359求调剂 +3 123456789qw 2026-03-31 4/200 2026-04-02 00:06 by 义文wang
[考研] 298求调剂 +4 什么是胖头鱼 2026-03-30 6/300 2026-04-01 22:06 by 客尔美德
[考研] 265求调剂 +11 yelck 2026-04-01 12/600 2026-04-01 19:12 by 549790059
[考研] 08工科275求调剂,可跨考。 +5 AaAa7420 2026-03-31 5/250 2026-04-01 15:21 by 159357hjz
[考研] 311求调剂 +10 李芷新1 2026-03-31 10/500 2026-04-01 14:38 by chenqifeng666
[考研] 323分 食品与营养调剂 +3 嘿ooo 2026-03-31 3/150 2026-03-31 09:38 by longlotian
信息提示
请填处理意见