计算化学公社

标题: 求助：博文http://sobereva.com/121操作遇到问题 [打印本页]

作者
Author: liaorongbao 时间: 2020-11-27 17:32
标题: 求助：博文http://sobereva.com/121操作遇到问题
本帖最后由 liaorongbao 于 2020-11-27 18:37 编辑

http://sobereva.com/121是指导如何在.pdb文件中设好原子间键连方式然后用博文中的脚本读取后再让vmd软件原版显示出来由.pdb文件定义的键连方式图案的方法。我采用两个H2O分子作为练习的案例，两个水分子距离比较远，本来是没有氢键的。在练习中我假设需要用VMD软件显示第一个水分子的O与第二个水分子的H之间的氢键。
用gview制作2H2O.pdb文件如下：
REMARK 1 File created by GaussView 6.0.16
HETATM 1  O          0    -1.259  -1.631  -1.251                      O
HETATM 2  H          0    -0.299  -1.631  -1.251                      H
HETATM 3  H          0    -1.580  -0.726  -1.251                      H
HETATM 4  O          0    -3.541  -4.720  -0.665                      O
HETATM 5  H          0    -3.165  -5.154  -1.434                      H
HETATM 6  H          0    -2.870  -4.179  -0.242                      H
END
CONECT 1 2 3 6
CONECT 2 1
CONECT 3 1
CONECT 4 5 6
CONECT 5 4
CONECT 6 4 1
原子个数有6个，所以博文中的set natm 6应该是没问题的。
博文中set nskipline 后边跟着的数字比natm大1，但容易发现CONECT显示在第9行。
第一个问题：当前情况采用set nskipline 7还是set nskipline 8呢？

现在我已经在2H2O.pdb文件中设定第一个水分子的O原子与第二个水分子的H原子有单键(就当作代表氢键吧)。
博文中有个Tcl脚本，按2H2O.pdb体系设定脚本的原子数natm=6和跳过的行数nskipline=7(是不是应该写成8？)，然后博文中.pdb文件的路径和名称改写成我的.pdb文件的路径和名称。
即为把原博文路径语句set rdpdbcon [open "D://CM//my program//Multiwfn//vtx.pdb" r]
改为我的2H2O.pdb文件的路径set rdpdbcon [open "E://Program Files (x86)//University of Illinois//VMD//2H2O.pdb" r]
这样改没问题吧？
第二个问题：目录符号是//而不是\，百度了一下编程语言中用//，是为了避免转义，是不是这样？

然后我就双击vmd.exe，把2H2O.pdb拖入 VMD main窗口，这个时候vmd.exe是不显示2H2O.pdb中设定的氢键的，因为Tcl脚本还没发力。
第三个问题：博文说：“下面这个Tcl脚本就是根据上述原理通过读取CONECT字段来设定VMD中原子间连接关系。拷贝到Tcl窗口下即可运行”。这句话也许对多数人不是问题但我却不知怎么操作
我的方法如下：
复制博文的脚本，按我的文件情况修改脚本。然后把修改好的脚本制作成connect.bat文件，把这个文件放到vmd.exe所在目录里。
在vmd.exe的黑色命令窗口执行source connect.bat，显示cannot read "self": no such variable
脚本语句[size=11.375px]scan [string range $line 6 84] "%d %d %d %d %d %d %d %d %d %d %d %d %d" self cn(1) cn(2) cn(3) cn(4) cn(5) cn(6) cn(7) cn(8) cn(9) cn(10) cn(11) cn(12)中有个self，是不是没有被赋值？也不知道怎么回事

如下是我的connect.bat文件
#Programmed by Sobereva, 2012-Jan-25
set natm 6
set nskipline 7
#itype=0 means reset connectivity by user-defined list, =1 means add self-defined connectivity list to the original one
set itype 1
set rdpdbcon [open "E://Program Files (x86)//University of Illinois//VMD//2H2O.pdb" r]
#set rdpdbcon [open "E:\Program Files (x86)\University of Illinois\VMD\2H2O.pdb" r]
#E:\Program Files (x86)\University of Illinois\VMD
#Skip other lines
for {set i 1} {$i<=$nskipline} {incr i} {
gets $rdpdbcon line
}

#Cycle each atom
set ird 1
for {set iatm 1} {$iatm<=$natm} {incr iatm} {

if {$ird==1} {
for {set i 1} {$i<=12} {incr i} {set cn($i) 0}
gets $rdpdbcon line
scan [string range $line 6 84] "%d %d %d %d %d %d %d %d %d %d %d %d %d" self cn(1) cn(2) cn(3) cn(4) cn(5) cn(6) cn(7) cn(8) cn(9) cn(10) cn(11) cn(12)

set tmplist {}
#Formation of connectivity list
for {set i 1} {$i<=12} {incr i} {
  if {$cn($i)==0} {break}
  lappend tmplist [expr $cn($i)-1]
}
}

if {$self==$iatm} {
#puts Atom\ serial:\ $iatm\ \ User-connectivity:\ $tmplist
set sel [atomselect top "serial $iatm"]
if {$itype==0} {
  $sel setbonds "{$tmplist}"
} else {
  set orglist [$sel getbonds]
  $sel setbonds "{[concat [lindex $orglist 0] $tmplist]}"
}
$sel delete
set ird 1
} else {
set ird 0
}

}
close $rdpdbcon

作者
Author: sobereva 时间: 2020-11-27 18:12
把问题搞得太复杂了
直接在VMD里用mouse - add/remove bonds，点击相应两个原子就能连上键，然后保存成mol文件，以后再载入mol时成键关系还是你想要的状态。仔细看
谈谈VMD可视化程序的连接关系的判断和设置问题
http://sobereva.com/534（http://bbs.keinsci.com/thread-16396-1-1.html）

作者
Author: liaorongbao 时间: 2020-11-27 18:34

sobereva 发表于 2020-11-27 18:12
把问题搞得太复杂了
直接在VMD里用mouse - add/remove bonds，点击相应两个原子就能连上键，然后保存成mol ...

博文http://sobereva.com/534中那个简单的处理个别个连线的方法我也读到了。实际上我面对的体系有几十个原子需要在VMD中人工连线。2H2O.pdb这个案例只是为了简化问题。我不想按http://sobereva.com/534费力地在VMD中手动连接几十个原子。现在就是想学学批量处理的http://sobereva.com/121。求社长再回复一下具体怎么操作

作者
Author: sobereva 时间: 2020-11-27 18:45
1 nskipline 8

2 tcl语言本来就那么规定的。用/或\\也可以

3 读取有问题，导致self没被赋值

作者
Author: liaorongbao 时间: 2020-11-27 18:53

sobereva 发表于 2020-11-27 18:12
把问题搞得太复杂了
直接在VMD里用mouse - add/remove bonds，点击相应两个原子就能连上键，然后保存成mol ...

把上边的7改成8似乎就出现了，谢谢社长关注！

欢迎光临计算化学公社 (http://bbs.keinsci.com/)