基于耦合簇与双杂化密度泛函方法的两点基组外推法等价性及稳健参数化研究

时间:2026年5月28日
来源:The Journal of Physical Chemistry A

编辑推荐:

基组外推(Basis-set extrapolation, BSE)至完全基组(complete basis set, CBS)极限是多项高精度热化学计算方案的基石。本研究评估了使用BSE时常作的一些近似处理。虽然外推参数通常使用耦合簇方法拟合,但在某些情况下

广告
   X   

基组外推(Basis-set extrapolation, BSE)至完全基组(complete basis set, CBS)极限是多项高精度热化学计算方案的基石。本研究评估了使用BSE时常作的一些近似处理。虽然外推参数通常使用耦合簇方法拟合,但在某些情况下,这些相同参数被应用于其他计算方法。尽管拟合参数可能存在显著差异,但当应用于两个基准数据集时,耦合簇方法与方法特定拟合参数之间的性能差异很小。尽管存在多种双基组外推至CBS极限的方案,但研究表明三种常见的外推方案——指数型、指数平方根型与反幂次型——均可简化为一种Schwenke风格的校正,即对较大基组的能量,使用两个连续基组能量缩放差值的校正(ECBS = En + f·(En – En-1))。此外,这种简化意味着不同方案之间的差异在数学上是无意义的,并且给定一种方案的外推参数,即可找到相应的Schwenke型外推因子f及其他两种方案的外推因子,反之亦然。最后,本文提供了耦合簇和选定双杂化密度泛函交换-相关泛函针对特定基组选择的Schwenke型两点外推参数。
基组近似(即使用一组数学函数描述波函数)是计算化学的基石,尤其对于波函数理论(wavefunction theory, WFT)和密度泛函理论(density functional theory, DFT,亦称Kohn-Sham理论)方法。理想情况下,研究者希望在无限基组下工作,因为这样计算将是精确的(在给定的计算方法内)。然而,计算机(及其用户)往往无法处理无限基组,因此使用有限基组。这为计算引入了固有误差,至少可以通过使用越来越大的基组来最小化误差。但即便如此也有其限制,因为计算成本往往随基组大小极不友好地增长(例如,Hartree-Fock方法计算量随基组函数数N按N4缩放,而CCSD(T)方法按N7缩放)。还存在其他近似方法,如恒等式分辨(resolution of the identity, RI)或局域化方法(例如Orca中的DLPNO近似),也可用于降低计算成本,但这超出了本研究的范围。这正是基组外推(BSE)发挥作用之处。研究者观察到,性质(如总电子能)随基组增大往往平滑且单调地收敛。研究者可以利用这一行为,使用一组随基组增大的少量计算来外推至无限基组,通常称为完全基组(CBS)极限。这构成了多项高精度热化学方案(high-accuracy thermochemical schemes, HATS)的基础,如Weizmann-n、HEAT、Gaussian-n和CBS方法(近期综述见参考文献),这些方法通常基于耦合簇或Møller-Plesset方法。本文后续章节将讨论多种已提出的BSE外推方案及其基础方程。一个常见选择是使用指数平方根外推(En ≈ ECBS + A·exp(−α√L))处理Hartree-Fock能量,以及使用反幂次外推(En ≈ ECBS + B·n−β)处理相关能。Neese和Valeev针对多种基组组合优化了α和β值用于双基组外推(two-basis-set extrapolation, 2BSE),这些参数通常用于大多数应用。虽然确定2BSE参数的大多数过程涉及拟合基准能量(如Neese和Valeev的研究),但Fishman等人采用了新颖方法,要求基组叠加误差(basis set superposition error, BSSE)在CBS极限处消失。然而,WFT方法随基组大小的陡峭缩放性始终会对可计算体系的大小施加限制。由于DFT计算效率高,它是大型体系计算的支柱。特别是双杂化(double-hybrid, DH)泛函已被证明在仅需一小部分成本的情况下,其精度接近成本高得多的CCSD(T)方法,尤其在结合各种线性化技术(如RI-MP2、RIJK、RIJCOSX)时。Chai和Head-Gordon使用反幂次外推(β = 3)对ωB97X(2) DH泛函的PT2部分进行CBS外推。Chuang和Chen随后考虑了B2PLYP DH泛函各种参数化方案的多种外推方案;在他们的研究中,考虑了本文研究的所有三种外推方案,并针对每种泛函优化了参数。Karton和Martin使用B2GP-PLYP评估了基组外推,包括常规和显式关联(explicitly correlated)方法(即B2GP-PLYP-F12),但其β外推参数取自CCSD-F12外推(具体为Hill等人);Mehta和Martin随后使用GMTKN55基准数据集进行了更广泛的比较。Kraus评估了选定混合和DH泛函的性能,并确定了适当的外推参数(此外还重新评估了Hartree-Fock CBS外推)。另一个例子是Dohm等人在修订其作者MOBH35基准集部分条目时,使用了BSE外推的PWPB95能量,这些修订基于使用CCSD(T)方法的W1型CBS外推。尽管BSE应用广泛,但仍存在一些未解决问题,本研究将予以解答:Neese和Valeev指出,他们通过最小化21个小分子总能量的平均绝对偏差(mean absolute deviation, MAD)来拟合α和β 2BSE参数。如3.3节所述,有理由认为这可能并不理想,因为Hartree-Fock(为避免与分子氟化氢(测试集成员之一)混淆,该能量将称为SCF,即自洽场(self-consistent field)的首字母缩写)能的外推应独立于相关处理,并且SCF能量比相关能大2-3个数量级。Schwenke提出了一种基于两个连续基组计算能量缩放差值的2BSE,其中SCF能量使用指数平方根外推,CCSD(T)相关能使用反幂次外推。这有何优势?如3.2节所示,各种2BSE方法均可简化为Schwenke型外推,使得不同方案之间的差异意义减弱。多项HATS将为CCSD(T)相关能拟合的2BSE参数用于其他相关方法。这是否合理?为CCSD(T)拟合的2BSE参数也用于DH泛函。这是否是合理的方法?本研究使用了多项WFT从头算方法:Hartree-Fock(记为SCF);耦合簇单双激发及准微扰三重态处理(coupled-cluster with single and double excitations and a quasiperturbative treatment of the triples contribution, CCSD(T));二次组态相互作用单双激发及三重态校正(quadratic configuration interaction with single and double excitations and a triples correction, QCISD(T));二阶Møller-Plesset微扰理论(second-order Møller-Plesset perturbation theory, MP2)。所有WFT方法均使用限制开壳Hartree-Fock(restricted open-shell Hartree-Fock, ROHF)参考。研究使用了选定的DH DFT交换-相关泛函:Grimme原始DH泛函,基于Becke88(B88)交换泛函和Lee-Yang-Parr(LYP)相关泛函(B2PLYP);Goerigk和Grimme的DH泛函,基于Perdew-Wang91(PW91)交换泛函和Becke95(B95)相关泛函(PWPB95);Schwabe和Grimme的DH泛函,基于改进的Perdew-Wang(mPW91)交换泛函(其中mPW91泛函使用内置的LibXC库)和LYP相关泛函(mPW2PLYP);Kozuch和Martin原始实现的色散校正DH泛函(具体使用Grimme及其合作者的第三版经验色散校正(D3)并采用Becke-Johnson阻尼(D3BJ)),具有自旋分量标度(spin-component-scaled, SCS)MP2类相关能,使用Perdew-Burke-Ernzerhof(PBE)交换泛函和Perdew-86相关泛函(DSD-PBEP86);Santra等人修订版(revDSD-PBEP86-D3BJ)及其使用Grimme及其合作者第四版经验色散校正(D4)的版本(revDSD-PBEP86-D4);Karton等人对B2PLYP的通用重新参数化版本(B2GP-PLYP)。DFT计算对开壳体系使用非限制Kohn-Sham(unrestricted Kohn-Sham, UKS)形式,对闭壳体系使用限制Kohn-Sham(restricted Kohn-Sham, RKS)形式。Molpro的“高”网格精度设置用于所有DFT计算。DSD双杂化泛函的Molpro输入文件取自Jan M. L. (Gershom) Martin教授的主页。研究使用了四个系列的基组:Dunning的相关一致基组(correlation-consistent basis sets, cc-pVnZ和aug-cc-pVnZ, n = D, T, Q, 5, 6),默认包含(aug-)cc-pV(n+d)Z基组(该基组在第二行主族元素(即Al–Ar)上包含额外的d极化函数),对于Mg使用的aug-cc-pV5Z基组取自Iron、Oren和Martin的核-价基组但移除了核-价函数;Neese和Valeev的原子自然轨道(atomic natural orbital, ANO)基组(ANO-pVnZ和aug-ANO-pVnZ, n = D, T, Q, 5);Weigend和Ahlrichs的def2基组系列,包括不带(即def2-SVP、def2-TZVPP和def2-QZVPP)和带(即def2-SVPD、def2-TZVPPD和def2-QZVPPD)弥散函数的版本;Jensen的极化一致基组(polarization consistent basis sets, pc-n和aug-pc-n, n = 1–4)。ANO-pVnZ和aug-ANO-pVnZ基组从基组交换(Basis Set Exchange, BSE)网站下载。注意BSE中的ANO-pVTZ版本(撰写时)以及Neese和Valeev论文支持信息中的ANO-pVTZ实际是ANO-pVDZ的副本(这是BSE中已记录的错误);按照Neese在错误报告中的建议,ANO-pVTZ基组通过从ANO-pVQZ基组中移除每个角动量最外层函数获得(最终基组以选定格式提供于支持信息S7节,并作为独立文件附加于Orca、Gaussian和Molpro格式)。此外,所有可用来源中氢的ANO-pV5Z和aug-ANO-pV5Z基组与其四ζ(Qζ)对应基组相同,但照常使用。BSEF21、S66和MB43-16的几何结构分别取自各自出版物;BSEF21部分几何结构和BSEF74集中的新几何结构(见3.1节和支持信息S2节)在与原始研究相同的理论水平(即使用Gaussian16中实现的Becke三参数混合泛函(B3LYP)与def2-TZVP基组)下重新优化。MB43-16中一个成员的生成能为负值;为简化统计分析,将该反应方向翻转使所有能量为正。参数拟合使用Microsoft Excel的Solver模块,默认“GRG非线性”引擎。Maple 2024用于协助部分代数推导。对于CCSD(T),研究评估了拟合过程:CCSD(T)是精确热化学计算公认的主力方法。Neese和Valeev指出,他们通过最小化BSEF21数据集(或其更准确版本)总CBS能量(即Etot = ESCF + Ecorr)的MAD,为多种基组对的标准CCSD(T)计算优化了α和β。这些值常用于各种外推方法。需要解决的三个主要问题是:(i) 同时拟合SCF和相关能的α和β参数是否合理?(ii) 将CCSD和(T)相关能视为一个整体是否合理?(iii) 所用的外推方法(即指数平方根型 vs. 反幂次型 vs. 指数型)有何影响?当通过最小化Etot的MAD同时优化SCF和相关能参数时,SCF外推变得依赖于相关方法的选择(例如CCSD(T)、CCSDT、QCISD(T)、MP2或甚至无相关)。这是不合逻辑的,可能导致与体系物理不符的外推。此外,SCF参考能量与CCSD相关能之间,以及CCSD与(T)相关能之间存在2-3个数量级的差异(见表1)。SCF组分主导总能量,而相关贡献在不同分子间可能变化。同时优化参数将忽略每个组分的真实外推行为,而专注于获得最准确的总能量。这是不可取的。表1总结了使用三点指数外推(公式22)和aug-cc-pVnZ(n = Q, 5, 6)基组计算的BSEF21数据集中每个成员的CBS极限能量(SCF、总相关、CCSD相关和(T)相关,单位均为Ha)以及计算的相关能差异(定义为ΔECBScorr = ECBScorr – ECBSCCSD – ECBS(T),单位为mHa和cal/mol)。这些差异非常小——数量级为μHa(或cal/mol)。在高精度热化学应用中,可以为相关能的两个分量使用单一外推参数β,但在此类方法中,每一点差异都很重要,鉴于额外计算成本可忽略不计,建议使用独立的β参数。这在W1和W2等方法中尤为重要,这些方法使用比SCF和CCSD贡献小一号的基组来外推(T)贡献。从表1数据中,可以为2BSE拟合α、β和/或γ参数。对于如何拟合参数有两种思路。可以独立拟合每个组分的α或β至ECBS,或在单一拟合中同时优化两个参数以拟合总能量在完全基组极限的值(ECBStot = ECBSSCF + ECBSCCSD(T));在其研究中,Neese和Valeev指出他们选择了第二种方案。虽然后者由于拟合更灵活可能给出更低误差,但这可能通过引入误差补偿而引入非物理性(即SCF能量的主导地位允许各组分能量偏离最佳拟合以改善整体误差)。然而,这可能不是各能量组分的最佳拟合。如果绘制外推能量的MAD随αSCF和βCCSD(T)变化的网格图(针对{ T,Q }ζ外推,相对于表1中的ECBSSCF和ECBScorr),会得到一个合理的αSCF和βCCSD(T)值的低谷。低谷底部呈S形,每个臂有平台。全局最小值在αSCF = 1.581和βCCSD(T) = 16.733;后者与Helgaker等人最初提出的β = 3、Truhlar及其合作者提出的βMP2 = 1.91、βCCSD = 1.94、βCCSD(T) = 2.02,或Neese和Valeev的值β = 3.05都非常不同。实际上,全局最小值位于低谷的一个分支上,而这些其他值都位于另一个分支上。另一方面,如果独立考虑每个组分,则在αSCF = 5.771、βCCSD(T) = 3.339、βCCSD = 3.404和β(T) = 3.361处获得最小值。相比之下,使用相同aug-cc-pVnZ基组的{ T,Q }ζ外推,Neese和Valeev获得αSCF = 5.79和βCCSD(T) = 3.05。Ranasinghe和Petersson指出,使用Schwenke型形式(公式11)参数优化更平滑。类似的MAD随fSCF和fCCSD(T)变化的图现在显示出一个轮廓平滑且具有明显最小值的直线低谷。事实上,如果在网格上对每个fSCF点绘制使MAD最小的fCCSD(T),则会得到一条直线。此外,最小MAD轮廓平滑且有明显最小值。因此,存在一组明确的最优值,在扫描范围内全局最小值为fSCF = 1.899和fCCSD(T) = 0.007,对应αSCF ≈ 1.579和βCCSD(T) ≈ 17.272。独立绘制每个能量组分的MAD随f变化的图呈现V形曲线,全局最小值在fSCF = 0.271和fCCSD(T) = 0.603,对应αSCF ≈ 5.768和βCCSD(T) ≈ 3.399。同时拟合两个f参数再次向体系引入非物理性,这再次清楚地表明,建议独立考虑每个组分。从图3D(和图2D)的另一个观察是,(T)相关能对外推参数f(T)(或β(T))依赖性弱——MAD随外推参数变化的图要平坦得多。这解释了表1中当CCSD(T)相关能作为单一实体外推或两个组分独立外推时CBS值的微小差异(即ΔECBScorr)。结合非常接近的最小值,这也解释了CCSD(T)和CCSD曲线几乎完全重叠。因此,如上所述,原则上可以使用单一外推参数fCCSD(T),尽管使用独立外推参数不会增加任何额外计算成本。

