计算化学公社

标题: 调用Multiwfn计算分子间电荷转移积分的程序[测试版] [打印本页]

作者
Author: snljty 时间: 2021-5-30 14:00
标题: 调用Multiwfn计算分子间电荷转移积分的程序[测试版]
本帖最后由 snljty 于 2023-9-12 13:05 编辑

更新：加入了一个新功能，可以提供一个dimer的NBO的.47文件，只要里面包含Fock矩阵，程序可以从里面读取Fock矩阵，就不用折腾了。对于Gaussian用户，读取上一次dimer单点的chk文件后，写Guess=Read Population=NBORead，结尾空一行写$NBO ARCHIVE FILE= MYDIMER $END，重跑一遍，就会在工作目录下生成MYDIMER.47文件，然后执行本程序的时候加上-d47选项和这个文件的名字即可。这样就不需要提供下文的Fock矩阵了。

更新：刚才改了几个typo，并加入了一个新功能。当数值误差原因导致系数矩阵不可逆的时候，可以用-dF选项提供一个含有完整的dimer的Fock矩阵的文件，将直接从里面读取Fock矩阵。这个文件第一行应该是注释，后面应该是紧跟一个下三角矩阵记录的dimer的Fock矩阵。对于Gaussian用户，可以在计算一次dimer的单点后，读取上一次单点计算的chk文件，关键词加上Guess=Read IOp(5/33=3)，然后找输出文件最后一次出现Fock matrix (alpha)或者Fock matrix (beta)的地方，就是Fock矩阵，带着这行标题拷贝出来到新文件即可。

顺便可以猜一下，这个程序为什么要用Fortran写。

写在最前面：
此程序在若干个体系测试过，但是未进行更广泛的测试。另外在某些基函数有接近线性依赖，但是还没有被量化软件判断为线性依赖的情况下，可能由于Multiwfn输出的数据精度有限，读取导致的一些数值误差问题，出现系数矩阵不满秩进而不可逆的情况。

正文：
论坛中已经有很多计算电荷转移积分的程序了，比如http://bbs.keinsci.com/thread-11369-1-1.html以及http://bbs.keinsci.com/thread-14311-1-1.html。这些程序基本都是基于Gaussian的。考虑到Multiwfn极为丰富的波函数读取与操作功能，这里（2个月前写好的，因为测试较少没有发上来）写了一个利用卢老师的Multiwfn程序（http://sobereva.com/multiwfn/）读取轨道信息，然后读取Multiwfn的输出文件并使用格点能修正法计算分子间电荷转移积分的程序。程序地址在https://github.com/snljty/charge_transfer_integral。这个程序是用Fortran写的，Fortran 95的语法。

输入文件是dimer以及两个monomer各自包含GTO和MO信息的文件，诸如Gaussian的fch，ORCA的molden等文件。使用方法可以交互输入，也可以用命令行参数比如./Calc_coupling.x -d dimer.fch -m1 monomer1.fch -m2 monomer2.fch。注意monome1r的坐标应该是dimer坐标的前面部分，monomer2的坐标应该是dimer坐标的后半部分，两者坐标要和其在dimer中的坐标完全相同，而且顺序不可改变。另外注意Gaussian要加上NoSymmetry关键词。如果可能有基函数线性依赖关系，量化软件去掉这些成分后会导致系数矩阵为非方阵，不存在逆矩阵从而无法用这个方法计算。即使量化程序用某些关键词指定不去除线性依赖的基函数成分，后面对系数矩阵求逆也会报错。

bin目录里.exe拓展名的是Windows下的可执行文件，.x是WSL下编译的Linux版可执行文件。

使用前请务必确认Multiwfn已经正确安装且配置好，已经将Multiwfn可执行文件所在目录加入环境变量PATH，使得命令行输入Multiwfn即可启动Multiwfn，并且尽量用较新版本Multiwfn。如果要自行编译，还需要链接包含LAPACK和BLAS内函数的数学库。因为这个程序通常没什么计算量，图省事我直接链接LAPACK编译的。

这里放一下README.md文件和Windows下的可执行文件，还有几个测试文件。其余内容包括源代码在上面的github地址。

欢迎测试和反馈。

