计算化学公社
标题: 使用Multiwfn超级方便地计算出概念密度泛函理论中定义的各种量 [打印本页]
作者Author: sobereva 时间: 2019-5-21 12:15
标题: 使用Multiwfn超级方便地计算出概念密度泛函理论中定义的各种量
注：轨道权重福井函数/双描述符特别适合前线轨道简并或近简并体系，也是使用本文提到的功能来计算，但本文没涉及，在另一篇博文里专门做了介绍《通过轨道权重福井函数和轨道权重双描述符预测亲核和亲电反应位点》（http://sobereva.com/533），非常建议看过本文后阅读。
注：本文对应2024-Apr-17及以后更新的Multiwfn及手册，不要用老版本！
使用Multiwfn超级方便地计算出概念密度泛函理论中定义的各种量
Using Multiwfn to easily calculate various quantities defined in the concept density functional theory
文/Sobereva@北京科音 First release: 2019-May-21 Last update: 2024-Apr-17
1 前言
由Parr最初开始发展的概念密度泛函理论（conceptual density functional theory，CDFT）也称密度泛函反应性理论（DFRT），是量子化学和波函数分析领域中的一个重要组成，专门用来预测和解释化学物质的反应活性、反应位点等问题，对于研究化学反应有重要意义。相关介绍以及笔者写过的相关博文见《概念密度泛函综述和重要文献合集》（http://bbs.keinsci.com/thread-384-1-1.html）。在“量子化学波函数分析与Multiwfn程序培训班”（http://www.keinsci.com/workshop/WFN_content.html）里专门有一节由笔者全面讲授反应位点的预测与反应活性分析，对本文涉及的各种量做了非常详细具体的介绍并给了诸多应用例子，推荐感兴趣者参加。
CDFT框架里有大量概念，包含很多实空间函数、局部指数和全局指数，它们经常出现在各种文献里用来讨论反应问题。只要有N、N+1、N-1电子态的波函数和能量，就可以计算它们中的绝大部分。N指的是体系的电子数，比如某个分子最稳定状态是中性的，则N就是指中性时的电子数，N+1就是指带一个负电荷的阴离子状态的电子数，N-1就是指带一个正电荷的阳离子态的电子数。其实计算CDFT里的那些量并不复杂，手动算起来没什么难度，然而如果把所有量全都一一算出来，还是比较费事，而且容易弄错数据。另外，网上也经常有人问我怎么去算那些量，由于初学者的量子化学基础知识往往很薄弱，经常解释起来很费劲，甚至怎么说也说不明白。于是，笔者决定在Multiwfn里直接加入一个功能，通过极简单的步骤，就可以一次性把所有CDFT框架内常见的量全都算出来，过程完全傻瓜化。笔者相信这个功能有极高的实用性，经常要计算CDFT里涉及的量的读者务必要重视此功能。
这个功能是Multiwfn的主功能22。Multiwfn最新版可以在官网http://sobereva.com/multiwfn免费下载。相关入门知识看《Multiwfn入门tips》（http://sobereva.com/167）、《Multiwfn FAQ》（http://sobereva.com/452）。
使用本文介绍的Multiwfn的功能发表文章时，除了要按照Multiwfn程序启动时的提示引用Multiwfn的原文外，也请同时引用下面的书籍章节，其中专门介绍了Multiwfn程序中的概念密度泛函理论分析的功能和实现：
Tian Lu, Qinxue Chen, Realization of Conceptual Density Functional Theory and Information-Theoretic Approach in Multiwfn Program. In Conceptual Density Functional Theory, WILEY-VCH GmbH: Weinheim (2022); pp 631-647. DOI: 10.1002/9783527829941.ch31
2 能计算的量
这个功能能直接计算的量主要包含下面这些，分为三大类，具体定义在这里就不写了，因为在Multiwfn手册3.25节都已经明确给出了。
• 全局指数（即对整个体系而言的数值）
垂直电离能（Vertical ionization potential，VIP）
垂直电子亲和能（Vertical electron affinity，VEA）
Mulliken电负性（Mulliken electronegativity）
化学势（Chemical potential）
电子的硬度（Hardness）
电子的软度（Softness）
亲电指数ω（Electrophilicity index）
亲电指数ωcubic（在ω之后提出，形式更严格）
亲电描述符ε（Electrophilic descriptor）
亲核指数（Nucleophilicity index）
• 实空间函数（即变量是三维空间坐标的函数）
福井函数（Fukui function）
双描述符（Dual descriptor）
局部软度（Local softness）
局部超软度（Local hyper-softness）
局部亲电指数（Local electrophilicity index）
局部亲核指数（Local nucleophilicity index）
• 原子指数（即每个原子各有一个数值）
简缩福井函数（Condensed Fukui function）
简缩双描述符（Condensed dual descriptor）
简缩局部软度（Condensed local softness）简缩局部超软度（Condensed local hyper-softness）
相对亲电指数（Relative electrophilicity index）
相对亲核指数（Relative nucleophilicity index）
简缩局部亲电指数（Condensed local electrophilicity index）
简缩局部亲核指数（Condensed local nucleophilicity index）
简缩局部ωcubic亲电指数
亲核和亲电超离域度（superdelocalizabilities）
上述原子指数都可以写为基于原子电荷计算的形式。根据一些文章的测试，Hirshfeld电荷非常适合这种目的（不建议用NPA等其它电荷），所以Multiwfn自动计算它们都是用的Hirshfeld电荷，相关讨论见《TCA上的一篇对比不同原子电荷预测反应位点、亲电/亲核性的文章》（http://bbs.keinsci.com/thread-15512-1-1.html）。
除上述外，Multiwfn还可以计算在反应位点预测方面比双描述符更严格的双描述符势（dual descriptor potential），笔者专门有文章介绍：《使用Multiwfn计算双描述符势预测反应位点》（http://sobereva.com/708）。
3 实例：对苯酚计算CDFT定义的各种量
这个功能的用法和细节在手册3.25节都介绍了，本文就不再累述，下面直接看一个具体使用例子，对苯酚计算各种量。看完例子后建议再看一下手册这一节。
3.1 准备.wfn文件
用本文介绍的功能计算各种量之前需要先在当前目录下给出N.wfn、N+1.wfn和N-1.wfn，分别记录当前体系三个电子态的波函数信息。可以自行提供它们，做法见《详谈Multiwfn支持的输入文件类型、产生方法以及相互转换》（http://sobereva.com/379）；但如果你是Gaussian用户，按照下面的流程操作省事得多。
计算前述所有量要用的几何结构都是相同的，都是N电子态的极小点结构。对于苯酚这个例子而言，我们需要对苯酚中性状态进行几何优化，这一步是用户自己手动做的，用什么程序、什么级别无所谓，靠谱就行。普通有机体系用常用的B3LYP-D3(BJ)/6-31G*级别一般就够了，如果想更便宜一些，用半经验、GFN-xTB之类也可以。在B3LYP/6-31G*下优化中性苯酚得到的结构已经直接给了，是Multiwfn目录下的examples目录中的phenol.xyz。
启动Multiwfn，依次输入
examples\phenol.xyz //一开始载入的文件只要含有中性苯酚优化后的结构即可，其它格式如pdb/fch/gjf/mol/mol2等也都可以
22 //计算CDFT里定义的各种量
1 //产生波函数文件
[直接按回车] //这一步让你输入gjf文件里对应的计算单点任务的关键词。按回车代表用默认的B3LYP/6-31G*，此级别虽然很便宜，但对于考察CDFT定义的量来说一般基本已可以满足需要
[直接按回车] //这一步让你输入计算N、N+1和N-1态用的电荷和自旋多重度。按回车代表用对三个态分别用默认的0 1、-1 2和1 2
现在当前目录下已经有了N.gjf、N+1.gjf和N-1.gjf，是用于产生.wfn文件的单点任务的Gaussian输入文件。如果你不自己改gjf里的.wfn文件路径的话，把这些gjf都用Linux下的Gaussian执行后，N.wfn、N+1.wfn和N-1.wfn会产生在当前目录下，而如果用Windows下的Gaussian执行则会产生在临时文件目录下（比如D:\study\G16W\Scratch下）。如果你当前机子里已经装了Gaussian，而且之前Multiwfn的settings.ini里的gaupath已经被你设为了实际的Gaussian可执行文件路径，那么建议此时直接在Multiwfn窗口里输入y，这样Multiwfn就会自动调用Gaussian计算这个三个gjf文件，并在当前目录下产生相应的三个wfn文件，然后gjf和out文件会自动被删掉。显然，如果你当前的机子里没Gaussian，应当把gjf放到有Gaussian的机子上算，算完之后再把得到的.wfn文件拷回来。
注：让Multiwfn顺利调用Windows版Gaussian需要做额外的设置，见Multiwfn手册Appendix 1。
计算这些gjf文件过程中，特别是计算N+1.gjf（通常对应阴离子）时，可能出现SCF不收敛问题，届时按照此文加上帮助SCF收敛的关键词即可：《解决SCF不收敛问题的方法》（http://sobereva.com/61）。
三个态的波函数文件也不是必须叫做N.wfn、N+1.wfn和N-1.wfn并处在当前目录下。用选项2或3开始计算时，若Multiwfn发现当前目录下缺少这些文件的一个或多个，程序就会让用户直接输入相应电子态的.wfn或.wfx或.fch或.mwfn文件的路径。因此，比如你想用自行算出来的N.fch、N+1.fch、N-1.fch做后面的计算也完全可以。
另外，如果你没有Gaussian而只有ORCA，也可以选择-2 Choose the quantum chemistry program used in option 1之后再选2。之后再选选项1产生输入文件时，Multiwfn产生的就是三个态的ORCA的输入文件N.inp、N+1.inp和N-1.inp，手动算完了之后就会在当前目录下产生N.wfn、N+1.wfn和N-1.wfn（其它产生的文件可删掉）。如果你把settings.ini里的orcapath设为了实际ORCA可执行文件路径，那么Multiwfn在产生输入文件后会问你是否自动调用ORCA去计算它们并在计算结束后删掉其它文件。
3.2 计算各种指数
现在当前目录下已经有了N.wfn、N+1.wfn和N-1.wfn，可以开始算了。在当前模块里选择2 Calculate various quantitative indices，然后程序就会依次载入wfn文件，读取其中的能量信息和波函数信息，自动计算Hirshfeld电荷，最后给出各种指数。对于6-31G*基组下的苯酚，在Intel 4核机子下不出10秒钟就都算完了。算完的数据被输出到了当前目录下的CDFT.txt中，此文件内容如下：
Note: the E(HOMO) of TCE used for evaluating nucleophilicity index is the value evaluated at B3LYP/6-31G* level
Hirshfeld charges, condensed Fukui functions and condensed dual descriptors
Units used below are "e" (elementary charge)
Atom q(N) q(N+1) q(N-1) f- f+ f0 CDD
1(C ) -0.0587 -0.1185 0.0852 0.1439 0.0598 0.1018 -0.0841
2(C ) -0.0390 -0.1674 0.0268 0.0658 0.1284 0.0971 0.0626
3(C ) -0.0597 -0.1873 0.0319 0.0916 0.1276 0.1096 0.0360
4(C ) 0.0737 0.0216 0.1739 0.1001 0.0522 0.0762 -0.0480
5(C ) -0.0731 -0.1956 0.0090 0.0821 0.1225 0.1023 0.0404
6(C ) -0.0415 -0.1730 0.0336 0.0751 0.1315 0.1033 0.0563
7(H ) 0.0417 -0.0048 0.0993 0.0576 0.0464 0.0520 -0.0112
8(H ) 0.0450 -0.0193 0.0904 0.0455 0.0642 0.0548 0.0187
9(H ) 0.0473 -0.0160 0.0953 0.0480 0.0633 0.0557 0.0153
10(H ) 0.0387 -0.0230 0.0854 0.0467 0.0617 0.0542 0.0150
11(H ) 0.0444 -0.0207 0.0914 0.0470 0.0651 0.0561 0.0181
12(O ) -0.1977 -0.2443 -0.0517 0.1460 0.0465 0.0963 -0.0995
13(H ) 0.1789 0.1482 0.2294 0.0505 0.0307 0.0406 -0.0197
Condensed local electrophilicity/nucleophilicity index (e*eV)
Atom Electrophilicity Nucleophilicity
1(C ) 0.02576 0.45535
2(C ) 0.05533 0.20827
3(C ) 0.05496 0.28982
4(C ) 0.02249 0.31688
5(C ) 0.05280 0.25977
6(C ) 0.05665 0.23773
7(H ) 0.02001 0.18234
8(H ) 0.02767 0.14394
9(H ) 0.02729 0.15195
10(H ) 0.02659 0.14779
11(H ) 0.02807 0.14875
12(O ) 0.02005 0.46203
13(H ) 0.01325 0.15965
Condensed local softness (e/Hartree), relative electrophilicity/nucleophilicity (dimensionless) and condensed local hyper-softness (e/Hartree^2)
Atom s- s+ s0 s+/s- s-/s+ s(2)
1(C ) 0.3761 0.1562 0.2661 0.4154 2.4075 -0.5746
2(C ) 0.1720 0.3355 0.2538 1.9501 0.5128 0.4271
3(C ) 0.2392 0.3333 0.2863 1.3933 0.7177 0.2459
4(C ) 0.2617 0.1363 0.1990 0.5210 1.9195 -0.3276
5(C ) 0.2145 0.3202 0.2674 1.4927 0.6699 0.2762
6(C ) 0.1964 0.3436 0.2700 1.7498 0.5715 0.3847
7(H ) 0.1506 0.1213 0.1360 0.8058 1.2411 -0.0764
8(H ) 0.1189 0.1678 0.1433 1.4116 0.7084 0.1279
9(H ) 0.1255 0.1655 0.1455 1.3188 0.7583 0.1045
10(H ) 0.1221 0.1612 0.1416 1.3210 0.7570 0.1024
11(H ) 0.1228 0.1703 0.1465 1.3860 0.7215 0.1239
12(O ) 0.3816 0.1216 0.2516 0.3186 3.1387 -0.6794
13(H ) 0.1319 0.0804 0.1061 0.6095 1.6408 -0.1346
E(N): -307.464860 Hartree
E(N+1): -307.383614 Hartree
E(N-1): -307.163438 Hartree
E_HOMO(N): -0.218913 Hartree, -5.9569 eV
E_HOMO(N+1): 0.161297 Hartree, 4.3891 eV
E_HOMO(N-1): -0.464864 Hartree, -12.6496 eV
Vertical IP: 0.301421 Hartree, 8.2021 eV
Vertical EA: -0.081246 Hartree, -2.2108 eV
Mulliken electronegativity: 0.110088 Hartree, 2.9956 eV
Chemical potential: -0.110088 Hartree, -2.9956 eV
Hardness (=fundamental gap): 0.382667 Hartree, 10.4129 eV
Softness: 2.613235 Hartree^-1, 0.0960 eV^-1
Electrophilicity index: 0.015835 Hartree, 0.4309 eV
Nucleophilicity index: 0.116285 Hartree, 3.1643 eV
上面的数据基本没有什么需要特别解释的，信息都非常易懂，也可以去对照手册3.25节的式子。唯一值得特别一提的是亲核指数（以及局部亲核指数），这个量是依赖于TCE（四氰基乙烯）的HOMO能量定义的，TCE的HOMO能量按理说是应当使用和当前计算级别完全一样的级别来计算的，这样得到的全局和局部亲核指数才是严格的，但Multiwfn在给出的时候是自动用笔者事先算好的B3LYP/6-31G*下的TCE的HOMO值-0.335198 Hartree算的。如果你当前用的不是B3LYP/6-31G*级别且追求更严格的结果，应当自行用当前级别去算TCE的HOMO值然后手动根据定义得到全局和局部亲核指数。
在Multiwfn手册里也有手动计算简缩福井函数和简缩双描述符的例子，见4.7.3节，用的也是Hirshfeld电荷和B3LYP/6-31G*级别，因此结果和上面自动算出来的完全一样。相比之下，用本文介绍的功能实在方便太多了，所有重要的量一口气全都输出了！
关于ωcubic和亲电描述符（ε）的计算
前面例子中输出信息里的亲电指数是Parr早年提出的形式（ω），也是目前用的最多的。后来在J. Phys. Chem. A, 124, 2090 (2020)又有人提出了新的而且更严格的形式，称为ωcubic。文中发现对于R-X...NH3型卤键二聚体（R为不同基团），卤原子上的简缩局部形式的ωcubic值与结合能有很不错的线性关系，R^2达到0.94，比ω的相关性明显更好，因此ωcubic很有实际意义，可以用于预测弱相互作用强度、解释弱相互作用内在特征。亲电描述符ε于Int. J. Quantum Chem., 124, e27366 (2024)中提出，文中通过35个有机分子测试发现它与Mayr亲电性参数的相关性远好于ω。
由于ωcubic和ε都涉及到比ω的更高阶项，即超硬度（hyperhardness），因此计算时还额外依赖于N-2态的能量，所以默认情况下Multiwfn不计算它们。如果想计算它们的话，对于本例应如下操作。载入输入文件后，依次输入
22 //计算CDFT里定义的各种量
-1 //要求计算ωcubic和ε
1 //产生波函数文件
[直接按回车] //用默认的B3LYP/6-31G*
[直接按回车] //对N、N+1、N-1、N-2态用默认的电荷和自旋多重度（对于中性闭壳层分子这是适合的）
此时当前目录下就有了N.gjf、N+1.gjf、N-1.gjf、N-2.gjf，你可以让Multiwfn直接调用Gaussian计算也可以自行计算。算完后就有了这四个态的.wfn文件。之后进入选项2，在输出的CDFT.txt中可见不仅我们之前看到的那些量都给出来了，还给出了体系整体的ωcubic、各个原子的简缩ωcubic，以及ε。与此同时，计算它们过程中用到的第二垂直电离能也会顺带输出出来。
3.3 计算福井函数和双描述符
接下来计算福井函数和双描述符。在当前模块里选择3 Calculate grid data of Fukui function, dual descriptor and related functions，然后选择一个合适的格点设定，对于当前这样小体系选择Medium quality grid就够了（大体系应当用High quality grid或其它选项，详见http://sobereva.com/452中的Q39的讨论），然后程序会依次载入当前目录下的N.wfn、N+1.wfn和N-1.wfn并计算电子密度格点数据，之后会看到一个菜单。通过选择相应选项，可以把各种类型的福井函数以及双描述符直接显示成等值面图。比如此例我们依次选择选项1、2、3、4来把f+、f-、f0福井函数和双描述符都依次绘制出来，下图是把图像手动合并到一起的图，等值面数值都是用的0.007：
(, 下载次数 Times of downloads: 416)
通过主功能5手动计算福井函数和双描述符的方法在Multiwfn手册4.5.4节已经详细介绍了，得到的图像和上图完全一样。使用本文介绍的步骤，明显远比手动操作方便得多得多得多！
值得一提的是，如果你想要更好的显示质量，可以借助VMD。比如我们选择6 Export grid data of f- as f-.cub in current folder，然后把当前目录下产生的f-.cub按照《在VMD里将cube文件瞬间绘制成效果极佳的等值面图的方法》（http://sobereva.com/483）里介绍的方法载入VMD并通过Tachyon渲染器渲染，只需要很简单几步就可以得到下图，可见效果极好！
(, 下载次数 Times of downloads: 399)
当前功能也可以对局部软度、局部亲电/亲核指数观看等值面和导出cube文件，因为它们的定义都是对特定类型的福井函数乘上一个数值，而Multiwfn当前功能里可以设乘上的系数。比如局部软度s-定义为全局软度乘以福井函数f-，我们根据前面算出的数据已知苯酚的局部软度是2.613235 Hartree^-1，因此若我们想导出局部软度cube文件，就选-1 Set the scale factor to various grid data，输入2.613235，之后选6 Export grid data of scaled f- as f-.cub in current folder导出的cube文件就对应于s-的格点数据了，若选2 Visualize isosurface of scaled f-则看到的是s-的等值面。
另外，在J. Math. Chem., 62, 461 (2024)中还介绍了超软度的概念，对应于双描述符乘上全局软度的平方，使得双描述符有尺寸一致性。显然在当前界面里把要乘的数值设为全局软度的平方，再选8 Export grid data of scaled dual descriptor as DD.cub in current folder，产生的cube文件就对应于超软度了。
4 总结
本文介绍了Multiwfn的超级便利的一口气算出来概念密度泛函理论里涉及到的各种量的功能，这使得初学者也可以非常顺利、快速且不犯错误地计算这些量。相信此功能对于将概念密度泛函理论广泛地应用于实际化学问题的研究有积极的促进作用，也鼓励大家在日常研究中充分使用此功能，以从繁复的手动操作中解脱出来。
经常有用户问本文介绍的方法是否只能用于闭壳层中性体系、用于其它体系时怎么对N、N+1和N-1状态在Multiwfn里输入净电荷和自旋多重度。在这里我就再次着重强调一下：显然本文的方法可以用于非中性、非闭壳层体系。例如，若当前体系是个中性自由基（并假定是普通有机类体系，没什么复杂的情况），N状态的净电荷和自旋多重度是0 2，N+1时是-1 1，N-1时是1 1。原因很容易理解：N+1和N-1的时候分别比N状态多一个和少一个电子，净电荷自然分别是N状态的减1和加1。N状态原本有一个未配对电子，多一个电子时必然令其配对，少一个电子时它必然丢失，体系自然就成了闭壳层状态，故N+1和N-1时的自旋多重度为1。注意这里假设N+1和N-1状态时不涉及到分子轨道简并因而变成高自旋的可能，一般有机体系都如此，而碰上过渡金属配合物体系时，拿不准的话应当对N+1和N-1态分别计算最低几种可能的自旋态的能量看哪种最低以判断它们应当是什么自旋多重度。
类似地，对于普通有机类体系，其它几种情况的净电荷和自旋多重度为：
• 若当前体系是个+1电荷的阳离子自由基，N状态是1 2，N+1状态是0 1，N-1状态是2 1。
• 若当前体系是个-1电荷的阴离子自由基，N状态是-1 2，N+1状态是-2 1，N-1状态是0 1。
• 若当前体系是个+1电荷的闭壳层阳离子，N状态是1 1，N+1状态是0 2，N-1状态是2 2。
• 若当前体系是个-1电荷的闭壳层阴离子，N状态是-1 1，N+1状态是-2 2，N-1状态是0 2。
还值得一提的是，Multiwfn亦能够计算亲核和亲电超离域度（superdelocalizabilities），类似于能量权重的基于轨道计算的福井函数。这个量最早由Schüürmann在Environ. Toxicof. Chem., 9, 417 (1990)和Quant. Struct.-Act. Relat., 9, 326 (1990)中提出，被广泛用于构建QSAR方程。但原文的计算形式只适合基于半经验方法进行计算，而Multiwfn使用的计算形式是笔者对之进行广义化后的，也能够用于HF和DFT波函数。具体细节见Multiwfn手册3.25.5节，用户只需要提供分子原本状态的fch、molden、mwfn之类含有基函数信息的文件，进入主功能22后选择选项8即可得到各个原子的超离域度。
作者Author: nianbin 时间: 2019-9-18 16:06
很好
作者Author: Villars 时间: 2020-2-12 10:26
社长，您好！我按照您上文中3.1节的操作步骤得到了N.gjf、N+1.gjf和N-1.gjf三个文件，然后将这三个文件放到超算平台上计算得到N.wfn、N+1.wfn和N-1.wfn。但是我将wfn文件放回到原来的路径下，使用2 Calculate various quantitative indices功能没有得到CDFT.txt文件。请问社长是因为Win版的Multiwfn不识别linux下计算出的wfn文件吗？
作者Author: sobereva 时间: 2020-2-12 10:45
看屏幕上的提示。
所有平台上Gaussian算出来的wfn文件内容都是相同的
作者Author: Villars 时间: 2020-2-12 11:21
本帖最后由 Villars 于 2020-2-13 00:17 编辑
感谢社长的回复！我在“2 Calculate various quantitative indices”这步的后面输入2以后，Multiwfn就直接弹出了，然后我在包含Multiwfn.exe文件的路径下确实没有找到您提到的CDFT.txt文件，您看这种情况可能是哪里出的问题？
作者Author: Villars 时间: 2020-2-12 13:52
本帖最后由 Villars 于 2020-2-12 13:59 编辑
社长您好！我刚才用您给的例子phenol.xyz，按照上文中的步骤重新计算了一遍。电荷和自旋多重度的设置都使用了默认设置，生成了N.gjf、N+1.gjf和N-1.gjf三个文件。然后在win版的G16下计算，生成了N.wfn、N+1.wfn和N-1.wfn三个文件。这次计算Multiwfn没有直接弹出，但是还是没有得到您给的CDFT文件。具体屏幕显示和CDFT文件报错如图所示。我看到上面显示没有找到系统能量，可是我打开wfn文件，上面SCF Done后面有显示单点能。您看这个可能是什么原因？
作者Author: Villars 时间: 2020-2-13 00:02
社长您好！问题终于解决了，原来是我在超算上把输出文件命名为wfn文件了，就是写成 g16 <xxx.gjf> xxx.wfn了，然后又把这个xxx.wfn放到multiwfn下计算才会出现这种报错，真是大错特错！幸好看了您的“详谈Multiwfn支持的输入文件类型、产生方法以及相互转换” 这篇文章才及时醒悟，谢谢社长！
作者Author: sobereva 时间: 2020-3-13 18:17
今日更新了Multiwfn，支持了J. Phys. Chem. A, 124, 2090 (2020)提出的ωcubic指数，对于考察弱相互作用很有益，在本文新增了相关介绍。
作者Author: 齐帆 时间: 2020-3-15 20:07
社长，您好。看有文献上说是电离电势IP的值和HOMO轨道的能量值是互为相反数的关系，那为什么给出的实例中两者不满足这种关系。自己算的分子中也不满足，一直很疑惑。希望社长帮助解答一下啦
作者Author: snljty 时间: 2020-3-15 21:11
本帖最后由 snljty 于 2020-3-15 21:22 编辑
先借问卢老师一下，上一次基础班讲义的260面的第2张PPT，编号90的那张的图，能麻烦您发一下么？就是那个循环伏安法测的性的氧化/还原势只能说是和理论计算的HOMO、LUMO能量有一定内在联系...那张，带有好几个箭头那张，想存一下那张图片。谢谢您！
在论坛上搜一下 电离势 sob老师以前贴过好几次那几张讲义中的图。摘自sob老师的讲义。
(, 下载次数 Times of downloads: 326)
(, 下载次数 Times of downloads: 319)
作者Author: sobereva 时间: 2020-3-17 13:31
那叫Koopmans定理。Google一搜资料一堆。对于精确的交换相关泛函是满足的，但是对于近似的交换相关泛函都不精确满足，只不过有的满足好，有的满足不好。QTP系列泛函算是满足最好的之一。
作者Author: sobereva 时间: 2020-3-17 14:17
(, 下载次数 Times of downloads: 1172)
作者Author: snljty 时间: 2020-3-17 14:27
谢谢老师！
作者Author: sobereva 时间: 2020-10-14 07:25
新更新
从2020-Oct-14更新的Multiwfn开始，不再强制用当前目录下的N.wfn、N+1.wfn和N-1.wfn作为后续计算的输入文件。用选项2或3开始计算时，若Multiwfn发现当前目录下缺少这些文件的一个或多个，程序就会让用户自行输入相应电子态的.wfn或.wfx或.fch或.mwfn文件的路径。因此，比如你想用自行算出来的N.fch、N+1.fch、N-1.fch做后面的计算也完全可以。
作者Author: sobereva 时间: 2021-6-21 03:15
今日更新的Multiwfn新支持了计算轨道超离域度，在本文文末增加了下面文字
还值得一提的是，Multiwfn亦能够计算亲核和亲电超离域度（superdelocalizabilities），类似于能量权重的基于轨道计算的福井函数。这个量最早由Schüürmann在Environ. Toxicof. Chem., 9, 417 (1990)和Quant. Struct.-Act. Relat., 9, 326 (1990)中提出，被广泛用于构建QSAR方程。但原文的计算形式只适合基于半经验方法进行计算，而Multiwfn使用的计算形式是笔者对之进行广义化后的，也能够用于HF和DFT波函数。具体细节见Multiwfn手册3.100.16.5节，用户只需要提供分子原本状态的fch、molden、mwfn之类含有基函数信息的文件，进入主功能100的子功能16后选择选项8即可得到各个原子的超离域度。
作者Author: lei234 时间: 2021-7-7 20:54
我是从fchk转化到wfn，现在有了N.N+1,N-1三个的wfn文件，并放在了一个文件夹中，请问怎么开始计算？您文中说的从当前模块选2 calculate various quantitative indices在哪？求指点，感谢！
作者Author: sobereva 时间: 2021-7-8 12:25
把此文完整、通读一遍
使用Multiwfn超级方便地计算出概念密度泛函理论中定义的各种量
http://sobereva.com/484（http://bbs.keinsci.com/thread-13328-1-1.html）
当然是指当前计算概念密度泛函相关的量的模块
(, 下载次数 Times of downloads: 300)
一个字一个字阅读
作者Author: sobereva 时间: 2021-7-15 20:30
今日更新的版本此功能已支持直接调用ORCA进行计算，修改了相应文字
另外，如果你没有Gaussian而只有ORCA，也可以选择-2 Choose the quantum chemistry program used in option 1之后再选2。之后再选选项1产生输入文件时，Multiwfn产生的就是三个态的ORCA的输入文件N.inp、N+1.inp和N-1.inp，手动算完了之后就会在当前目录下产生N.wfn、N+1.wfn和N-1.wfn（其它产生的文件可删掉）。如果你把settings.ini里的orcapath设为了实际ORCA可执行文件路径，那么Multiwfn在产生输入文件后会问你是否自动调用ORCA去计算它们并在计算结束后删掉其它文件。
作者Author: 斯斯 时间: 2021-11-10 00:13
普通有机体系用常用的B3LYP-D3(BJ)/6-31G*级别一般就够了，如果想更便宜一些，用半经验、GFN-xTB之类也可以。
您好，想问一下社长，为什么在算cdft的时候，对基组没什么要求呢，有什么文献可以支持这种说法么。（审稿人要求我用更高的机组对比一下，我用的就是B3LYP/6-31G*算的有机化合物）
作者Author: sobereva 时间: 2021-11-10 02:06
我从来没说“如果想更便宜一些，用半经验、GFN-xTB之类也可以”
GFN-xTB做个预优化还行，实际的优化，但凡算得动的都不要直接用GFN-xTB。对于得到波函数、算能量，更是绝对不能用GFN-xTB，实在太糙
没说对基组没什么要求，只是对基组要求不高，6-31G*就行了，但是如果你降到6-31G也是绝对不能接受的。
碰到了外行审稿人让你对比你就找一个体系对比一下就完了，结果肯定和6-31G*没有定性差异。
作者Author: 斯斯 时间: 2021-11-10 09:43
好的，谢谢社长
作者Author: 斯斯 时间: 2021-11-11 21:24
社长，想请教您一下，我按着这个方法进行方法的对比，来回复审稿人，可以么。
载入B3LYP-D3（BJ)/6-311G**优化后的文件mol2
100 //其它功能（Part 1）
16 //计算CDFT里定义的各种量的模块
1 //产生波函数文件
[直接按回车]
[直接按回车]
现在当前目录下已经有了N.gjf、N+1.gjf和N-1.gjf，更改里面的基组一个是B3LYP/6-31g*对比方法为M06-2x/6-31g*
然后高斯计算得到wfn文件。对比两个方法得到的数据
作者Author: sobereva 时间: 2021-11-12 00:16
“更改里面的基组一个是B3LYP/6-31g*对比方法为M06-2x/6-31g*” 明显病句。明明这俩计算级别的基组都是6-31G*。
gjf里的计算级别想自己改可以自己改
作者Author: 斯斯 时间: 2021-11-12 10:39
对的，应该是计算方法的不同，是我说错了，
所以我在B3LYP-D3（BJ)/6-311G**优化后的分子基础上，两个方法的比对应该是B3LYP/6-311G**和M062X/6-311G**下分别计算，进行比对，是这样么,
看您的另一片博文上说的与研究内容、目的有关。比如要对比研究一系列泛函A、B、C、D对于某类分子的极化率、NMR、偶极矩之类问题的计算精度，那就没必要用这四种泛函分别去优化结构和计算性质了，只需要用一个对当前体系合理的泛函去优化即可，之后计算性质的时候再切换成不同泛函去分别计算
作者Author: sobereva 时间: 2021-11-13 03:35
这么对比是可以的
作者Author: 斯斯 时间: 2021-11-13 18:43
好的，谢谢社长大人
作者Author: Hear 时间: 2022-7-19 00:07
老师，您好。我想问下我在算简缩福井函数的时候，优化结构的文件里除了有基组和方法外，可以加scrf=(smd,solvent=water)，还有关于经验色散矫正的这种命令嘛？如果可以加的话，算单点的时候是不是也得加上和优化的输入文件保持一致啊。
作者Author: sobereva 时间: 2022-7-19 14:15
溶剂模型可以加
色散校正如果是诸如DFT-D3这样只校正能量，而且校正量只依赖于几何结构的，加不加结果都完全一样
作者Author: Hear 时间: 2022-7-19 15:27
老师，您好！我现在还想算一下酶环境下的福井函数。我是在turbomole上用的QM/MM方法，1.我想的是优化完结构以后，在算单点能的时候改一下电荷和自旋多重度，分别算N、N+1、N-1状态下的单点，算完之后产生相应的wfn文件，之后再用Multiwfn软件去做。2.在用QM/MM时，是不是应该把QM区域设成只有底物呢？
作者Author: sobereva 时间: 2022-7-19 15:50
是
作者Author: Hear 时间: 2022-7-19 16:45
老师，您好！我发现我把QM区域改为只有底物以后，单点算不了了。我可以把优化好的底物结构，放到高斯里算N,N+1,N-1状态下的单点嘛，这样算出来的是不是就不算是在酶环境下的结果啦
作者Author: sobereva 时间: 2022-7-21 11:34
这只能体现酶环境对分子构象的影响，但体现不了酶环境的势场对底物电子结构的影响。要体现后者，至少要把酶的原子用背景电荷表现出来
作者Author: sobereva 时间: 2024-2-24 06:18
给本文增加了计算局部超软度（local hyper-softness）的说明
作者Author: LJZZZZ 时间: 2024-2-29 22:35
社长您好，我在生成了3个N，N+1,N_1文件，按照步骤设置了环境变量，结果我一输入字母,快速过了一些代码然后自动关闭，我在相应的路径下就没有看到其他文件的生成
作者Author: sobereva 时间: 2024-2-29 23:25
Multiwfn FAQ
http://sobereva.com/452
上文里搜 闪退
作者Author: LJZZZZ 时间: 2024-3-1 00:28
社长你好，我检查了我的gif文件，最后一行也加上了生成wfn文件的路径，但在muti中还是生成不出来。
作者Author: sobereva 时间: 2024-3-1 00:40
分清楚gif和gjf
没有叫muti的东西
自己看Gaussian输出文件内容判断，看不懂就把输入输出文件上传
作者Author: LJZZZZ 时间: 2024-3-1 00:45
不好意思，muti wfn他不让我连续打出来，我用高斯输出wfn文件时总是卡在这个I302，而且还是新换的电脑，麻烦您帮我看一下
作者Author: LJZZZZ 时间: 2024-3-1 00:47
这是我的gjf文件
作者Author: sobereva 时间: 2024-3-1 01:32
因为那叫Multiwfn，不叫mutiwfn。论坛不允许乱拼程序名的情况出现
默写10遍Multiwfn，并且认真看下文
计算化学中的一些常见不良写法和用词
http://sobereva.com/298（http://bbs.keinsci.com/thread-1358-1-1.html）
作者Author: sobereva 时间: 2024-3-1 01:34
gjf文件内容本身没问题
远程说不清楚，Multiwfn直接调用windows下Gaussian如果没法正常运行，就自己手动用Gaussian一个一个算gjf文件
作者Author: LJZZZZ 时间: 2024-3-1 13:09
收到谢谢社长！
作者Author: LJZZZZ 时间: 2024-3-1 13:10
我是用Gaussian一个一个算，结果卡在I302.exe了，一直没有办法解决。
作者Author: LJZZZZ 时间: 2024-3-1 13:17
可以看一眼我的out文件，一直卡在了I302.exe，没有办法解决，我用的是09W这个版本，好像是说cpu跟软件不兼容，这个我也不太懂，不知道换16版本能不能解决
作者Author: sobereva 时间: 2024-3-2 06:30
Gaussian的安装方法及运行时的相关问题
http://sobereva.com/439（http://bbs.keinsci.com/thread-10814-1-1.html）
这里该说的都说了
用64bit windows版Gaussian或者用linux版
作者Author: LJZZZZ 时间: 2024-3-2 22:19
社长您好，今天发现Maximum Displacement 0.003030 0.001800 NO，没有收敛，请问有什么好的解决方法吗。
作者Author: sobereva 时间: 2024-3-2 23:12
跟阈值足够接近，可以接受
下文也说了
Gaussian中几何优化收敛后Freq时出现NO或虚频的原因和解决方法
http://sobereva.com/278（http://bbs.keinsci.com/thread-633-1-1.html）
作者Author: LJZZZZ 时间: 2024-3-2 23:57
但听说投稿的没有4个yes会被拒稿，是真的嘛社长
作者Author: zako 时间: 2024-3-3 01:32
Why Multiwfn doesn't support Interacting Quantum Atoms (IQA) analysis ? It's s very interesting feature.
作者Author: sobereva 时间: 2024-3-3 03:56
zako 发表于 2024-3-3 01:32
Why Multiwfn doesn't support Interacting Quantum Atoms (IQA) analysis ? It's s very interesting feat ...
It is extremely expensive, while I don't find it is very useful
作者Author: wzkchem5 时间: 2024-3-3 17:49
没有这码事。那是有的人误以为只有四个YES高斯才认为算收敛，又想给别人强调收敛的重要性，所以故意把问题说得更严重、吓唬人的。我不信有人发现Maximum Displacement虽然是NO、但离YES这么近，就会拒稿，顶多提出小修，要求作者把结构完全收敛。
一定要注意，高斯判断收敛的条件不是四个YES，如果有的NO足够接近收敛了，且其他的YES收敛得比收敛限还要严好多，那么这些YES可以弥补那些NO的不足，导致高斯仍然认为计算收敛了。关键在于看高斯在输出文件里怎么说，有没有说找到了稳定点，而不是看有几个YES。如果高斯说收敛了，那么即便没有四个YES也能用；如果高斯说没收敛，且你的结果对结构精度比较敏感，那再按sob老师说的方法消掉NO。
说实话，我觉得没做构象搜索被拒稿的概率可能还大一些。很多人花了很大气力证明自己找到的局部极小值点结构确实是局部极小值点，但可能也只比没完全收敛准0.1kcal/mol而已，而不做构象搜索导致用错误的局部极小值点结构做计算，很容易导致好几个kcal/mol的误差。
作者Author: LJZZZZ 时间: 2024-3-4 17:09
谢谢谢谢，感谢，我知道啦
作者Author: LJZZZZ 时间: 2024-3-6 00:53
社长请问那一段生成N.wfn文件的操作是怎么弄呀，您看一下我的out文件截图
作者Author: sobereva 时间: 2024-3-6 03:26
明显N.wfn已经生成了
如果当前目录找不到，去Gaussian的scratch目录下找
作者Author: LJZZZZ 时间: 2024-3-6 19:50
他给我保存到16w目录下了
作者Author: LJZZZZ 时间: 2024-3-12 22:10
我想知道一下，Gaussian16w会不会存在程序一直跑，但不报错，我看log文件一直没有收敛，已经跑了一个多星期了，文件大小也有1万多kb了，因为我是做蛋白的，用的是4肽进行DFT计算，之前一个4肽花了4天就好了，最后log文件也才5000多kb,这次4肽跑了这么久却没有好，我想问是不是出了什么问题。
作者Author: wzkchem5 时间: 2024-3-12 22:30
是SCF不收敛还是结构优化不收敛？
两个分子结构非常类似，但是一个收敛需要的步数是另一个的好几倍，也是常有的事情
作者Author: LJZZZZ 时间: 2024-3-12 22:54
out文件太大传不上去，这是我的gjf文件
作者Author: LJZZZZ 时间: 2024-3-12 22:56
已经收敛200次左右了，一开始还有两个yes,现在是不是没希望了
作者Author: LJZZZZ 时间: 2024-3-13 11:38
社长请问一下，算法的精度会不会影响收敛，比如我用一个简单的算法收敛了200次左右都没成功，换一个精度更高的算法会不会好一点
作者Author: wzkchem5 时间: 2024-3-13 16:27
displacement太大了，可以考虑缩小最大步长
方法参见http://sobereva.com/164
作者Author: sobereva 时间: 2024-3-14 10:15
“算法”指代不明
直接说清楚用的关键词
那不叫“收敛了200次”，对于SCF那叫“迭代了200次没收敛”，对于几何优化那叫“跑了200步没收敛”，用词必须严谨
SCF不收敛、几何优化不收敛该怎么解决怎么解决，仔细看下文，换计算级别绝对不是首先要考虑的做法
解决SCF不收敛问题的方法
http://sobereva.com/61
量子化学计算中帮助几何优化收敛的常用方法
http://sobereva.com/164
每次说“收敛”的时候必须说清楚是对于什么而言的
作者Author: LJZZZZ 时间: 2024-3-14 15:55
好的社长，我还想请问一下Gaussian16W可以同时进行多个DFT计算吗
作者Author: LJZZZZ 时间: 2024-3-14 21:56
社长，我的Gaussian16w卡在这5个多小时，是出了什么问题吗，一直不动
作者Author: wzkchem5 时间: 2024-3-15 05:20
上传完整输入文件。由截图不能唯一确定原因
作者Author: LJZZZZ 时间: 2024-3-15 14:18
完整的太大了，上传不了，有办法解决嘛
作者Author: wzkchem5 时间: 2024-3-15 16:37
这就是为什么我说的是输入文件而非输出文件。
输入文件总不至于大到上传不了吧？如果连输入文件都传不上来，至少给出输入文件的前10行
作者Author: sobereva 时间: 2024-3-15 17:36
可以同时跑多个
作者Author: sobereva 时间: 2024-3-15 17:36
输入不可能大到上传不了，仔细看置顶的社员必读贴了解上传文件常识
作者Author: LJZZZZ 时间: 2024-3-15 21:51
社长请问可以分享一下Gaussian16的Linux版本安装包嘛
作者Author: LJZZZZ 时间: 2024-3-15 23:31
请问有Gaussiian16linux版本安装包吗，可以分享一下吗
作者Author: wzkchem5 时间: 2024-3-15 23:34
要么从官方渠道买正版，要么用ORCA等免费软件。世界上不是只有高斯一个量化软件
作者Author: sobereva 时间: 2024-3-16 07:02
论坛置顶的社员必读贴明确说了不得公开求商业程序。下次再求封号处理
不得在本论坛公开发帖求商业性、收费程序，此类帖子会被删帖，并且记违规。
作者Author: LJZZZZ 时间: 2024-3-22 22:34
社长我想请问一下，对同一物质进行相同的多次DFT计算（设置的参数都一样），得到的HOMO轨道和LUMO轨道的能量值可能不一样嘛，每次结果。
作者Author: sobereva 时间: 2024-3-23 09:21
有可能有极轻微不同，见下文
数值误差对计算化学结果重现性的影响
http://sobereva.com/88
作者Author: LJZZZZ 时间: 2024-3-24 11:55
社长，请问一下，我看到很多人通过计算HOMO轨道和LUMO轨道能量差，预测其抗氧化活性，随后在体外进行ABTS和DPPH自由基清除试验，发现其预测结果和体外实验结果一直，请问这个有依据吗？或者说可不可能存在预测结果与体外实验结果不一致的情况，还有ABTS和DPPH是两个不同溶液体系，那我进行DFT计算时，是不是得在不同溶液体系都计算一遍
作者Author: sobereva 时间: 2024-3-25 11:06
这本质上相当于fundamental gap（等价于电子硬度）的近似，远不如用fundamental gap说事来得严格
一致或不一致都有可能，取决于影响结果的机理。对比的体系越相似、其它影响反应的因素越少（位阻、扩散等）越可能有较好相关性
不同溶液需要分别算
作者Author: LJZZZZ 时间: 2024-4-8 01:43
社长您好，我想请问描述一种多肽分子具有5.24eV的高能隙，具有高度稳定性，这是作为抗氧化剂的理想特性。这种描述有道理吗，我在一篇文献中看到的
作者Author: sobereva 时间: 2024-4-8 02:45
光从这上面不能说明什么
下文有更多讨论
正确地认识分子的能隙(gap)、HOMO和LUMO
http://sobereva.com/543（http://bbs.keinsci.com/thread-16758-1-1.html）
作者Author: 九曜 时间: 2024-4-9 09:48
请问老师，例如计算福井函数，对DFT的计算级别要求高么，大基组下获得的结果会更准确吗。
作者Author: sobereva 时间: 2024-4-9 10:12
对于普通中性有机分子，通常B3LYP/6-31G*就够
作者Author: LJZZZZ 时间: 2024-4-15 01:06
老师请问，方法设定中的Spin参数怎么计算，像我想计算多肽，那我该如何去计算它的自旋多重度
作者Author: sobereva 时间: 2024-4-15 06:19
自旋多重度，这是量化计算最基本的常识。不具备这种常识性知识的话强烈建议参加北京科音初级量子化学培训班（http://www.keinsci.com/workshop/KEQC_content.html）扫盲，里面还专门讲了自旋多重度的判断。
对于普通有机、生物分子体系，中性非自由基状态下自旋多重度都是1
作者Author: LJZZZZ 时间: 2024-4-15 09:52
但您看这里好像最低都是2，那我该如何设置呀
作者Author: sobereva 时间: 2024-4-16 00:58
拿具体体系说事，起码摆出来结构
作者Author: sobereva 时间: 2024-4-17 07:40
今日更新的Multiwfn的CDFT模块已经支持了亲电描述符ε的计算，并相应地更新了本文。此描述符于Int. J. Quantum Chem., 124, e27366 (2024)中提出，文中通过35个有机分子测试发现它与Mayr亲电性参数的相关性远好于Parr的亲电指数ω。
作者Author: LJZZZZ 时间: 2024-5-7 20:22
这是我的out文件，突然报2070错误，社长能不能帮我看看怎么解决
作者Author: sobereva 时间: 2024-5-8 20:19
内存不够
Gaussian计算怎么设置内存，仔细看
Gaussian的安装方法及运行时的相关问题
http://sobereva.com/439（http://bbs.keinsci.com/thread-10814-1-1.html）
作者Author: LJZZZZ 时间: 2024-5-8 20:26
是设置程序时候分配的内存嘛，我之前是1GB，现在弄了2GB重新跑
作者Author: sobereva 时间: 2024-6-1 16:07
Multiwfn现在还可以计算在反应位点预测方面比双描述符更严格的双描述符势（dual descriptor potential），见《使用Multiwfn计算双描述符势预测反应位点》（http://sobereva.com/708）。
作者Author: LJZZZZ 时间: 2024-6-26 22:41
社长您好，我想请问一下如何在out文件中找到我需要的键长和每个原子所带电荷的数据呀
作者Author: LJZZZZ 时间: 2024-6-26 23:15
我想问一下右边标注每个原子的图是用什么软件画的呀
作者Author: sobereva 时间: 2024-6-27 02:54
估计是自己ps的
作者Author: sobereva 时间: 2024-6-27 02:55
Multiwfn的概念密度泛函理论计算模块不输出键长信息
想得到所有键的键长，用主功能0图形界面上方的Tools - Export all internal coordinates
Multiwfn的概念密度泛函理论计算功能输出的CDFT.txt直接就有各个原子的Hirshfeld电荷。如果要算其它原子电荷，用主功能7里的相应选项计算相应原子电荷，如ADCH、EEM等
作者Author: tangdu 时间: 2025-6-15 13:28
sobereva老师您好，请问若当前体系是个+1电荷的闭壳层阳离子，可以算双描述符势嘛？双描述符势可以定量看到分子中每个原子的情况嘛
作者Author: sobereva 时间: 2025-6-15 17:11
可以
不能
作者Author: 121089576 时间: 2025-7-7 19:30
本帖最后由 121089576 于 2025-7-8 15:10 编辑
有一个问题，Multiwfn生成的N、N+1、N-1文件中未包含色散校正和溶剂模型，如果实际是溶剂环境，是否应该手动在文件中添加溶剂模型？难道CDFT指数的i计算不需要考虑色散校正和溶剂模型？
-----------更新
关于溶剂模型问题已解决，搜到了sob老师在别的帖子中的回复，在此贴出
考察溶剂环境下的情况就带上。是否更准确不一定，看具体情况，但至少原理上更合理
但是关于色散校正问题依然不理解，博文中sob老师也是直接用的B3LYP/6-31G*做的演示，这个泛函不加色散校正不是完全不能用吗？
-----------更新
经实测，是否添加色散校正对CDFT计算无影响，得到的值都一样。
作者Author: sobereva 时间: 2025-7-9 02:48
得明确是什么色散校正。有的色散校正诸如VV10是依赖于电子密度的，这会影响诸如福井函数的计算结果（尽管影响甚微）。而最常用的DFT-D3只依赖于几何结构，N、N-1、N+1态用的几何结构是相同的，所以D3校正不会产生丝毫影响
| 欢迎光临 计算化学公社 (http://bbs.keinsci.com/) |
Powered by Discuz! X3.3 |