关于双杂化泛函,已确定MP2通常不是精确热化学计算的首选方法,甚至被许多DFT泛函超越(参见GMTKN55基准集的WTMAD2值:MP2为6.91 kcal/mol,SCS-MP2为5.35 kcal/mol,而B97M-V、M06-2X、ωB97M-V、DSD-PBEP86和ωB97M(2)分别为6.37、4.79、3.29、3.10和2.19 kcal/mol)。另一方面,包含MP2类相关能组分的双杂化(DH)泛函已被证明其精度甚至接近耦合簇方法。在某些情况下,CCSD-T优化的参数被用于外推DH泛函的能量,例如Dohm等人修订作者MOBH35基准集的工作。因此,决定评估典型的DH泛函DSD-PBEP86在这些基准测试中的表现。技术说明:在基组外推过程中省略了D3BJ经验色散校正,因为它仅依赖于几何结构,与使用的基组无关。对于外推双杂化泛函至基组极限,可以考虑两种方法:分别外推Kohn-Sham(EKSDH)和MP2类(EMP2DH)组分,类似于标准MP2计算;或外推最终总能量(EtotDH)。先验来看,没有理由选择一种方案而非另一种。BSEF74*集的CBS能量通过与WFT方法相同的方式确定(即使用公式22和Dunning的aug-cc-pVnZ基组进行{ Q,5,6 }ζ外推),同时针对EtotDH及其两个组分EKSDH和EMP2DH。CBS能量列于表S12,相应的拟合两级(Schwenke)CBS外推参数及相关MAD列于表5。箱线图显示误差随基组增大而减小,{ Q,5 }ζ外推的所有误差均低于2 kcal/mol。同样,三点外推比相应的两点外推具有更大的误差分布。对于S66和MB16-43数据集,使用双杂化泛函相比MP2有显著改进。使用两个参数(表6,对于{ T,Q }ζ CBS外推ΔMAD < 1 kcal/mol)有非常微小的优势。鉴于DFT泛函(包括双杂化)的(半)经验性质,使用单一参数非常合理;然而,如前所述,使用两个参数将在零额外计算成本下提高一些精度。同样,使用DH优化的参数与使用CCSD-T优化参数的优势惊人地可以忽略(表7)。这些CBS拟合针对其他选定的DH泛函重复进行:B2PLYP、B2GP-PLYP、mPW2PLYP、PWPB95、revDSD-PBEP86-D3BJ和revDSD-PBEP86-D4。所得外推参数列于表5(aug-cc-pVnZ),其他基组系列的参数列于表S19-S25,而f及相关α、β和γ参数列于表S26。尽管基组外推是常规做法,但关于其使用仍存在一些未解决的问题,本研究试图予以解答。研究发现,同时为Hartree-Fock(SCF)和CCSD(T)相关能贡献优化参数会导致与独立优化每个参数时观察到的行为不一致的值。同样,对CCSD(T)相关能使用单一参数而非为每个组分(即CCSD和(T))使用一个参数的做法是合理的,因为能量差异非常小。同样,将为CCSD(T)优化的外推参数用于其他相关方案(如QCISD(T)或MP2)不会引入显著误差。常用的外推方案通常取其参数来自Neese和Valeev的研究;尽管这些可能不是理想值,但针对S66和MB16-43基准集的评估显示,由此产生的差异很小。最后,建议使用Schwenke型外推,特别是在拟合外推参数时;本研究考虑的其他三种外推方案(指数平方根型、反幂次型和指数型)均可表达为Schwenke型外推。此外,研究表明,对于两个连续基组的外推,由于所有外推方案都简化为Schwenke型外推,不同外推方案在数学上是等价的,给定一个外推参数即可找到另外两个。双杂化DFT泛函也被考虑;已使用的CCSD-T优化参数仍给出合理的外推,并且通常为总能量优化单一参数或为Kohn-Sham和MP2类能量分别优化两个参数的优势很小。BSEF74基准集主要由双原子分子组成。虽然它们覆盖了一系列电子态,但人们不禁要问,更大的参考集是否会给出不同结果,特别是考虑到外推方法之间的差异可能很小。如何扩展参考集的问题很棘手。如果添加S66测试集,则会担心拟合偏向闭壳有机体系。可以添加其他数据集,但总是担心参考集是否平衡——直到达到GMTKN55。问题随后变成计算资源限制:更大的参考集意味着需要计算更多能量,更大的分子意味着更昂贵的计算(并且这些计算并非线性缩放!),并且在某些点aug-cc-pV6Z能量计算将变得过于昂贵。在这种情况下,将被迫依赖使用aug-cc-pV5Z作为外推中的较大基组(直到它也变得太昂贵)或诉诸分辨恒等式(resolution of the identity, RI)等方法使计算更可行(但这引入了全新的复杂性维度——观察到的小差异是来自外推、RI近似还是潜在的两者?)。在所有这些额外成本和努力之后,最可能的情景是本研究的结论不会改变。最终,使用BSEF74参考集可能提供了一个良好的成本与收益平衡。总之,尽管常用外推参数可能不是理想的,但由此产生的误差非常小,意味着没有真正理由重新评估当前基准值。也就是说,向前推进时,建议考虑使用本研究修订的参数。当开始考虑高精度热化学(其中每一点小误差都很重要)时,这变得更加重要。

生物通微信公众号
微信
新浪微博


生物通 版权所有