(, 下载次数 Times of downloads: 370)
(, 下载次数 Times of downloads: 242)
(, 下载次数 Times of downloads: 242)

作者
Author: pyyu 时间: 2021-6-7 15:46
作者能把linux版本也放上来吗

github总是进不去

作者
Author: snljty 时间: 2021-6-7 16:34
本帖最后由 snljty 于 2021-6-7 16:45 编辑

pyyu 发表于 2021-6-7 15:46
作者能把linux版本也放上来吗github总是进不去

没有问题，链接已经发给你了。但是我的Linux版gcc缺一些静态库，有些动态库和你电脑上可能不太一样。如果有问题欢迎再反馈。

作者
Author: zzffzz33 时间: 2021-7-6 15:04
楼主你好，我使用的时候输入dim 以及 mon后，运行了三秒左右直接闪退了，请问应该是什么原因啊？是没有正确调用multiwfn吗？

作者
Author: snljty 时间: 2021-7-6 15:10

zzffzz33 发表于 2021-7-6 15:04
楼主你好，我使用的时候输入dim 以及 mon后，运行了三秒左右直接闪退了，请问应该是什么原因啊？是没有正确 ...

终端里跑应该不会闪退吧，至少终端不会退出。麻烦给我截图看下屏幕上的全部输出。如果方便的话，可以把文件发给我测试一下。

作者
Author: zzffzz33 时间: 2021-7-6 16:43

snljty 发表于 2021-7-6 15:10
终端里跑应该不会闪退吧，至少终端不会退出。麻烦给我截图看下屏幕上的全部输出。如果方便的话，可以把文 ...

感谢您的回复，我用的就是test文件。把test里的三个out文件放在程序目录下，然后按提示敲了下三次回车之后闪退。我是不是遗漏了什么步骤，本人比较小白，

作者
Author: snljty 时间: 2021-7-6 17:21

zzffzz33 发表于 2021-7-6 16:43
感谢您的回复，我用的就是test文件。把test里的三个out文件放在程序目录下，然后按提示敲了下三次回车之 ...

启动一个cmd.exe，在cmd.exe里面cd到文件所以目录，执行calc_coupling.exe，然后照常操作，然后麻烦把输出信息截屏或者拷贝下来发给我看。

作者
Author: indec 时间: 2021-7-7 12:10
请问可以支持 HOMO-1 和 HOMO-1 之间的电荷转移积分吗？

作者
Author: snljty 时间: 2021-7-7 12:47

indec 发表于 2021-7-7 12:10
请问可以支持 HOMO-1 和 HOMO-1 之间的电荷转移积分吗？

可以，有个选项是输出所有轨道之间的转移积分。你可以用calc_coupling.exe --help查看帮助菜单，然后看一下有个--full选项，然后从输出文件中找就行。只输出指定轨道间的功能没有做。

作者
Author: MAX-0309 时间: 2021-7-17 10:38
楼主，你好，我是用您开发的程序进行电荷转移积分计算，两个monomer之间的HOMO电荷转移积分可以算，LUMO电荷转移计算出现的是这个结果？请问是为什么呢？谢谢

作者
Author: snljty 时间: 2021-7-17 15:21

MAX-0309 发表于 2021-7-17 10:38
楼主，你好，我是用您开发的程序进行电荷转移积分计算，两个monomer之间的HOMO电荷转移积分可以算，LUMO电 ...

Windows下算的么？方便把文件发给我测试一下么？谢谢

作者
Author: MAX-0309 时间: 2021-7-17 16:48

snljty 发表于 2021-7-17 15:21
Windows下算的么？方便把文件发给我测试一下么？谢谢

好的，麻烦您啦

作者
Author: xmdu 时间: 2021-8-31 15:13
有使用教程吗

作者
Author: snljty 时间: 2021-9-13 11:50

xmdu 发表于 2021-8-31 15:13
有使用教程吗

本文一楼已经写明白了全部需要的操作指南，要注意的容易忽视的点还特意用红字标出来了，还贴上了接近1000字的README.md作为手册，程序启动直接用--help参数也可以查看好几十行不能更详细的使用说明。连测试文件都一起提供了，我真不知道你连着问我两遍是还想要什么样的教程。

