计算化学公社

标题: cp2k跑分子动力学应该怎么分析结果？ [打印本页]

作者
Author: zhangs 时间: 2024-9-24 11:00
标题: cp2k跑分子动力学应该怎么分析结果？
各位老师，我使用cp2k跑了一个分子动力学，不知道怎么分析结果，目前只能分析距离随时间变化，我想得到一些关于电荷/电子转移还有成键的相关信息，应该如何做？

作者
Author: sobereva 时间: 2024-9-24 21:42
结合Multiwfn做波函数分析，有的是方法
Multiwfn支持的分析化学键的方法一览
http://sobereva.com/471（http://bbs.keinsci.com/thread-12488-1-1.html）
谈谈怎么考察、计算、分析化学体系的电子密度
http://sobereva.com/715（http://bbs.keinsci.com/thread-46380-1-1.html）

可以每隔一定帧数导出一个molden文件，按下文批量用multiwfn进行分析
详谈Multiwfn的命令行方式运行和批量运行的方法
http://sobereva.com/612（http://bbs.keinsci.com/thread-24929-1-1.html）
上面的文章里就有现成的分析ORCA做AIMD过程中电子结构变化的例子

作者
Author: ljh123 时间: 2024-11-13 13:57

sobereva 发表于 2024-9-24 21:42
结合Multiwfn做波函数分析，有的是方法
Multiwfn支持的分析化学键的方法一览
http://sobereva.com/471（h ...

sob老师请问一下，我已经用CP2K跑完了AIMD，现在应该已经不能再每隔一定帧数导出一个molden文件了。那我该怎么获得multiwfn需要的输入文件呢？是把xyz轨迹文件里的每一帧坐标都拆成单一文件计算单点能获得波函数吗？这样的话计算量似乎有点高

作者
Author: sobereva 时间: 2024-11-13 23:40

ljh123 发表于 2024-11-13 13:57
sob老师请问一下，我已经用CP2K跑完了AIMD，现在应该已经不能再每隔一定帧数导出一个molden文件了。那我 ...

是
没必要每步都产生，诸如可以每10帧产生一次

作者
Author: ljh123 时间: 2024-11-14 13:56

sobereva 发表于 2024-11-13 23:40
是
没必要每步都产生，诸如可以每10帧产生一次

非常感谢sob老师非常细致认真的回答，总能解决我很多疑问。

作者
Author: Uus/pMeC6H4-/キ 时间: 2024-11-14 15:43
本帖最后由 Uus/pMeC6H4-/キ于 2024-12-15 15:24 编辑

sobereva 发表于 2024-11-13 23:40
是
没必要每步都产生，诸如可以每10帧产生一次

社长您好，我也在考虑利用之前您提到的REFTRAJ字段来这样重算部分MD轨迹，感觉似乎可以在重算中提高计算级别。比如一个KCl+100H2O构成的水溶液体系（编辑2：应实验室要求确认此为合理的简化模型），我打算先在便宜的PBE-D3(BJ)/DZVP-MOLOPT-SR-GTH级别下用800 Ry截断能跑分子动力学，每5帧或10帧记录一个坐标的snapshot，然后在更贵的PBE0-D3(BJ)/TZVP-MOLOPT-PBE0-GTH级别下用稍高一些如1000 Ry截断能结合REFTRAJ重新计算每个snapshot的单点能和原子受力。坐标、单点能和原子受力是用于训练NNP的。

虽然高计算级别意味着时间加倍，但是保留用于高级别单点能计算的snapshot数量只有原先分子动力学的几分之一，如果估计同样计算资源（集群单个节点112核，单任务计算时长上限一周）单点能计算的总耗时和之前分子动力学模拟的总耗时相当的话还是可以接受的。请指点一下这个设计是否合理。

作者
Author: sobereva 时间: 2024-11-15 04:56

Uus/pMeC6H4-/キ发表于 2024-11-14 15:43
社长您好，我也在考虑利用之前您提到的REFTRAJ字段来这样重算部分MD轨迹，感觉似乎可以在重算中提高计算 ...

