计算化学公社
标题: 使用Multiwfn做表面活性剂aNCI分析遇到的几个小问题 [打印本页]
作者Author: 牧生 时间: 2021-9-14 17:37
标题: 使用Multiwfn做表面活性剂aNCI分析遇到的几个小问题
本帖最后由 牧生 于 2021-9-15 14:34 编辑
在Centos下进行的操作，使用amber力场，使用在线服务器得到表面活性剂的GAFF力场的itp。
一、向盒子中加入5个十二烷基苯磺酸钠阴离子SDBS，加满水，能量最小以后，跑npt。再冻结其中一个SDBS （resid 4），去掉控压，跑nvt，总共得到1001帧。
二、使用VMD观察轨迹运动时先载入nvt.gro，删除仅有的一帧，再载入nvt.xtc，看到resid 4没有动，其余的水分子和SDBS都在动。
三、待轨迹运动完了以后，将进度条拖回到第0帧，使用语句resid 4 or same resid as within 5 of resid 4，选择resid 4周围5 A的分子。图像如下：
(, 下载次数 Times of downloads: 20)
第一个需要请教的：在这里，有怪异的线条，应该是水分子从盒子另一头穿过来了，这个会影响结果吗？
然后save ，选择resid 4 or same resid as within 5 of resid 4，存为aNCI.xyz
四、运行gmx editconf–f nvt.tpr –o nvt.pdb，然后vmd打开nvt.pdb，同样将resid 4 or same resid as within 5 of resid 4保存为aNCI.pdb。
(, 下载次数 Times of downloads: 19)
这个pdb看起来就很正常。 打开aNCI.pdb看到resid 4的原子编号为52-102
五、使用命令Multiwfn aNCI.xyz，依次选择20 3 1,1001 11 52-102 2.5 A 0.1 （希望得到的图形光滑，所以选择0.1）， 计算结束以后选择6，在当前目录输出得到avgRDG.cub和avgsl2r.cub。再Multiwfn avgRDG.cub，选择13, 13, 1.5 100， 2 选择52-102，-2，重新输出avgRDG.cub。
六、将作图脚本avgRDG.vmd拷到当前目录，然后打开vmd以后，source avgRDG.vmd，得到图像如下 （我发现这里无论在当前目录，还是VMD安装目录下操作，没有差别）
(, 下载次数 Times of downloads: 18)
需要请教的第二个问题：此前两个图中，水分子都非常靠近SDBS，但是这个图中，水分子就跑开了，这样对吗？
第三个问题，按数字0查看此时等值面内的分子，resid和residue 就都变成1了。。这个编号的变化，是正常的吗？
(, 下载次数 Times of downloads: 24)
第四个问题，图中左边蓝色圆圈1中的浅蓝色，可以归结为氧原子与周围水的作用，但我直觉上觉得这个蓝色有点浅，。。上边圆圈2内的红色面可以归结为分子弯曲时，原子间的斥力，但是也同时出现了明显的蓝色，那么，如何解释这个蓝色呢？
我知道其中一定有不正确的操作，但是又不是太清楚错在哪
2021.9.15 编辑，重新使用了整个轨迹去计算，得到结果如下
(, 下载次数 Times of downloads: 18)
现在看起来就比较正常了，端基磺酸根上带负电氧原子与周围的原子有强烈的作用，颜色有蓝色和红色。
但是，如果针对大体系，如何才能保存部分位置，以及特定分子及其周围的水分子
作者Author: sobereva 时间: 2021-9-15 04:31
1 没有影响。实际上你可以把体系用多帧叠加方式显示，只要SDBS附近一圈都包满了水即可
2 播放一下aNCI.xyz轨迹，看看一开始那些水是否后来已经分散开，不再紧密包裹SDBS，如果是这样，aNCI结果没法用。这种情况必须考虑所有的水，和下文的情况一样，这样在轨迹的xyz文件里SDBS附近第一配位层才始终有水
使用Multiwfn研究分子动力学中的弱相互作用
http://sobereva.com/186
3 不用管这个
4 先把前面的2解决了
作者Author: 牧生 时间: 2021-9-15 08:57
现在情况就是，播放aNCI.xyz轨迹，水分子就会跑散开。
我在另存xyz轨迹之前已经选中了update selection every frame，播放时，SDBS周围的水分子并没有明显跑散开。但是，另存为xyz以后，播放该xyz轨迹时，水分子就散开了。
作者Author: sobereva 时间: 2021-9-15 09:38
所以当前轨迹根本不能用
update selection every frame对于保存轨迹不生效
作者Author: 牧生 时间: 2021-9-15 09:46
本帖最后由 牧生 于 2021-9-15 16:58 编辑
那么，有没有办法，在只冻住SDBS的情况下，保存SDBS及附近水分子为xyz格式的时候，保持SDBS附近一直都有水分子。
有个异想天开的办法，比如一个很大的盒子（原子数有数十万个），足够长时间的npt平衡以后，SDBS聚集形态已经稳定，只有水分子还在到处跑。。此时选择某个SDBS附近15A分子，另存为pdb格式的小盒子，冻结这个SDBS后，再对这个盒子进行nvt（原子数有几千个左右，虽然体系还是有点大，但是服务器也能算的动），结束后，对这个小盒子整体的轨迹进行分析。这样，就能实现考察动态过程中的表面活性剂-水间的相互作用。。
这个，从原则上讲，是否可行？？
作者Author: sobereva 时间: 2021-9-16 06:39
没有直接办法实现。
但可以自己写脚本，循环每一帧，每一帧都用within语句选择SDBS附近一定距离的水，导出成结构文件。但同样的within距离内每一帧的水的数目可能不同，而xyz文件要求必须每一帧原子数相同，因此对每一帧还需要做个循环，逐渐改变within的阈值，当选中的水分子数和第一帧时候选中的水分子数相同时结束循环。
你说的那个异想天开的做法可以。但是你没法恰好保证小盒子里的密度合适，所以取出这么一个子体系后，应当再做一下NPT，等盒子尺寸自发稳定之后再NVT
作者Author: 牧生 时间: 2021-10-3 14:53
本帖最后由 牧生 于 2021-10-3 16:00 编辑
我已经试过了，从大盒子里，取出一个小盒子，再进行NPT去调节密度的过程中，分子都会移动。即使冻结了需要分析的分子，其余的分子都会移动，造成目标分子附近的环境发生变化。如果直接NVT，分子照样会到处移动，液体形状也变化，有可能变成球形。所以，我觉得从大盒子里挖出小盒子，再想办法aNCI分析，此法行不通。再者，我现在觉得，即使强行用大体系去计算其中一个长条形分子的等值面，得到一个毛毛虫形状的等值面图，可能意义也不大。我们更多的是关注亲水头基与附近分子的作用，其余部分，有一些就行了，不必过多去计较。
照着博文http://bbs.keinsci.com/thread-9472-1-7.html 的方法，将我在1楼中，第四步中取出的pdb进行IGM分析
- 四、运行gmx editconf–f nvt.tpr –o nvt.pdb，然后vmd打开nvt.pdb，同样将resid 4 or same resid as within 5 of resid 4保存为aNCI.pdb。
得到的图形如下，两个SDBS分子仍然靠的近，水分子围绕在SDBS附近，没有散开。
(, 下载次数 Times of downloads: 15)
这个看起来就比较正常了，SDBS亲水头基附近出现了明显的等值面，颜色为蓝色的部分，解释为头基的氧原子与水作用，且作用较强。SDBS碳链附近有部分离散的绿色等值面，可解释为碳链与水的作用，或者SDBS之间的弱相互作用，且该种作用较弱。SDBS附近没有毛毛虫形状的等值面，这个应该是比较简洁好看，且好分析解释。
此外，在http://bbs.keinsci.com/thread-21826-1-1.html中，sob也提到了
- 笔者之前在《通过独立梯度模型(IGM)考察分子间弱相互作用》（http://sobereva.com/407）里还介绍过IGM分析，比起NCI分析的好处是可以自定义片段，只让片段间的相互作用显示出来，就免去了做格点数据屏蔽的步骤，而且哪怕用比aNCI明显更大的格点间距，等值面边缘的锯齿也很不明显。笔者在aNCI的思想上提出了aIGM，即平均的IGM，专用于分析动态环境中的弱相互作用，此分析将另文说明。
作者Author: sobereva 时间: 2021-10-3 23:01
IGM一个缺点是等值面太厚，有时候影响考察
aIGM看Multiwfn手册3.23.9节。用法很简单，选主功能20的子功能12照提示操作即可
作者Author: 牧生 时间: 2021-10-7 12:51
本帖最后由 牧生 于 2021-10-7 15:37 编辑
已经试过了使用aIGM方法，就以教程中苯酚-水体系，只用了前10帧，大致看了下效果，那个等值面太厚了，像个馒头，能不能解决一下啊。
(, 下载次数 Times of downloads: 15)
如果使用aIGM来分析SDBS与水作用，得到的图如下，效果似乎还没aNCI的那种薄层好看(为了省时间，只用了10帧。如果用1000帧，可能这个等值面会全部包裹住SDBS的头基，导致什么都看不出来。)
(, 下载次数 Times of downloads: 12)
这个等值面太厚，即使在isosurface里面去改，仍然是很厚，等值面甚至比原子还大，有时候还会覆盖原子本身。
其次，有个提醒新手的地方，如果使用aIGM，要把IGM_inter.vmd脚本中的cub改成当前实际得到的cub文件名。
作者Author: sobereva 时间: 2021-10-7 23:18
IGM就是如此，本来就是不如NCI。这是为什么我后来特意提出了IGMH解决这个问题，不过IGMH依赖于波函数，不支持轨迹平均
| 欢迎光临 计算化学公社 (http://bbs.keinsci.com/) |
Powered by Discuz! X3.3 |