作者
Author: snljty 时间: 2021-9-13 11:52
加了个功能，可以直接从包含dimer的Fock矩阵的NBO .47矩阵读dimer的Fock矩阵，这个比之前读Fock矩阵的方式更方便，也没有逆解Fock矩阵可能遇到的数值问题。对于Gaussian用户的使用方法上面已经提到，也可以用--help手册查看可执行程序或者看README.md里的说明。

作者
Author: zzffzz33 时间: 2021-10-24 20:56

snljty 发表于 2021-7-6 17:21
启动一个cmd.exe，在cmd.exe里面cd到文件所以目录，执行calc_coupling.exe，然后照常操作，然后麻烦把输 ...

非常抱歉老师。忘记回复您了，我用CMD指令打开可以正常算，就是直接点击软件运行一段后会闪退。不过是可以算了，非常感谢

作者
Author: zzffzz33 时间: 2021-10-25 11:14
请问下楼主，对于激子解离和复合的coupling能否实现？

作者
Author: pizaipizai@5989 时间: 2021-11-2 14:01
坐标不改变，在计算的时候怎么实现呢？我试过nosym和给坐标加$方法，程序都会自动把分子重新换个位置然后去scf

作者
Author: snljty 时间: 2021-11-2 15:15

zzffzz33 发表于 2021-10-25 11:14
请问下楼主，对于激子解离和复合的coupling能否实现？

等我有时间看一下，可能至少得一年...

作者
Author: snljty 时间: 2021-11-2 15:15

pizaipizai@5989 发表于 2021-11-2 14:01
坐标不改变，在计算的时候怎么实现呢？我试过nosym和给坐标加$方法，程序都会自动把分子重新换个位置然后去 ...

请问你用的软件是...?

作者
Author: pizaipizai@5989 时间: 2021-11-2 16:56

snljty 发表于 2021-11-2 15:15
请问你用的软件是...?

我用的orca，我同学跟我讲高斯可以nosym一下貌似就不会变坐标，不过我还没试，想着先来请教一下撸主怎么弄的

作者
Author: snljty 时间: 2021-11-2 18:03

pizaipizai@5989 发表于 2021-11-2 16:56
我用的orca，我同学跟我讲高斯可以nosym一下貌似就不会变坐标，不过我还没试，想着先来请教一下撸主怎么 ...

我记得ORCA不会自己平移旋转坐标啊。两个monomer的任务把另外一个分子的坐标直接删去试试？

作者
Author: pizaipizai@5989 时间: 2021-11-2 21:40

snljty 发表于 2021-11-2 18:03
我记得ORCA不会自己平移旋转坐标啊。两个monomer的任务把另外一个分子的坐标直接删去试试？

我也很迷惑，我目前是这么解决的，比如用dimer算，dimer.inp是一套坐标，dimer.out又是一套坐标；然后算monomer的时候，再从dimer.inp里面把单体坐标抠出来，然后算完了观察一下monomer.out，发现跟dimer.out能对上。但还没有用撸主的程序算积分，我先把scf的bug给de了就去试

作者
Author: pizaipizai@5989 时间: 2021-11-3 00:59
撸主，我看上面有一哥们直接运行exe会闪退，我也直接运行了一下，程序走了下面的步骤然后也闪退了。留给你debug用
# Input the name of file contains MO infomation of the dimer.
# If press <Enter> directly, "dimer.fchk" will be used.

