计算化学公社

标题: 巨大体系的范德华表面静电势图的快速绘制方法 [打印本页]

作者
Author:
sobereva    时间: 2019-5-7 16:11
标题: 巨大体系的范德华表面静电势图的快速绘制方法
重要提示下文是很早以前写的,如今已经没意义了!如今使用Multiwfn官网上的最新版本,结合下文提到的xtb产生的molden文件作为输入文件,直接用《使用Multiwfn+VMD快速地绘制静电势着色的分子范德华表面图和分子间穿透图》(http://sobereva.com/443)介绍的ESPiso脚本绘制才是最好的做法。因为后来得益于《Multiwfn使用的高效的静电势算法的介绍文章已于PCCP期刊发表!》(http://sobereva.com/614)中介绍的Multiwfn计算静电势速度的巨大提升,以及《巨幅降低Multiwfn结合VMD绘制分子表面静电势图耗时的一个关键技巧》(http://sobereva.com/602)所描述的巨幅降低耗时的技巧,用443号文中的做法比下文的做法不仅更省事、省时,还不需要调用收费的Gaussian里的cubegen。

提示:如果你要绘制1000原子以上的特大体系的表面静电势图,应当使用《基于原子电荷极快速绘制超大体系的分子表面静电势图》(http://sobereva.com/639)里的做法,用Multiwfn基于原子坐标和原子电荷来计算,耗时超级低!

巨大体系的范德华表面静电势图的快速绘制方法
Rapidly plotting electrostatic potential colored van der Waals surface for huge system

文/Sobereva@北京科音
First release: 2019-May-7  Last update: 2020-Apr-26


有人问我,怎么绘制他的含有多达336个原子的分子的分子表面静电势图。对这么大体系,绘制这种图看似极其困难,很难吃得消,但其实只要方法得当,哪怕只有台普通四核机子,也完全没有压力,而且绘制的图像效果极佳,本文就来介绍一下流程。绘制用的步骤和《使用Multiwfn+VMD快速地绘制静电势着色的分子范德华表面图和分子间穿透图》(http://sobereva.com/443)里介绍的基本一样,但牵扯到一些其它过程和技巧,没读过此文的读者一定要先阅读一下。本文用到的Multiwfn可以在其官网http://sobereva.com/multiwfn上免费下载,读者用的Multiwfn版本应是2019年4月或以后更新的。

本文笔者的条件是:CPU为Intel i7-2630QM(4核),操作系统是Win7-64bit,Multiwfn是2019-May-5版,xtb是2019-Apr-16版,Gaussian是G16 A.03版,VMD是1.9.3版,这几个程序里只有Gaussian是收费的。笔者跑xtb时是在VMware装的CentOS 7.4虚拟机下跑的,因为此程序目前只有Linux版,不会装Linux的话看《在VMware 15中安装CentOS 7.6的完整过程视频演示》(http://sobereva.com/454)。

要绘制分子表面静电势图,不管通过什么途径,肯定得先产生波函数,至少得算一次单点任务。但如果你只有一台四核机子,遇见300多个原子的体系,用Gaussian算DFT就甭想了,哪怕在ORCA里用纯泛函开RI,照样很难算得动。虽然用半经验方法如PM6算几百个原子很快,但是一般的程序都不支持基于半经验方法得到的波函数去算静电势。碰到这种情况,最佳的出路就是用Grimme的xtb程序。此程序的安装和基本用法在《将Gaussian与Grimme的xtb程序联用搜索过渡态、产生IRC、做振动分析》(http://sobereva.com/421)已经介绍了,不会用xtb的人一定要看一下。

先把那个336原子的大分子的结构文件载入Multiwfn(支持mol、mol2、pdb、xyz、gjf等众多格式),进入主功能100的子功能2,选择相应选项导出为xyz文件,比如叫336.xyz。然后把此文件随便放到一个空目录下,在此目录下运行此命令通过xtb做单点计算:xtb 336.xyz --molden(如果这个结构是未经优化过的,可再加上--opt让xtb优化此结构)。在笔者的4核机子上只花了20秒钟就算完了,在当前目录下出现了molden.input,这是Molden输入文件,可以被Multiwfn读取并做波函数分析。

把Multiwfn目录下的settings.ini里的nthreads设为机子里实际的CPU物理核心数,使得CPU运算能力能在接下来的计算中充分发挥。

下面产生此体系的电子密度格点数据。启动Multiwfn,载入molden.input,然后依次输入以下命令:
5  //计算格点数据
1  //计算电子密度
3  //高质量格点
计算过程花费仅不到1分钟,然后选2导出格点数据,此时当前目录下出现了density.cub,将之改名为density1.cub备用。

然后我们要计算静电势格点数据,用Multiwfn调用Gaussian目录下的cubegen速度会很快,但是这要求必须以fch作为输入文件,所以我们要把当前的波函数导出为fch文件。接着在Multiwfn里输入
0  //返回主菜单
100  //其它功能part 1
2  //导出文件
7  //导出fch文件
336.fch  //输出文件名
马上当前目录下就有了336.fch。现在退出Multiwfn。

我们修改Multiwfn目录下的settings.ini,将其中的cubegenpath路径设为你机子里实际的Gaussian目录下的cubegen可执行文件的路径,比如D:\study\G16W\cubegen.exe。关于这点,在http://sobereva.com/435里有更多说明。

当前这个体系很大,当Multiwfn调用cubegen去处理这个文件时,会由于G16的cubegen默认可用内存上限不够而报错。为解决这个问题,对于Windows用户,进入操作系统的控制面板,选择“系统”-“高级系统设置”-“高级”-“环境变量”,然后在用户变量里点击“新建”,变量名填GAUSS_MEMDEF,变量值填1400MB,然后点“确定”。之后cubegen处理fch文件时最大可用内存就为1400MB了,对当前体系完全够了(如果你的Gaussian是64bit版本,GAUSS_MEMDEF也可以设更大,比如可以填4GB。而对于32bit版Gaussian,GAUSS_MEMDEF通常最多只能设1400MB,再大往往就无法运行了)。如果你用的是Linux版Gaussian,应当运行比如export GAUSS_MEMDEF=4GB命令来让cubegen可以最大调用4GB内存(建议把这行命令也加入到~/.bashrc文件中,否则每次开启终端后都得再运行一次此命令)。

启动Multiwfn,依次输入
336.fch
5  //计算格点数据
12  //计算静电势
3  //高质量格点
现在会看到Multiwfn正在调用cubegen计算静电势。这一步是整个过程最耗时的一步,花了50分钟(而在笔者的2*E5-2696v3 36核服务器上三分多钟就能算完)。然后选2导出格点数据,此时当前目录下出现了totesp.cub,将之改名为ESP1.cub。

现在把density1.cub和ESP1.cub都拷到VMD目录下,把Multiwfn目录下的examples\drawESP目录下的ESPiso.vmd也拷到VMD目录下。启动VMD,在文本窗口输入source ESPiso.vmd并回车,马上就可以看到下图:
(, 下载次数 Times of downloads: 283)

这图虽然效果已经很好了,但看着还有点乱,因此我们进入Graphics - Materials,选里面的EdgyGlass材质(这是用于显示当前分子表面用的材质),然后适当调节各个参数令效果尽可能好,最后改成下面这样:
(, 下载次数 Times of downloads: 260)

然后进入Graphics - Representation,点击对应显示分子结构的rep,把CPK显示方式切换为Licorice,并让键的粗细从默认的0.3改为0.2。然后点击对应Isosurface显示方式的那个rep,选择Trajectory标签页,把Color Scale Data Range范围设成比默认-0.03~0.03更宽的-0.04~0.04使得颜色变化更柔和些。最后图像效果如下,两个视角都给了,可见效果相当令人满意!(可能有人觉得这个图还没一开始看到的酷炫,但此时的图更容易分辨分子表面上不同区域静电势的差异)
(, 下载次数 Times of downloads: 252)

图中红色和蓝色分别对应于静电势为正和为负的区域。如果你想把色彩刻度轴加上去也很容易,做法在http://sobereva.com/443文中已经介绍了。

本文的做法绘制甚至四五百个原子的体系的静电势图也是完全没问题的,哪怕只有一台四核机子,也完全可以算得动,只不过是时间长一些而已。本文给出的静电势分布其实算不上很准确,毕竟xtb程序做的GFN-xTB计算只相当于半经验层面的DFT方法,用的是极小基,不过对于几百个原子的静电势图的绘制,谁也不会那么关注定量大小,只要定性正确就足够了,当前绘制的图完全可以满足这个需要。本文的例子也充分体现出,对于好几个百原子的精度要求不很高的波函数分析,xtb+Multiwfn可谓是黄金组合。

下面是xtb产生的波函数与精度理想的B3LYP/def2-TZVP波函数绘制的分子表面静电势图的对比,可见差异不是很大,基于xtb的结果已经完全可以满足需要了。
(, 下载次数 Times of downloads: 247)

不过,如果你想基于xtb的波函数按照《使用Multiwfn的定量分子表面分析功能预测反应位点、分析分子间相互作用》(http://sobereva.com/159)去定量考察分子表面静电势特征,那么xtb的波函数质量还是差点意思。比如对于上面这个多巴胺的分子表面静电势最大、最小点数值,基于xtb算的结果与基于B3LYP/def2-TZVP波函数算的结果甚至能差出20 kcal/mol的程度。而对于计算拟合静电势电荷,比如CHELPG,基于xtb波函数算的结果与基于B3LYP/def2-TZVP算的结果差异最大的可达0.23,这也不是个小数目。所以,xtb波函数对应的静电势可以满足肉眼观看、定性分析需要,但定量层面上还明显达不到准确的标准。
作者
Author:
exity    时间: 2019-5-7 16:20
...",花了50分钟"...
作者
Author:
ShuangfeiZhu    时间: 2019-5-7 19:37
点赞!我以前是用orca开RI算一个单点分析过上百个原子,,贼慢,这个方法快的太多了
作者
Author:
Novice    时间: 2019-5-7 21:58
其实我已经用过这个方法了,哈哈,就是没有对比过xTB的结果和Gaussian算的差异是否显著,不知道社长知道不
作者
Author:
sobereva    时间: 2019-5-8 09:35
exity 发表于 2019-5-7 16:20
...",花了50分钟"...

若有个像样的服务器,整个绘制过程从头到尾不会超过10分钟。文中刻意用4核机子举例,主要目的是说明,哪怕是不是专门做计算的组,没有服务器,用个普通PC机按此文做法也能算得动。
作者
Author:
sobereva    时间: 2019-5-8 10:00
Novice 发表于 2019-5-7 21:58
其实我已经用过这个方法了,哈哈,就是没有对比过xTB的结果和Gaussian算的差异是否显著,不知道社长知道不

文末加入了对比图
作者
Author:
haowei    时间: 2019-5-10 21:25
请问老师,本文所使用的xtb需要与gaussian联用吗
作者
Author:
sobereva    时间: 2019-5-10 23:09
haowei 发表于 2019-5-10 21:25
请问老师,本文所使用的xtb需要与gaussian联用吗


作者
Author:
hzamber    时间: 2019-5-12 18:41
xtb有个--esp calculate electrostatic potential on VdW-grid选项,要是xtb可以计算cube格式电势静数据,再结合xtb生成的电子密度cube数据,是不是可以绘制更大分子的表面静电势图?



作者
Author:
sobereva    时间: 2019-5-13 00:17
hzamber 发表于 2019-5-12 18:41
xtb有个--esp calculate electrostatic potential on VdW-grid选项,要是xtb可以计算cube格式电势静数据, ...

这个选项只是生成分子表面上的点的静电势,没法实现本文的效果
作者
Author:
sobereva    时间: 2019-5-13 00:29
在本文末尾添加了一段,强调了xtb给出的波函数用于定量分析还精度不够
作者
Author:
Frank    时间: 2020-6-14 13:44
请教社长个问题,我想做一个带有5个embedding charge下-2价体系的静电势,从https://xtb-docs.readthedocs.io/en/latest/pcem.html这个网页可以看到怎么考虑embedding,但是不太清楚怎么写xtb的输入命令,例如xtb water_4.coord -I pcem.input这个命令应该怎么在哪个位置添加--charge --UHF以及--molden? 另外考虑embedding charge用的是*.coord和常规xtb计算需要的*.xyz文件也不同,应该怎么处理?
作者
Author:
mol    时间: 2020-6-27 16:30
sob老师您好,请教大概3500原子的聚集体系能算动不。。。
作者
Author:
sobereva    时间: 2020-6-27 17:55
mol 发表于 2020-6-27 16:30
sob老师您好,请教大概3500原子的聚集体系能算动不。。。


如果愣算的话,可能需要给很大内存cubegen才能工作
如果你非要算的话,建议先用比如1700原子左右的情况试试水
如果你的体系是一堆分子构成的,且对静电势要求不是特别高,我建议先用Multiwfn对每个单分子产生拟合静电势电荷得到单分子的chg文件,然后拼成整个体系的chg文件,让Multiwfn载入后计算基于原子电荷算静电势,精度基本够了(载入chg文件,然后选择计算的函数那一步选择基于原子电荷的静电势),这样的话几千原子的静电势计算也不是难事。

作者
Author:
mol    时间: 2020-6-28 21:28
sobereva 发表于 2020-6-27 17:55

如果愣算的话,可能需要给很大内存cubegen才能工作
如果你非要算的话,建议先用比如1700原子左右的 ...

多谢sob老师,用了第二种方法,电荷导入参考了http://sobereva.com/365
作者
Author:
beidiwu    时间: 2021-4-9 08:48
社长您好,我在尝试用该方法计算离子态有机物,得出的结果似乎丢失了电荷信息。按步骤生成xyz文件,fch文件似乎均按照charge 0, spin 2再进行下一步。请问是哪里出了问题?
作者
Author:
sobereva    时间: 2021-4-12 09:06
beidiwu 发表于 2021-4-9 08:48
社长您好,我在尝试用该方法计算离子态有机物,得出的结果似乎丢失了电荷信息。按步骤生成xyz文件,fch文件 ...

看不懂你的意思
确保xtb计算时候明确设定了实际的电荷和自旋多重度,并且确保之后Multiwfn载入xtb产生的molden.input文件后屏幕上正确提示了此体系的电荷和自旋多重度
作者
Author:
JZS12345    时间: 2021-4-24 09:46
sob老师,您好!请问一下按这种方法(主功能5)计算出来的静电势和主功能12计算出来的表面静电势的异同点在哪里。PS:我现在想从静电势分析两个物质结合的难易程度
作者
Author:
sobereva    时间: 2021-4-24 11:52
JZS12345 发表于 2021-4-24 09:46
sob老师,您好!请问一下按这种方法(主功能5)计算出来的静电势和主功能12计算出来的表面静电势的异同点在 ...

主功能5本身是计算静电势格点数据的,即计算一个矩形区域内均匀分布的各个点的静电势。而主功能12才是计算分子表面上各个点的静电势。目的和用处完全不同。
作者
Author:
JZS12345    时间: 2021-4-24 20:03
sobereva 发表于 2021-4-24 11:52
主功能5本身是计算静电势格点数据的,即计算一个矩形区域内均匀分布的各个点的静电势。而主功能12才是计 ...

好的,谢谢sob老师
作者
Author:
mol    时间: 2022-4-11 10:20
sobereva 发表于 2020-6-27 17:55

如果愣算的话,可能需要给很大内存cubegen才能工作
如果你非要算的话,建议先用比如1700原子左右的 ...

sob老师好,最近想用chg文件结合功能5计算分子表面各个静电势区域的面积分布,但是体系较大,2000多原子,出现如下错误,请问有什么好的建议么,多谢。
Spatial range of grid data:
X is from   50.1911  to  125.0999 Bohr
Y is from   39.0417  to  148.3246 Bohr
Z is from    4.2330  to   96.7351 Bohr
The number of point in x,y,z:   301   438   371  Total:  48911898

Calculating grid data of electron density...
Note: Promolecular approximation is used
Note: All exponential functions exp(x) with x< -40.000 will be ignored
Progress: [##################################################]  100.00 %     /
Calculation of grid data took up wall clock time        37 s

The number of boundary cubes:            926527
The number of corners of boundary cubes:   1849560
Generating isosurface by Marching Tetrahedra algorithm, please wait...
The number of surface vertices (V):    2745878
The number of surface edges (E):       8237514
The number of triangular facets (F):   5491676
V+F-E-2=        38
Warning: V+F-E-2=0 is violated! Probably grid spacing is too large or the isosurface is not closed
Volume enclosed by the isosurface:170683.80325 Bohr^3  (********** Angstrom^3)
Among all surface vertices:
Min-X:     27.2773  Max-X:   65.7392 Angstrom
Min-Y:     21.4109  Max-Y:   77.8648 Angstrom
Min-Z:      2.8069  Max-Z:   50.7219 Angstrom

Eliminating redundant surface vertices...
After elimination, V=   830465,  E=  2491276,  F=  1660850,  V+F-E-2=    37
Isosurface area: 38573.07375 Bohr^2  (********** Angstrom^2)
Sphericity:  0.3858

Calculating electrostatic potential at surface vertices, please wait patiently

Initializing LIBRETA library (fast version) for ESP evaluation ...

Intel MKL ERROR: Parameter 8 was incorrect on entry to DGEMM .
forrtl: severe (41): insufficient virtual memory
作者
Author:
sobereva    时间: 2022-4-11 18:46
mol 发表于 2022-4-11 10:20
sob老师好,最近想用chg文件结合功能5计算分子表面各个静电势区域的面积分布,但是体系较大,2000多原子 ...

说明具体用的输入文件,以及怎么操作的
作者
Author:
mol    时间: 2022-4-12 16:44
sobereva 发表于 2022-4-11 18:46
说明具体用的输入文件,以及怎么操作的

sob老师好,输入文件是chg文件,2786个原子,包括四列,原子名,x,y,z,charge。
muliwfn3.8中依次输入:
12
0
产生以上错误信息
作者
Author:
sobereva    时间: 2022-4-12 16:58
mol 发表于 2022-4-12 16:44
sob老师好,输入文件是chg文件,2786个原子,包括四列,原子名,x,y,z,charge。
muliwfn3.8中依次输入: ...

chg文件不包含波函数信息,所以12-0之后没法正常产生分子表面
作者
Author:
mol    时间: 2022-4-12 17:00
sobereva 发表于 2022-4-12 16:58
chg文件不包含波函数信息,所以12-0之后没法正常产生分子表面

酱紫,那像这种情况有比较好的处理方法吗
作者
Author:
sobereva    时间: 2022-4-13 04:45
mol 发表于 2022-4-12 17:00
酱紫,那像这种情况有比较好的处理方法吗

看此文
基于原子电荷极快速绘制超大体系的分子表面静电势图
http://sobereva.com/639http://bbs.keinsci.com/thread-28913-1-1.html
作者
Author:
mol    时间: 2022-4-13 08:38
sobereva 发表于 2022-4-13 04:45
看此文
基于原子电荷极快速绘制超大体系的分子表面静电势图
http://sobereva.com/639(http://bbs.kein ...

非常感谢sob老师!
作者
Author:
wilely    时间: 2025-1-21 15:55
用xtb里面的“esp”命令直接算可以吗?
作者
Author:
sobereva    时间: 2025-1-22 06:22
wilely 发表于 2025-1-21 15:55
用xtb里面的“esp”命令直接算可以吗?

--esp选项是将范德华表面顶点的坐标和静电势输出为xtb_esp.dat和xtb_esp.cosmo文件,和本文介绍的绘制表面静电势图的全流程没有可比性
作者
Author:
wilely    时间: 2025-1-22 09:52
sobereva 发表于 2025-1-22 06:22
--esp选项是将范德华表面顶点的坐标和静电势输出为xtb_esp.dat和xtb_esp.cosmo文件,和本文介绍的绘制表 ...

我体系600多原子,用该方法算的好慢。
作者
Author:
sobereva    时间: 2025-1-22 13:03
wilely 发表于 2025-1-22 09:52
我体系600多原子,用该方法算的好慢。

用更好的CPU所有核并行、linux版Multiwfn、确保计算方式得当
作者
Author:
xxx127    时间: 2025-3-10 15:13
sob老师好,我想请教一下,按照这种方法算出来的精度是不是比按照这篇文章《基于原子电荷极快速绘制超大体系的分子表面静电势图》(http://sobereva.com/639)算出来的要高很多呢
作者
Author:
sobereva    时间: 2025-3-10 22:49
xxx127 发表于 2025-3-10 15:13
sob老师好,我想请教一下,按照这种方法算出来的精度是不是比按照这篇文章《基于原子电荷极快速绘制超大体 ...

更高
但还比不上基于DFT得到的静电势
作者
Author:
xxx127    时间: 2025-3-14 15:04
sobereva 发表于 2025-3-10 22:49
更高
但还比不上基于DFT得到的静电势

非常感谢老师的回复!
作者
Author:
xu123123    时间: 2025-6-2 21:41
老师好,请问我在分析体系的静电势的时候,寻找表面的最大和最小值,发现只能看见最小值,无法看见最大值,最大值显示的是*****这一串字符。这个是超出了程序每列显示的最大位数吗,请问可以修改最大位数吗?
作者
Author:
sobereva    时间: 2025-6-2 22:18
xu123123 发表于 2025-6-2 21:41
老师好,请问我在分析体系的静电势的时候,寻找表面的最大和最小值,发现只能看见最小值,无法看见最大值, ...

输入文件用得不对,认真看当前帖子以及下文
认真看
详谈Multiwfn支持的输入文件类型、产生方法以及相互转换
http://sobereva.com/379http://bbs.keinsci.com/thread-6020-1-1.html


作者
Author:
xu123123    时间: 2025-6-3 15:51
sobereva 发表于 2025-6-2 22:18
输入文件用得不对,认真看当前帖子以及下文
认真看
详谈Multiwfn支持的输入文件类型、产生方法以及相互 ...

谢谢老师回答,我用的输入文件是.fch文件,这个文件是由.molden文件在multiwfn(3.8_dev)软件中生成出来的.fch文件,老师您的意思是把.molden手动把原子序数那一列改成实际核电荷数,这样才能正常分析体系表面的极值吗?
作者
Author:
sobereva    时间: 2025-6-3 22:39
xu123123 发表于 2025-6-3 15:51
谢谢老师回答,我用的输入文件是.fch文件,这个文件是由.molden文件在multiwfn(3.8_dev)软件中生成出来的 ...

molden文件能提供的信息和fch是一样的,直接把molden文件载入Multiwfn算静电势就完了,干嘛转成fch,毫无意义,完全多此一举
并且说清楚molden文件是什么体系、用什么程序以什么输入文件计算产生的。我根本没说那个意思
作者
Author:
xu123123    时间: 2025-6-4 00:00
好的老师,molden文件是金团簇的中性体系,molden是通过orca产生的,那我再用molden文件算一算静电势。
作者
Author:
xu123123    时间: 2025-6-4 09:52
sobereva 发表于 2025-6-3 22:39
molden文件能提供的信息和fch是一样的,直接把molden文件载入Multiwfn算静电势就完了,干嘛转成fch,毫无 ...

老师,我把ocra生成的molden文件直接放入multiwfn里计算后,发现还是无法找到表面静电势的最大值,请问老师我哪里个地方没有做对。
作者
Author:
sobereva    时间: 2025-6-4 09:56
xu123123 发表于 2025-6-4 09:52
老师,我把ocra生成的molden文件直接放入multiwfn里计算后,发现还是无法找到表面静电势的最大值,请问老 ...

molden文件压缩后传网盘贴链接,说明Multiwfn里输入的所有命令
那叫ORCA

倘若用了赝势,需要恰当修改molden文件,下文明确说了,否则静电势算得没意义
详谈Multiwfn支持的输入文件类型、产生方法以及相互转换
http://sobereva.com/379http://bbs.keinsci.com/thread-6020-1-1.html
懒得改就用ORCA和Multiwfn的最新版本,不需要自己修改,上文说了


作者
Author:
xu123123    时间: 2025-6-4 15:05
sobereva 发表于 2025-6-4 09:56
molden文件压缩后传网盘贴链接,说明Multiwfn里输入的所有命令
那叫ORCA

老师这个是体系molden的网盘链接:https://pan.baidu.com/s/1_-Kc6Fwr01fV9aWYE9MUjg?pwd=zysm

在multiwfn的指令是
12
0
作者
Author:
sobereva    时间: 2025-6-4 17:07
xu123123 发表于 2025-6-4 15:05
老师这个是体系molden的网盘链接:https://pan.baidu.com/s/1_-Kc6Fwr01fV9aWYE9MUjg?pwd=zysm

在mul ...

按我前面说的改完了molden文件再计算
Multiwfn都提示得如此明确了:
(, 下载次数 Times of downloads: 116)

作者
Author:
xu123123    时间: 2025-6-4 22:16
好的谢谢老师
作者
Author:
xu123123    时间: 2025-6-5 19:22
sobereva 发表于 2025-6-4 17:07
按我前面说的改完了molden文件再计算
Multiwfn都提示得如此明确了:

我按老师说的用orca6.0和multiwfn的两个软件的最新版已经成功将静电势计算的问题解决了,谢谢老师的指点。

老师我还有一个疑问是我分别用orca5.0和orca6.0生成的同一个结构的molden文件,我发现两个molden的Au的原子索引都是79,C和H分别是6和1 都是一样的,两个文件最前面就只有xyz坐标不一样,请问老师这个是orca内部计算的原因吗?
作者
Author:
sobereva    时间: 2025-6-6 00:58
xu123123 发表于 2025-6-5 19:22
我按老师说的用orca6.0和multiwfn的两个软件的最新版已经成功将静电势计算的问题解决了,谢谢老师的指点 ...

自己根据计算流程搞清楚xyz坐标为什么不同,没细节别人没法判断
作者
Author:
xu123123    时间: 2025-6-6 10:53
本帖最后由 xu123123 于 2025-6-6 10:55 编辑

这个图是新旧OCRA两个版本分别生成的同一结构的molden文件。

旧版本的molden导入multiwfn会提醒一个警告”当您忘记将.molden文件[atoms]字段中的原子指数修改为实际核电荷时,可能会使用ECP”,并且测试发现旧版本无法正确计算静电势的值,只有ORCA6.0新版本能正确计算。

然后我打开了新旧版本生成的molden文件,发现旧版本和新版本的OCRA的原子指数是一样的,也就是核电荷数是一样的,

我看了老师的另一个帖子,是在旧版本生成的molden文件的做些修改即可(即更改原子指数)也能正常放入multiwfn里正常计算静电势值,但是两个版本OCRA生成的molden文件里面Au C H原子的原子指数都是一样的,一个新版本生成的molden能正常计算静电势,另一个旧版本不能正常计算静电势,二者只有xyz坐标不一样,为什么旧版本OCRA5.0不能正常计算,而新版本OCRA6.0能计算静电势,二者仅仅是坐标不一样吗?
作者
Author:
Uus/pMeC6H4-/キ    时间: 2025-6-6 11:13
xu123123 发表于 2025-6-6 10:53
这个图是新旧OCRA两个版本分别生成的同一结构的molden文件。

旧版本的molden导入multiwfn会提醒一个警告 ...

ORCA 6.0与5.0在输出的molden文件上有个[Pseudo]字段的重要区别,看http://bbs.keinsci.com/thread-48197-1-1.html这帖。可能这个字段在molden文件末尾,所以只看开头就会错过之。
至于坐标的差异,如果确信输入文件的结构与关键词均完全一致,可以把关键词贴出来看看有没有哪个会涉及坐标变换。
作者
Author:
xu123123    时间: 2025-6-6 11:47
Uus/pMeC6H4-/キ 发表于 2025-6-6 11:13
ORCA 6.0与5.0在输出的molden文件上有个字段的重要区别,看http://bbs.keinsci.com/thread-48197-1-1.htm ...

谢谢老师解答,我再看看帖子
作者
Author:
sobereva    时间: 2025-6-7 01:31
xu123123 发表于 2025-6-6 10:53
这个图是新旧OCRA两个版本分别生成的同一结构的molden文件。

旧版本的molden导入multiwfn会提醒一个警告 ...

别乱拼ORCA

不同大版本默认设置不同,导致同一个输入文件计算结果有差异非常正常,下文也说了,你不给具体输入文件没法判断
谈谈不同量子化学程序计算结果的差异问题
http://sobereva.com/573http://bbs.keinsci.com/thread-20122-1-1.html

如今没有任何理由还用ORCA 5


作者
Author:
xu123123    时间: 2025-6-18 14:54
老师好,我想请问这个静电势柱状图中的纵轴 area 的意思是,area 的值越大,这一点所占表面积就越大吗?

这个图是不考虑这个点的值是正值还是负值,只考虑这点的表面积的大小吗?
作者
Author:
sobereva    时间: 2025-6-19 05:54
xu123123 发表于 2025-6-18 14:54
老师好,我想请问这个静电势柱状图中的纵轴 area 的意思是,area 的值越大,这一点所占表面积就越大吗?

...

体现的是不同静电势范围对应的分子表面上的面积
你说的“点”意义不明
作者
Author:
xu123123    时间: 2025-6-19 15:05
sobereva 发表于 2025-6-19 05:54
体现的是不同静电势范围对应的分子表面上的面积
你说的“点”意义不明

就是一个静电势的极值点,

这个图画的是所有的静电势的点,还是所有的静电势极值点呢?

然后就是面积越大,说明这个对应静电势能量的极值点的数量就越多是吗?不知道我理解的对不对 老师。
作者
Author:
sobereva    时间: 2025-6-19 17:22
xu123123 发表于 2025-6-19 15:05
就是一个静电势的极值点,

这个图画的是所有的静电势的点,还是所有的静电势极值点呢?

这图跟极值点毫无关系
极值点是表面上一个具体的点,而当前的图统计的是范德华表面上不同静电势数值范围包含的面积
作者
Author:
xu123123    时间: 2025-6-19 19:30
sobereva 发表于 2025-6-19 17:22
这图跟极值点毫无关系
极值点是表面上一个具体的点,而当前的图统计的是范德华表面上不同静电势数值范围 ...

谢谢老师的解答。




欢迎光临 计算化学公社 (http://bbs.keinsci.com/) Powered by Discuz! X3.3