|
本帖最后由 ggdh 于 2022-3-22 18:54 编辑
类似的部分程序(待补充):
Turbomole: 支持ri加速的scs-cc2, 但是要收费
eT: 有multilevel-cc加速, 但不支持三重态计算
orca: ADC不支持计算三重态, DLPNO-STEOM-CCSD 速度很快, 但容易出现收敛问题
MRCC: 待测试
1. 程序安装
1.1 安装anaconda或者miniconda
1.2 安装pyscf和adcc
- conda create -n pyscf -c pyscf pyscf
- conda activate pyscf
- conda install -c conda-forge adcc
复制代码 另外也可以选择别的方法比如pip, 参考:
https://adc-connect.org/v0.15.13/installation.html
https://pyscf.org/install.html
2. 输入文件的准备
2.1 准备一个优化好的xyz文件, 手边没有现成的可以用附件中的xyz文件, 其他格式的结构文件可以用Multiwfn (主功能100里的子功能2)转成xyz
结构图所示,该分子具有一定的MR-TADF结构特征
2.2 创建pyscf的输入文件BN.py, 其中内容如下
- from pyscf import gto, scf, lib
- import adcc
- lib.num_threads(12)
- # Run SCF in pyscf
- mol = gto.M(atom=gto.mole.fromfile("BN.xyz"), basis='cc-pvdz')
- scfres = scf.RHF(mol)
- scfres.kernel()
- adcc.set_n_threads(12)
- state = adcc.adc2(scfres, n_triplets=1)
- print(state.describe())
- print(state.describe_amplitudes())
复制代码 其中需要注意的是:
- 12为openmp并行核数, 这里并不支持mpi并行, 要分别为pyscf和adcc设
- BN.xyz为xyz文件名
- cc-pvdz为基组名
- n_triplets=1的意思是算一个三重态, 如果要算单重态改成n_singlets=5就是算5个单重态
2.3运行程序
- export PYSCF_TMPDIR=/tmp
- conda activate pyscf
- python BN.py | tee BN.out
复制代码 其中第一行是设缓存目录, 大体系很耗硬盘, 所以如果这里尽量设一个体积大的硬盘位置
然后部分输出结果如下:
- converged SCF energy = -539.461545133644
- Starting adc2 triplet Jacobi-Davidson ...
- Niter n_ss max_residual time Ritz values
- 1 4 0.32167 23.3s [0.29368883]
- 2 8 0.039105 12.7s [0.13900202]
- 3 12 0.0017295 21.4s [0.1199706]
- 4 16 0.00093827 24.2s [0.11816422]
- 5 20 0.00016175 24.8s [0.11744624]
- === Converged ===
- Number of matrix applies: 20
- Total solver time: 106s 519ms
- +------------------------------------------------------------------------+
- | adc2 triplet , converged |
- +------------------------------------------------------------------------+
- | # excitation energy osc str |v1|^2 |v2|^2 |
- | (au) (eV) |
- | 0 0.1174462 3.195875 0.0000 0.9119 0.08811 |
- +------------------------------------------------------------------------+
复制代码 可以看到其中有adc2的收敛情况, 激发态能级, 振子强度, 单激发态和双激发的比例
3. 加速的技巧
加速通常意味着损失精度或者有产生误差风险, 但是postHF往往会面对体系太大无法计算的情况(用这个程序必须保证内存必须要大)
所以往往不得不放弃部分精度来完成计算. 这其中涉及到的道德,IF,良心,电费等因素还需用户自行权衡.
首先是可以改变的是基组的大小和体系中非共轭部分的裁剪, 这就不多说了.
然后影响速度的参数主要是如下几个参数:
state = adcc.adc2(scfres, n_triplets=1, frozen_core=0,frozen_virtual=0, n_guesses=4, conv_tol=1e-6)
其中:
n_triplets/n_singlets: 所需计算的态的数量
frozen_core: 冻结内层轨道的数量, 从低往高数
frozen_virtual: 冻结空轨道的数量, 从高往低数
n_guesses: 计算的猜的数量, 默认是激发态的数量的2倍, 但不会小于4
conv_tol: 收敛标准, 默认是SCF收敛标准的10倍但不会小于1e-6, 上面的例子中这个值大概是在1e-4
分别对这几个参数做了测试, 耗时和精度下降的情况入下图所示
其中time列是耗时(秒), eV列是激发态跃迁能, ΔeV列代表精度的损耗,
假设道德底线是0.01 eV, 假设其他体系也有类似的规律, 那么从图中可得得结论如下
收敛标准conv_tol对耗时影响较大, 但对精度影响不大1e-4(默认)都能接受,保险的话设1e-5(如果要提高到1e-5, 注意必须得另外把scf的精度调高到1e-6)
猜数量n_guesses对耗时影响大, 但对精度影响不大(这个应该要看体系, 改小这个参数一定要小心)
冻结内层电子frozen_core对耗时影响大, 如果只冻结内层电子则对精度影响不大(这个体系的内层轨道--1S轨道--正好是14个)
冻结高层空轨道frozen_virtual对耗时影响没那么大, 但由于空轨道数量多, 总体对速度影响还是大, 但是对精度影响较大(0.01eV级别)
综上所述,conv_tol用默认就够, n_guesses可以不放心的稍微改小, frozen_core可以比较放心的冻结内层轨道, frozen_virtual调大则需要耗费道德值.
4. 自动设加速参数的输入文件
根据第三节的测试, 我们发现conv_tol和n_guesses比较好设, 但是frozen_core和frozen_virtual比较让人纠结
所以我写了根据占据轨道和空轨道能量的高低来自动设frozen_core和frozen_virtual的输入文件如下
- from pyscf import gto, scf, lib
- import adcc
- occ_threash=3 # orbitals 3 Ha lower than HOMO will be frozen
- virtual_threash=5 # orbitals 3 Ha higher than LUMO will be frozen
- lib.num_threads(12)
- # Run SCF in pyscf
- mol = gto.M(atom=gto.mole.fromfile("BN.xyz"), basis='cc-pvdz')
- scfres = scf.RHF(mol)
- scfres.kernel()
- for i,o in enumerate(scfres.mo_occ):
- if o == 0:
- LUMO_idx = i
- HOMO_idx = i-1
- break
- occ_threash = scfres.mo_energy[HOMO_idx] - occ_threash
- virtual_threash = scfres.mo_energy[LUMO_idx] + virtual_threash
- fo = len([i for i in scfres.mo_energy if i < occ_threash])
- fv = len([i for i in scfres.mo_energy if i > virtual_threash])
- print("{:d} occ orbitals lower than {:.3f} eV will be frozen".format(fo,occ_threash*27.212))
- print("{:d} virtual orbitals higher than {:.3f} eV will be frozen".format(fv,virtual_threash*27.212))
- # Run adcc
- adcc.set_n_threads(12)
- state = adcc.adc2(scfres, n_singlets=1, n_guesses=2, conv_tol=1e-4, frozen_core=fo, frozen_virtual=fv)
- print(state.describe())
- print(state.describe_amplitudes())
复制代码 大家可以自行修改文件开头的occ_threash和virtual_threash值来设置根据能量来固定多少轨道, 如果道德水平在0.01eV, 我估计可以设4,4 (比HOMO低4Ha 比LUMO高4Ha的轨道纳入冻结), 如果道德水平在0.1eV,我估计可以设3,1. 大家可以弄个小体系自行尝试.
5. 一个小的benchmark
一些讨论:- adc2得到的结果和用orca的STEOM-DLPNO-CCSD得到的结果近似
- 纯泛函和HF%小的泛函(B3LYP)会因为低估T1而产生大的ΔEST
- 高HF%的杂化泛函/范围分离泛函/双杂化泛函对T1的估算比较准, 但是会因为高估S1而产生大的ΔEST
- 双杂化泛函表现并不好, 但部分泛函加了SOS后情况得到了好转其中表现最好的是SOS-PBE-QIDH和SOS-RSX-QIDH
- 注意上述结论仅仅针对当前一个分子得到, 对其他体系是否适用还需大家自行测试
附件中包含了excel数据,以及分子结构xyz文件
|
-
-
BN.xyz
1.25 KB, 下载次数 Times of downloads: 9
-
-
mrtadf.xlsx
27.07 KB, 下载次数 Times of downloads: 4
评分 Rate
-
查看全部评分 View all ratings
|