# Neither "dimer.fchk" nor "dimer.fch" found!
# Input again.
dimer.molden.input
# Input the name of file contains MO infomation of the 1st monomer.
# If press <Enter> directly, "monomer1.fchk" will be used.
mono1.molden.input
# Input the name of file contains MO infomation of the 2nd monomer.
# If press <Enter> directly, "monomer2.fchk" will be used.
mono2.molden.input
# Using "dimer.molden.input" for dimer.
# Using "mono1.molden.input" for monomer1.
# Using "mono2.molden.input" for monomer2.
# Calling Multiwfn to obtain amount of orbitals and electrons of monomer1 ...
# Reading amount of orbitals and electrons of monomer1 ...
# Calling Multiwfn to obtain amount of orbitals and electrons of monomer2 ...
# Reading amount of orbitals and electrons of monomer2 ...
# Calling Multiwfn to obtain amount of orbitals and electrons of dimer ...
# Reading amount of orbitals and electrons of dimer ...
# number of electrons of monomer1/monomer2/dimer:
#                            126    126 252
# number of  orbitals of monomer1/monomer2/dimer:
#                         1206 1206  2412
# Calling Multiwfn to obtain amount of orbitals and electrons of monomer1 ...
# Reading coefficient matrix of monomer1 ...
# Calling Multiwfn to obtain amount of orbitals and electrons of monomer2 ...
# Reading coefficient matrix of monomer2 ...
# Calling Multiwfn to obtain amount of orbitals and electrons of dimer ..
# Reading overlap matrix of dimer ...
# Reading coefficient matrix of dimer ...
# Inversing coefficient matix of dimer ...
# Calculating Fock matrix ...

作者
Author: pizaipizai@5989 时间: 2021-11-3 22:10
撸主，我拿相同体系放到超算上跑了一遍，也报了个什么内存引用错误。留着给你debug用吧

# Using "dimer.molden.input" for dimer.
# Using "mono1.molden.input" for monomer1.
# Using "mono2.molden.input" for monomer2.
# Calling Multiwfn to obtain amount of orbitals and electrons of monomer1 ...
# Reading amount of orbitals and electrons of monomer1 ...
# Calling Multiwfn to obtain amount of orbitals and electrons of monomer2 ...
# Reading amount of orbitals and electrons of monomer2 ...
# Calling Multiwfn to obtain amount of orbitals and electrons of dimer ...
# Reading amount of orbitals and electrons of dimer ...
# number of electrons of monomer1/monomer2/dimer:
#                            126    126 252
# number of  orbitals of monomer1/monomer2/dimer:
#                         1206 1206  2412
# Calling Multiwfn to obtain amount of orbitals and electrons of monomer1 ...
# Reading coefficient matrix of monomer1 ...
# Calling Multiwfn to obtain amount of orbitals and electrons of monomer2 ...
# Reading coefficient matrix of monomer2 ...
# Calling Multiwfn to obtain amount of orbitals and electrons of dimer ...
# Reading overlap matrix of dimer ...
# Reading coefficient matrix of dimer ...
# Inversing coefficient matix of dimer ...

Program received signal SIGSEGV: Segmentation fault - invalid memory reference.

Backtrace for this error:
#0  0x4c3cba
#1  0x4c3653
#2  0x51a05f
#3  0x40e1f5
#4  0x40c8d5
#5  0x40bd8a
#6  0x402ebf
#7  0x40ba93
#8  0x40089e
#9  0x50c308
#10  0x4013e9
Segmentation fault (core dumped)

作者
Author: snljty 时间: 2021-11-4 13:47

pizaipizai@5989 发表于 2021-11-3 22:10
撸主，我拿相同体系放到超算上跑了一遍，也报了个什么内存引用错误。留着给你debug用吧
# Using "dimer ...

麻烦上传一下文件。另外是否用unlimit修改了程序堆内存的使用上限？

作者
Author: Larry2022 时间: 2022-2-25 20:14
本帖最后由 Larry2022 于 2022-2-25 20:18 编辑

你好，非常感谢您的分享，对我帮助很大，此外我想请问一下您这个程序使用的什么方法来计算电荷转移积分的呢？最后，在发表论文时该怎么引用您的程序呢？谢谢

作者
Author: snljty 时间: 2022-2-25 20:44

Larry2022 发表于 2022-2-25 20:14
你好，非常感谢您的分享，对我帮助很大，此外我想请问一下您这个程序使用的什么方法来计算电荷转移积分的呢 ...

格点能修正法(前面正文里有提到)，这个方法的文献可以找一下然后引用，Multiwfn建议引用，目前这个程序不用引用。

作者
Author: Larry2022 时间: 2022-2-25 20:59

snljty 发表于 2022-2-25 20:44
格点能修正法(前面正文里有提到)，这个方法的文献可以找一下然后引用，Multiwfn建议引用，目前这个程序不 ...

好的

作者
Author: kaiden 时间: 2022-4-20 22:29
本帖最后由 kaiden 于 2022-4-20 22:42 编辑