很合理
由于每一帧单独算单点+受力的时候没法像MD那样读取上一步收敛的波函数当初猜，所以应注意检查批量算单点+受力过程的SCF收敛情况

作者
Author: Uus/pMeC6H4-/キ 时间: 2024-11-15 11:41

sobereva 发表于 2024-11-15 04:56
很合理
由于每一帧单独算单点+受力的时候没法像MD那样读取上一步收敛的波函数当初猜，所以应注意检查批 ...

波函数收敛的确很重要呢。

其实从已有的水溶液研究文献J. Phys. Chem. Lett. 2017, 8, 7, 1545–1551和Nature Chemistry 16, 1838–1844 (2024)来看，在CP2K中利用ADMM技术时计算级别甚至可达revPBE0-D3(BJ)方法结合TZV2P档次基组。不过还有几个小细节要确认一下是否正确：

1.我看到博文《CP2K做杂化泛函计算的关键要点和简单例子》写了PBE0泛函的定义方法，用了手册所述&PBE字段书写。我理解的revPBE0和revPBE的关系也类似于PBE0与PBE，revPBE0由75%的revPBE交换项+25%的HF交换项+100%的revPBE相关项（实与原PBE一致）构成，那么书写时只要加一行PARAMETRIZATION REVPBE就行。

2.在CP2K的data/BASIS_MOLOPT_UZH中有专为PBE0优化的TZV2P-MOLOPT-PBE0-GTH系列基组，但是用revPBE0时应该用回data/BASIS_MOLOPT_UCL里的TZV2P-MOLOPT(-SR)-GTH。（毕竟BASIS_MOLOPT_UZH没写原始文献）

3.上述JPCL论文采用的ADMM辅助基组是cpFIT3，但如今最好改成上述博文推荐的admm-tzp。

作者
Author: sobereva 时间: 2024-11-15 11:48

Uus/pMeC6H4-/キ发表于 2024-11-15 11:41
波函数收敛的确很重要呢。

其实从已有的水溶液研究文献J. Phys. Chem. Lett. 2017, 8, 7, 1545–1551 ...

1 或许正确，但我没具体了解过revPBE0。若无特殊必要用，还是建议用流行得多的PBE0

2 只要是杂化泛函，最好都用给PBE0优化的基组，原理上比用对纯泛函优化的基组更理想

3 对

作者
Author: Uus/pMeC6H4-/キ 时间: 2024-11-15 13:23

sobereva 发表于 2024-11-15 11:48
1 或许正确，但我没具体了解过revPBE0。若无特殊必要用，还是建议用流行得多的PBE0

2 只要是杂化泛函 ...

用revPBE-D3做液态水体系AIMD的文章能看到好几篇（没准今年第一性原理泛函流行程度投票里面有因此投票给这个泛函的），感觉可能跟J. Phys. Chem. C 2014, 118, 50, 29401–29411比较若干泛函的结论有关：

It is found that, among these seven functionals, BLYP, BLYP-D2, and PBE-D3 functionals produce slower dynamics of spectral diffusion, while BLYP-D3, revPBE-D3, DRSLL-PBE, and DRSLL-optB88 functionals produce faster dynamics. The latter results are comparable to each other, and they are in good agreement with available experimental results.

想想也是，如果其他类型的体系也要低级别跑完MD后高级别重算单点能和受力，通用性更强的PBE0会更合适。

作者
Author: ljh123 时间: 2024-12-13 20:52

Uus/pMeC6H4-/キ发表于 2024-11-14 15:43
社长您好，我也在考虑利用之前您提到的REFTRAJ字段来这样重算部分MD轨迹，感觉似乎可以在重算中提高计算 ...

