仓库源文站点原文


layout: post title: Pi-Pi堆积距离和堆积角度的计算 categories:


Pi-Pi堆积的概念和意义

Pi-Pi堆积的表征量

表征Pi-Pi堆积的物理量一般使用两个, 距离和角度. 距离为某一组原子质心或中心到另一组原子所在平面的距离, 角度为两组原子所在平面之间的夹角. 计算这两个量就需要知道如何计算一组原子所在平面的方程. 对于完全刚性的平面组原子, 三个非共线原子就可以确定其所在的平面方程(平面的三点式方程). 对于近似刚性的原子组, 则需要通过多元线性拟合来确定平面方程.

有了平面方程 $ax+by+cz+d=0, \vec n=(a,b,c)$, 可计算空间一点 $\vec r(x,y,z)$ 到平面的距离

$$d={\vec r \cdot \vec n +d\over \abs {\vec n} }$$

也可计算此点到平面的垂足

$$\vec v=\vec r - {\vec n \cdot \vec r+d \over \abs{\vec n}^2} \vec n$$

Pi-Pi堆积量的计算

gmx analyze可以进行多元线性拟合得到平面方程. 通过脚本调用即可, 虽然麻烦但不困难, 就不再示例了.

这里我们用vmd加上tcl数学库来进行计算. vmd并没有自带这个库, 须自行安装. 我们只需要其中的线性代数库linalg.tcl, 其中包含了多元线性拟合的奇异值分解算法leastSquaresSVD, 使用倒也方便. 这里多说一句, vmd加上tcl的数学库后基本可以进行各式分析了, 虽然速度未必佳.

既然是决定用vmd进行分析了, 那顺便也学习一下vmd的绘图方法吧. 将计算结果实时显示出来, 也更容易确定计算结果是否正确.

代码

gmxtool

使用方法

得到了轨迹以后, 先对轨迹进行处理: 分子完整化, , 然后就可以计算Pi-Pi堆积量了.

  1. 运行MD模拟, 获得轨迹
  2. 对轨迹进行处理: 先gmx trjconv -pbc whole完整化分子, 再gmx trjconv -center -fit rot+trans对分子进行居中叠合
  3. vmd载入轨迹: vmd conf.gro traj.xtc
  4. 修改pistack.tcl中的两个原子组
  5. vmd命令窗口执行source pistack.tcl
  6. 分析得到的输出文件pistack.xvg, 也可以播放轨迹, 查看每一帧计算结果

示例

两个苯分子的模拟

输出文件pistack.xvg中每行数据会依次列出

存在问题

  1. 轨迹过大时速度不佳
  2. 考虑到周期性边界条件, 所得距离和角度可能存在不连续性
  3. vmd每帧绘制时draw material Transparent似乎不起作用? 必须存在错误语句才可以?
  4. tcl是totally command language, 也是toy command language, 太(Tai)差(Cha)了(Le)