楼主您好，我最近尝试使用您的程序算电荷转移积分，但是过程遇到些问题，请教下面报错会是什么问题呢？谢谢！（window版程序，已经正确设置Multiwfn环境变量，用楼主的sample测试正常。我的体系中不加-d47 时显示系数矩阵不可逆，加了-d47时出错如下）# Using "dimer.fchk" for dimer.
# Using "monomer1.fchk" for monomer1.
# Using "monomer2.fchk" for monomer2.
# Using "MYDIMER.47" for Fock matrix of dimer.
# Calling Multiwfn to obtain amount of orbitals and electrons of monomer1 ...
# Reading amount of orbitals and electrons of monomer1 ...
# Calling Multiwfn to obtain amount of orbitals and electrons of monomer2 ...
# Reading amount of orbitals and electrons of monomer2 ...
# Calling Multiwfn to obtain amount of orbitals and electrons of dimer ...
# Reading amount of orbitals and electrons of dimer ...
# number of electrons of monomer1/monomer2/dimer:
#                            406    406 812
# number of  orbitals of monomer1/monomer2/dimer:
#                         1280 1280  2560
# Calling Multiwfn to obtain amount of orbitals and electrons of monomer1 ...
# Reading coefficient matrix of monomer1 ...
# Calling Multiwfn to obtain amount of orbitals and electrons of monomer2 ...
# Reading coefficient matrix of monomer2 ...
# Calling Multiwfn to obtain amount of orbitals and electrons of dimer ...
# Reading overlap matrix of dimer ...
# Reading Fock matrix of dimer from an NBO .47 file ...
At line 834 of file calc_coupling.f90
Fortran runtime error: Bad real number in item 1 of list input

Error termination. Backtrace:

Could not print backtrace: libbacktrace could not find executable to open
#0  0xffffffff
#1  0xffffffff
#2  0xffffffff
#3  0xffffffff
#4  0xffffffff
#5  0xffffffff
#6  0xffffffff
#7  0xffffffff
#8  0xffffffff
#9  0xffffffff
#10  0xffffffff
#11  0xffffffff
#12  0xffffffff
#13  0xffffffff
#14  0xffffffff
#15  0xffffffff

作者
Author: snljty 时间: 2022-4-21 00:07

kaiden 发表于 2022-4-20 22:29
楼主您好，我最近尝试使用您的程序算电荷转移积分，但是过程遇到些问题，请教下面报错会是什么问题呢？谢谢 ...

麻烦把文件发给我，我来测一下。

作者
Author: kaiden 时间: 2022-4-21 20:59

snljty 发表于 2022-4-21 00:07
麻烦把文件发给我，我来测一下。

压缩文件有点大，我用wetransfer传给您。下载链接已经私信发给您了，麻烦啦！

作者
Author: milanox 时间: 2022-6-29 15:02
你好，我在使用您给的测试文件计算时。出现图中的报错，显示 “Cannot open file 'Multiwfn_in.txt'”。但是我在multiwfn文件夹下也没有找到这个文件，想请教下需要怎么处理这种情况。

作者
Author: snljty2 时间: 2022-6-29 15:28

milanox 发表于 2022-6-29 15:02
你好，我在使用您给的测试文件计算时。出现图中的报错，显示 “Cannot open file 'Multiwfn_in.txt'”。 ...

在什么目录下运行的？看起来是没有权限，换个普通点的目录试试。

作者
Author: 未名|梧桐雨 时间: 2022-9-11 14:06
大佬你好，我用了您发的Multiwfn计算转移积分的程序，发现计算我的4个结构单元的寡聚物体系不收敛出现报错，提示（4700，4700）不收敛，请问下是超过4700个轨道的分子就无法计算嘛？谢谢~

作者
Author: 未名|梧桐雨 时间: 2022-9-11 14:18
报错如图所示，麻烦大佬了！

作者
Author: seahurt 时间: 2022-9-12 13:37
用Fortran的原因除了性能还有其它吗？

作者
Author: 兰兰爱学习 时间: 2022-9-14 15:27
楼主想请教一下这样读取数值对吗？