你好，我在用REFTRAJ重算能量的时候遇到一些问题。因为我是用NPT跑的轨迹，因此每一帧的盒子尺寸都在波动，而使用REFTRAJ是似乎只能用一个盒子尺寸，轨迹只能用xyz格式，这样重跑的话对每一帧的能量似乎有一定影响。请问这种情况你遇到过吗，是不是只能单独提取结构算单点能了呢

作者
Author: Uus/pMeC6H4-/キ 时间: 2024-12-13 21:02

ljh123 发表于 2024-12-13 20:52
你好，我在用REFTRAJ重算能量的时候遇到一些问题。因为我是用NPT跑的轨迹，因此每一帧的盒子尺寸都在波动 ...

我自己没做过这种情况的重跑，但手册说可以添加CELL_FILE_NAME字段指定记录盒子信息的.cell文件来重跑，应该用先前&MOTION/&PRINT/&CELL控制输出的文件就行。当然要保证先前&MOTION/&PRINT/&CELL/EACH和&MOTION/&PRINT/&TRAJECTORY/EACH的记录频率一致，盒子和轨迹的记录能对上。

在&FORCE_EVAL/&SUBSYS那里也需要定义一个盒子和结构，我感觉用重跑第一帧的就行。

作者
Author: ljh123 时间: 2024-12-13 21:13

Uus/pMeC6H4-/キ发表于 2024-12-13 21:02
我自己没做过这种情况的重跑，但手册说可以添加CELL_FILE_NAME字段指定记录盒子信息的.cell文件来重跑， ...

我尝试过CELL_FILE_NAME字段，结果他的盒子信息默认用.cell文件最后一帧的盒子尺寸，还是无法实现NPT的重算

作者
Author: Uus/pMeC6H4-/キ 时间: 2024-12-13 21:30

ljh123 发表于 2024-12-13 21:13
我尝试过CELL_FILE_NAME字段，结果他的盒子信息默认用.cell文件最后一帧的盒子尺寸，还是无法实现NPT的重 ...

读取盒子信息还需要把默认为F的VARIABLE_VOLUME设成T哦。

作者
Author: ljh123 时间: 2024-12-13 21:34

Uus/pMeC6H4-/キ发表于 2024-12-13 21:30
读取盒子信息还需要把默认为F的VARIABLE_VOLUME设成T哦。

原来如此，十分感谢！还有就是想请教下，通过REFTRAJ如何输出动能呢？我用REFTRAJ得到的.ener文件里只有势能，其他全为0

作者
Author: Uus/pMeC6H4-/キ 时间: 2024-12-13 21:59
本帖最后由 Uus/pMeC6H4-/キ于 2024-12-13 22:03 编辑

ljh123 发表于 2024-12-13 21:34
原来如此，十分感谢！还有就是想请教下，通过REFTRAJ如何输出动能呢？我用REFTRAJ得到的.ener文件里只有 ...

轨迹重算应该是读取各结构计算单点能（电子的势能）和原子受力的，又不涉及按分子动力学走步算法随时间推移原子变化结构的计算，那自然没有（原子核）动能和温度的输出。

编辑：动能这个问题似乎在http://bbs.keinsci.com/thread-50101-1-1.html和http://bbs.keinsci.com/thread-50466-1-1.html这两帖都讨论过，但是我也不太清楚分子动力学逐帧“自由能”或其热平均具体怎么处理。

作者
Author: ljh123 时间: 2024-12-13 22:03

Uus/pMeC6H4-/キ发表于 2024-12-13 21:59
轨迹重算应该是读取各结构计算单点能（电子的势能）和原子受力的，又不涉及按分子动力学走步算法随时间推 ...

感谢你的解答。我只是好奇REFTRAJ为什么不输出动能，倒不是需要算他。我有大概的猜测，但是没人讨论就有点没底hhh

作者
Author: ljh123 时间: 2024-12-13 22:43

Uus/pMeC6H4-/キ发表于 2024-12-13 21:59
轨迹重算应该是读取各结构计算单点能（电子的势能）和原子受力的，又不涉及按分子动力学走步算法随时间推 ...

