计算化学公社

 找回密码 Forget password
 注册 Register
Views: 7906|回复 Reply: 6

[Multiwfn资源与经验] 使用Multiwfn计算分子的长宽高以及显示分子的主轴

[复制链接 Copy URL]

4万

帖子

99

威望

4万

eV
积分
89946

管理员

公社社长+计算化学玩家

发表于 Post on 2018-6-10 11:05:31 | 显示全部楼层 Show all |阅读模式 Reading model
使用Multiwfn计算分子的长宽高以及显示分子的主轴

文/Sobereva @北京科音
First release: 2018-Jun-10  Last update: 2021-Oct-23

0 前言

经常有人问怎么计算分子的长宽高,之前笔者在《谈谈分子半径的计算和分子形状的描述》(http://sobereva.com/190)一文中也曾讨论分子的长宽高怎么算,但那个博文用的蒽的例子过于简单,而实际分子形状往往比较复杂,而且结构文件里的分子的朝向可能歪斜着,没法直接用那篇博文里的做法。本文介绍通过Multiwfn程序非常方便、准确地计算长宽高的做法。读者请务必使用最新版本Multiwfn,可以在官网http://sobereva.com/multiwfn下载。如果读者对Multiwfn不了解,看《Multiwfn入门tips》(http://sobereva.com/167)和《Multiwfn FAQ》(http://sobereva.com/452)。没看过《谈谈分子半径的计算和分子形状的描述》的读者应先看一下了解相关基本知识。


1 计算分子长宽高的原理

计算长宽高的第一步,是先把分子“摆正”,让分子最长的轴冲着某个笛卡尔轴。手动在可视化程序里旋转不准确,而且比较麻烦。一个自动化的做法是计算分子的惯性矩张量,如下所示:
1.png
式中x,y,z是原子坐标,m是原子质量。这个3*3矩阵的三个本征矢对应分子的三个主轴,彼此正交,相应的本征值是绕着主轴的转动惯量。最小的本征值对应的主轴是分子的最长轴,最大的本征值对应的是分子的最短轴。因此,如果我们旋转分子,让分子的三个主轴朝向正好对应三个笛卡尔轴,分子就被摆正了。代码实现其实很容易,就是利用惯性矩张量的本征矢矩阵来对各个原子的坐标做个线性变换而已。

摆正之后,寻找x坐标最小的原子,减去其范德华半径,记为xmin;寻找x坐标最大的原子,加上其范德华半径,记为xmax,则x方向尺寸就是xmax-xmin。类似地,对y、z方向也都这么计算。这样分子的长、宽、高就都有了。

当然,能用长宽高来描述一个分子,前提是这个分子必须能够通过矩形来比较好地描述,诸如一个分子是月亮形,或者比如碳纳米管侧面接出来一个长链,或者是一个四面体形的团簇,此时长宽高这种描述就没什么意义,且有误导性。而且,注意很多大分子有很多不同构象,不同构象的长宽高算出来可能差异很大,这种情况也不适合用长宽高这种概念说事。另外,还有的分子用球形、柱形、椭球形等形式描述更恰当,此时也不要局限于用长宽高来描述。


2 计算长宽高实例

例如当前有一个pdb文件,用可视化程序打开后,看到其朝向和位置如下图所示。我们想计算其长宽高
2.png
.mol、.xyz、.pdb、.fch等含有结构信息的Multiwfn可以识别的格式都可以直接载入到Multiwfn里面计算长宽高。什么格式含有什么信息看《详谈Multiwfn支持的输入文件类型、产生方法以及相互转换》(http://sobereva.com/379)。启动Multiwfn后,载入此pdb文件,然后依次输入
100  //主功能100
21   //输出各种描述体系结构的信息
size  //计算分子尺寸
此时程序自动把体系平移和旋转,各个主轴都分别冲着各个笛卡尔轴了。屏幕上显示以下信息
Farthest distance:   21(H )  ---   45(H ):    10.250 Angstrom
vdW radius of   21(H ): 1.200 Angstrom
vdW radius of   45(H ): 1.200 Angstrom
Diameter of the system:    12.650 Angstrom
Radius of the system:     6.325 Angstrom
Length of the three sides:    12.650     7.513     7.619 Angstrom


从中可以看到,体系中相距最远的两个原子是H21和H45,二者的Bondi范德华半径也输出了,将俩原子的范德华半径加到它们之间的距离上,就是屏幕上输出的体系的直径,屏幕上输出的半径是其一半。可见利用Multiwfn计算分子半径/直径很方便。Length of the three sides后面的三个值就是体系的三个边长,可以被视为长宽高。

如果想把长宽高图形化地展现一下,选择屏幕上的选项1 Visualize the new orientation and molecular box,会看到下图,图中蓝色方框的三个边长正对应于分子的长宽高。下图把Ratio of atomic size设为了4.0,此时图上原子球的半径正对应于Bondi范德华半径,可见蓝色方框正好卡着范德华表面。

3.png

如果想把坐标已经被Multiwfn旋转后的当前体系导出成pdb文件,从而能在VMD等程序中显示,可以选择2 Export the geometry in new orientation as new.pdb in current folder,此时当前目录下就会产生new.pdb。此文件中有个CRYST1字段,记录了此体系的长宽高。如果用VMD程序(http://www.ks.uiuc.edu/Research/vmd/)载入此pdb文件,就会读入盒子信息,并可以绘制出来。

把new.pdb拖到VMD里载入,然后在命令行窗口输入pbc box,盒子就被显示出来了。将背景改为白色,选Display-Orthographic用正交视角,在Graphics-Representation里把体系设置以VDW风格显示,然后把长宽高手动ps到图上去,就看到了下图的效果。
4.png

顺带一提,利用VMD的绘图命令,还可以把主轴绘制出来。把new.pdb载入到VMD后,把以下内容复制到VMD的命令行窗口运行(最后要按一次回车确保最后一条命令运行了)
  1. set lena [molinfo top get a]
  2. set lenb [molinfo top get b]
  3. set lenc [molinfo top get c]
  4. set cenx [expr $lena/2]
  5. set ceny [expr $lenb/2]
  6. set cenz [expr $lenc/2]
  7. draw color red
  8. draw cylinder "$cenx $ceny 0" "$cenx $ceny $lenc" radius 0.15 resolution 20
  9. draw color green
  10. draw cylinder "0 $ceny $cenz" "$lena $ceny $cenz" radius 0.15 resolution 20
  11. draw color blue
  12. draw cylinder "$cenx 0 $cenz" "$cenx $lenb $cenz" radius 0.15 resolution 20
复制代码

调节显示方式为CPK,然后就会看到下面的图像,红绿蓝三个圆柱对应于三个主轴
5.png

下面是两篇已发表的文章中使用Multiwfn的这个功能的实例
6.png

下图中对分子不同片段表面分别着色突出显示的做法见《使用Multiwfn和VMD计算分子表面积和片段表面积》(http://sobereva.com/487)。

7.png

评分 Rate

参与人数
Participants 2
eV +10 收起 理由
Reason
kulaomega + 5 好物!
captain + 5 好物!

查看全部评分 View all ratings

北京科音自然科学研究中心http://www.keinsci.com)致力于计算化学的发展和传播,长期开办最高质量的各种计算化学类培训:初级量子化学培训班基础(中级)量子化学培训班分子动力学与GROMACS培训班量子化学波函数分析与Multiwfn程序培训班,内容介绍以及往届资料购买请点击链接查看。这些培训是计算化学从零快速入门以及进一步全面系统性提升研究水平的最佳途径。培训各种常见问题见《北京科音办的培训班FAQ》
欢迎加入“北京科音”微信公众号获取培训最新消息、避免错过网上最有价值的计算化学文章!
欢迎加入人气最高、水准最高的综合性理论与计算化学交流QQ群“思想家公社QQ群”:1号:18616395,2号:466017436,3号:764390338,搜索群号能搜到哪个说明目前哪个能加,合计9000人,讨论范畴相同
思想家公社的门口Blog:http://sobereva.com(发布大量原创计算化学相关博文)
Multiwfn主页:http://sobereva.com/multiwfn(最强大的量子化学波函数分析程序)
ResearchGate:https://www.researchgate.net/profile/Tian_Lu
Money and papers are rubbish, get a real life!

239

帖子

0

威望

2264

eV
积分
2503

Level 5 (御坂)

发表于 Post on 2018-6-10 16:06:49 | 显示全部楼层 Show all
上次我做材料的室友让我帮他算一下分子长度,我实在是为难了

275

帖子

0

威望

4494

eV
积分
4769

Level 6 (一方通行)

发表于 Post on 2019-4-8 10:57:47 | 显示全部楼层 Show all
sob老师您好,可否把本文用到的例子的pdb文件也上传一下,我想测试一下^_^

1020

帖子

0

威望

3510

eV
积分
4530

Level 6 (一方通行)

发表于 Post on 2019-4-8 20:08:11 来自手机 | 显示全部楼层 Show all
个人认为,如果计算纯几何的主轴,应该把质量权重去掉,得到纯几何的坐标张量

57

帖子

0

威望

957

eV
积分
1014

Level 4 (黑子)

发表于 Post on 2020-10-4 14:15:24 | 显示全部楼层 Show all
做了一点小改动,可以在VMD显示输入以下代码,这样显示出来的红绿蓝,分别对应,new.pdb中的x,y, z, 当三者数值相差不大时,这样比较便于指认,不容易弄错
  1. set lena [molinfo top get a]
  2. set lenb [molinfo top get b]
  3. set lenc [molinfo top get c]
  4. set cenx [expr $lena/2]
  5. set ceny [expr $lenb/2]
  6. set cenz [expr $lenc/2]
  7. draw color red
  8. draw cylinder "0 $ceny $cenz" "$lena $ceny $cenz" radius 0.15 resolution 20
  9. draw color green
  10. draw cylinder "$cenx 0 $cenz" "$cenx $lenb $cenz" radius 0.15 resolution 20
  11. draw color blue
  12. draw cylinder "$cenx $ceny 0" "$cenx $ceny $lenc" radius 0.15 resolution 20
复制代码


4万

帖子

99

威望

4万

eV
积分
89946

管理员

公社社长+计算化学玩家

 楼主 Author| 发表于 Post on 2021-10-24 02:46:10 | 显示全部楼层 Show all
给本文末尾又加入了一个文献中使用此功能的实例
北京科音自然科学研究中心http://www.keinsci.com)致力于计算化学的发展和传播,长期开办最高质量的各种计算化学类培训:初级量子化学培训班基础(中级)量子化学培训班分子动力学与GROMACS培训班量子化学波函数分析与Multiwfn程序培训班,内容介绍以及往届资料购买请点击链接查看。这些培训是计算化学从零快速入门以及进一步全面系统性提升研究水平的最佳途径。培训各种常见问题见《北京科音办的培训班FAQ》
欢迎加入“北京科音”微信公众号获取培训最新消息、避免错过网上最有价值的计算化学文章!
欢迎加入人气最高、水准最高的综合性理论与计算化学交流QQ群“思想家公社QQ群”:1号:18616395,2号:466017436,3号:764390338,搜索群号能搜到哪个说明目前哪个能加,合计9000人,讨论范畴相同
思想家公社的门口Blog:http://sobereva.com(发布大量原创计算化学相关博文)
Multiwfn主页:http://sobereva.com/multiwfn(最强大的量子化学波函数分析程序)
ResearchGate:https://www.researchgate.net/profile/Tian_Lu
Money and papers are rubbish, get a real life!

5

帖子

0

威望

161

eV
积分
166

Level 3 能力者

发表于 Post on 2022-10-9 13:15:48 | 显示全部楼层 Show all
请问一下,利用惯性矩张量的本征矢矩阵来对各个原子的坐标做个线性变换,是用旋转矩阵的逆矩阵处理原子坐标么

本版积分规则 Credits rule

手机版 Mobile version|北京科音自然科学研究中心 Beijing Kein Research Center for Natural Sciences|京公网安备 11010502035419号|计算化学公社 — 北京科音旗下高水平计算化学交流论坛 ( 京ICP备14038949号-1 )|网站地图

GMT+8, 2023-2-6 03:51 , Processed in 0.194002 second(s), 25 queries .

快速回复 返回顶部 返回列表 Return to list