作者
Author: 蛤蛤蛤蛤蛤蛤 时间: 2022-10-17 17:56
老师您好，我计算c12h26 c14h28等几种有机小分子电荷迁移模型，在cp2k中已经进行了结构优化用了，pbe泛函和tzv2p基组，从优化了的结构中我选择了几个二聚体，用cp2k molden文件计算二聚体转移积分，已经添加了盒子信息和有效电荷，出现overlap matrix 可能在哪里有问题呢

作者
Author: 蛤蛤蛤蛤蛤蛤 时间: 2022-10-17 18:15
这个是报错信息

作者
Author: 蛤蛤蛤蛤蛤蛤 时间: 2022-10-17 21:59
不好意思已经解决了，感谢

作者
Author: 蛤蛤蛤蛤蛤蛤 时间: 2022-10-17 22:02

未名|梧桐雨发表于 2022-9-11 14:18
报错如图所示，麻烦大佬了！

这个我也遇到了解决了，你控制打印让他打印出的和基函数一致就行了

作者
Author: snljty2 时间: 2022-10-17 22:10

兰兰爱学习发表于 2022-9-14 15:27
楼主想请教一下这样读取数值对吗？

读J_eff，后两组数据平时一般不用。读第一行和最后一行没问题。

作者
Author: snljty2 时间: 2022-10-17 22:13

seahurt 发表于 2022-9-12 13:37
用Fortran的原因除了性能还有其它吗？

读Gaussian里用Fortran输出的科学记数法双精度浮点数比较省事，就是那个带D的格式，因为有时候甚至不含那个D，比如指数超过两位，1.234567E+102会被输出成1.234567+102而不是1.234567D+102（这个虽然可以写程序的时候控制），别的语言读这个太麻烦了。另外就是Fortran调LAPACK比C方便一些，数组运算也方便很多，不用自己写一堆辅助函数。大概是这些原因。

作者
Author: 兰兰爱学习 时间: 2022-10-18 08:22

snljty2 发表于 2022-10-17 22:10
读J_eff，后两组数据平时一般不用。读第一行和最后一行没问题。

好的谢谢老师，第一行是空穴转移积分最后一行是电子转移积分这样对应是对的吧老师？

作者
Author: Windln 时间: 2022-12-17 12:35
“注意monome1r的坐标应该是dimer坐标的前面部分，monomer2的坐标应该是dimer坐标的后半部分，两者坐标要和其在dimer中的坐标完全相同，而且顺序不可改变”
在dimer中坐标不变，将monomer2即后面的坐标作为1，monomer1作为2，这样的操作是可以接受的吗？

作者
Author: Windln 时间: 2022-12-17 12:41

Windln 发表于 2022-12-17 12:35
“注意monome1r的坐标应该是dimer坐标的前面部分，monomer2的坐标应该是dimer坐标的后半部分，两者坐标要和 ...

还想问一个问题，单体的信息一定要在单体单独存在下重新优化得到，还是从二聚体中扣出单体结构只算单点能也行呢？

作者
Author: snljty 时间: 2022-12-21 09:39

Windln 发表于 2022-12-17 12:41
还想问一个问题，单体的信息一定要在单体单独存在下重新优化得到，还是从二聚体中扣出单体结构只算单点能 ...

不能重新优化。

作者
Author: 楠楠 时间: 2023-7-4 08:52
本帖最后由楠楠于 2023-7-10 08:59 编辑

老师好想请教一下，我先优化的dimer 然后基于优化好的dimer 保存了两个单体 dimer加上两个单体我都加上了关键词 NoSymmetry IOp(3/32=2,3/33=1) 进行单体计算
然后输入程序时出现如下报错请问老师是什么原因

作者
Author: 王OldMan 时间: 2024-11-10 22:14

snljty 发表于 2021-9-13 11:50
本文一楼已经写明白了全部需要的操作指南，要注意的容易忽视的点还特意用红字标出来了，还贴上了接近1000 ...

老师，您好。我也想查看这软件更详细的使用说明，可是我打开这个程序输入“help”后回车，程序一直弹出“File “help” not found！” 是我查看的方式不对吗？该如何查看使用说明？本人比较小白，请老师不吝赐教

欢迎光临计算化学公社 (http://bbs.keinsci.com/)