你好，我尝试了VARIABLE_VOLUME关键词，但一直出现报错，不知是否是我的关键词输入问题。输入文件部分如下：
&REFTRAJ
   CELL_FILE_NAME 1.cell
   VARIABLE_VOLUME .TRUE.
   TRAJ_FILE_NAME 1.xyz
   EVAL ENERGY_FORCES
   FIRST_SNAPSHOT 1
   LAST_SNAPSHOT 5000
   STRIDE 50
&END REFTRAJ
报错信息如下：
* ___                                                                      *
*  / \                                                                   *
* [ABORT]                                                                   *
*  \___/                            CPASSERT failed                         *
* |                                                                      *
*  O/|                                                                      *
* /| |                                                                      *
* / \                                              motion/integrator.F:1622 *
似乎是motion这节的问题，我把VARIABLE_VOLUME关键词删掉是可以正常运行的，但是加上就会报错。。。

作者
Author: Uus/pMeC6H4-/キ 时间: 2024-12-13 23:10
本帖最后由 Uus/pMeC6H4-/キ于 2024-12-13 23:12 编辑

ljh123 发表于 2024-12-13 22:43
你好，我尝试了VARIABLE_VOLUME关键词，但一直出现报错，不知是否是我的关键词输入问题。输入文件部分如 ...

源代码src/motion/integrator.F第1622行附近是这么写的：

! Check if we reached the end of the file and provide some info..
IF (my_end) THEN
IF (reftraj_env%isnap /= (simpar%nsteps - 1)) &
CALL cp_abort(__LOCATION__, &
"Reached the end of the Trajectory frames in the TRAJECTORY file. Number of "// &
"missing frames ("//cp_to_string((simpar%nsteps - 1) - reftraj_env%isnap)//").")
END IF

复制代码

我猜是&MD/&REFTRAJ里设置的重算起始帧、终止帧、跨步与&MD里STEPS设置的步数有所冲突了，但是我不知道这和VARIABLE_VOLUME有啥关系……

编辑：我自己设置跨步为1跑重算时发现重算设置确定的步数会覆盖掉STEPS确定的步数，干脆去掉了STEPS设定（此时用回默认值3），结合FIRST_SNAPSHOT 1和LAST_SNAPSHOT <总结构数>可以顺利重算完整个轨迹文件的所有结构。

作者
Author: Uus/pMeC6H4-/キ 时间: 2025-1-4 12:11

sobereva 发表于 2024-11-15 04:56
很合理
由于每一帧单独算单点+受力的时候没法像MD那样读取上一步收敛的波函数当初猜，所以应注意检查批 ...

悲报，SCF不收敛的问题果然极其高发，多个大小不一的平行体系都因此报错终止了……

对若干分子数和盒子尺寸不同的水溶液，按6楼的计算方法/基组做分子动力学和REFTRAJ重算，权衡精度和时长上限后截断能都用800 Ry。观察低级别纯泛函的MD轨迹没有明显结构问题，高级别杂化泛函重算时开头一小段轨迹也能顺利SCF收敛算出能量和受力，但某一步SCF就会突然开始震荡、严重偏离收敛限。

此时是否需要取.restart文件记录的最后一帧用纯泛函重新算单点产生波函数，给杂化泛函读取来继续重算呢？有没有杂化泛函重算一开始就能避免这种问题的设置（开IGNORE_CONVERGENCE_FAILURE不算）？感谢您的指点。

（我也知道有条件好的敢直接上杂化泛函跑MD，但我这计算资源支持纯泛函跑MD加杂化泛函重算的权宜之计都勉强，叹气）

作者
Author: sobereva 时间: 2025-1-5 01:28

Uus/pMeC6H4-/キ发表于 2025-1-4 12:11
悲报，SCF不收敛的问题果然极其高发，多个大小不一的平行体系都因此报错终止了……

对若干分子数和盒 ...

需要
没有。为减小中途SCF不收敛的可能，可以减小步长

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