量子化学研究中切换泛函应当注意的问题
量子化学研究中切换泛函应当注意的问题
Issues that should be paid attention to when switching DFT functionals in quantum chemistry research
文/Sobereva @北京科音
First release: 2018-May-2 Last update: 2024-Sep-24
量子化学研究者经常会用一种做法，即计算一种问题用A泛函，计算另一种问题则切换为B泛函。笔者在线上和线下答疑的时候，都常看到有人不恰当地切换泛函，故此文专门说说这个问题。
一定要明白，如果你要切换泛函，必须必须必须有非常正当的理由。一般有这三种正当理由：
(1) A泛函和B泛函适用的问题不同，因此在计算不同问题时用不同泛函。什么问题适合什么泛函，在此文里有汇总《简谈量子化学计算中DFT泛函的选择》（http://sobereva.com/272）。比如对某个有机大共轭体系，用PBE0优化基态结构，而用CAM-B3LYP计算电子光谱，这就是很合理的。因为PBE0对于优化有机体系基态的结构公认是适合的，在J. Chem. Theory Comput., 12, 459 (2016)的横测中也充分证明了这点。但PBE0对于大共轭体系的电子光谱计算往往不理想（HF成分偏低所致），而CAM-B3LYP则没这个问题，已被大量文献用于大共轭体系电子光谱的计算且得到了不错的结果。但对于基态几何结构优化问题，用CAM-B3LYP则不及PBE0，文献里也很少这么做。再加上CAM-B3LYP是范围分离泛函，耗时也比PBE0高，因此基态不用CAM-B3LYP去优化。
(2) A泛函的耗时明显比B泛函低，而精度也明显比B泛函差。因此用A泛函计算那些对计算级别要求低而且计算量较大的任务，而用B泛函较准确地计算那些对计算级别要求高的问题，即好钢用在刀刃上。最典型的情况就是用便宜的泛函优化，而用明显更好且明显更贵的泛函算单点能。例如计算弱相互作用，我们可以用B3LYP-D3优化结构，而用精度和耗时都显著高一个档次的双杂化泛函PWPB95-D3去算单点能。
(3) 与研究内容、目的有关。比如要对比研究一系列泛函A、B、C、D对于某类分子的极化率、NMR、偶极矩之类问题的计算精度，那就没必要用这四种泛函分别去优化结构和计算性质了，只需要用一个对当前体系合理的泛函去优化即可，之后计算性质的时候再切换成不同泛函去分别计算。
如果没有正当理由，绝对不要切换泛函，否则只会在发文章时给自己找麻烦！内行的审稿人看到你在文章中不同问题用了不同泛函，几乎总会思考一下你这么做的必要性，如果你切换泛函显得没有任何必要，就会产生疑点，被审稿人怀疑并询问这么做的目的。如果你有正当理由，就如实解释即可；但如果你没有正当理由，没法给出令审稿人信服的解释，这就尴尬了。你如果说我就乐意这么干，会把人家惹怒，最终可能被拒；如果你想方设法东拼西凑点理由来应付，可能有的审稿人会被你糊弄过去，也可能有的审稿人依然觉得你的做法不妥但懒得继续追究，而如果碰上特别较真的审稿人，那恐怕最后你只能把大量数据重算来使得所用泛函统一。因此，为了避免潜在的麻烦，如果不是有十分确切的原因、能确保被审稿人问及为什么换泛函时能给出有理有据的回复，就绝对不要切换泛函！之前在一次北京科音的量子化学培训班里，有个学员是做过渡金属配合物的，导师让他用PBE优化结构，然后用M06算能量，他让我提提看法。我问他为什么这两步用不同的泛函，他却说不出来理由。像这种情况，自己都不清楚为什么要换泛函，而且研究经验也不多，就更不能切换泛函了。本来不切换泛函，审稿人万一问起来为什么用那个泛函，可能回复起来都有点吃力；这回又切换泛函，吸引了审稿人的注意力，增加了被提问的几率，到时候还得解释那两个泛函分别用在那两个场合的原因，这无疑是自找麻烦。
哪怕你切换泛函有一定道理，但是理由并没有充分到特别有必要这样做，而且自己对理论方面的自信不足、没能耐回击审稿人的一些错误的看法，那么为了免得审稿人找碴，最好别切换泛函。考虑一种情况，B3LYP做优化和振动分析，M06-2X算单点，假设当前研究的是有机体系。这俩泛函这么搭配，确实有一定理由，因为M06-2X的耗时比起B3LYP要高一些（尤其是在振动分析上），M06-2X对积分格点要求比起B3LYP高很多，优化以及振动分析如果用M06-2X做的话，若不用较高积分格点精度（比如Gaussian里用ultrafine格点），出现虚频、几何优化难收敛的概率还是挺大的，而用较高积分格点精度则会进一步提升M06-2X的耗时。而且B3LYP结合DFT-D3后对大多数体系都能优化得不错，见《谈谈量子化学研究中什么时候用B3LYP泛函优化几何结构是适当的》（http://sobereva.com/557）。因此，B3LYP做opt freq，单点用M06-2X，这么切换泛函确实有实际意义。但如果希望减少被某些审稿人comment的可能性，则建议把泛函统一，即opt freq和单点计算全都用M06-2X。对有机体系，M06-2X做几何优化也是比B3LYP略微好一点点的，特别是存在分子内弱相互作用的时候（当然，对B3LYP可以通过加DFT-D校正弥补这点）以及存在大共轭的时候，何况M06-2X的计算耗时终究还是和B3LYP在一个档次，差距远没有普通泛函和双杂化泛函之间那么悬殊。当然绝不是说笔者不建议用B3LYP优化而M06-2X算单点这种组合，但这么做的时候应当在文章当中交代清楚切换泛函的理由，否则容易令某些缺乏经验的读者和审稿人起疑。至于前面提到的对过渡金属体系，PBE优化，M06算单点，这种做法则是不可取的，因为根本没有切换的意义，M06对这类体系的结果未必比PBE更合理，而且耗时还会提高，因此要用PBE就都用PBE，要用M06就都用M06，真要想结果更好，算单点时候应当用双杂化。
补充：有本文的读者表示，我看那知名的做有机反应计算的xxx大牛，人家就是B3LYP优化，M06L算单点，在高档次期刊发了很多文章。言下之意泛函这么切换才是正确、最佳的。在本文上面一段我刚说完B3LYP优化，单点切换成M06-2X是有理由的，而M06L对有机体系本来就不是什么好的选择，早就众所周知明显不如专为主族体系参数化的同门的M06-2X好，干嘛换成这个？这不是没事找事么？人家名气大，人家的做法就是对的？况且这大牛只不过是应用性计算做得多才出名，不是专门搞理论方法的，存在认识误区、有不恰当的做法是很正常的。而且人家声名在外，审稿人往往还是大牛的熟人，即便审稿人对其切换泛函有看法、觉得并不妥当，只要结果不存在明显、严重不合理性，审稿人看在作者的面子上一般也不会追究。而你和那位大牛不同，你如果盲目效仿这大牛的并不科学的做法，就会有极高的几率被审稿人盯上，还可能不依不饶，到时候你怎么对自己都说不清楚的切换泛函的做法予以合理解释？就算你把大牛这么做的一些文章给审稿人看，审稿人也未必会放过你。要牢记：在不懂基本原理的情况下盲目效仿他人的做法是极度危险的！
对于Gaussian，虽然支持密度拟合方法加速纯泛函的计算，但是由于这方面做得非常差，一般也没人去利用，所以在Gaussian里，纯泛函和杂化泛函的耗时基本没有区别。但是，在ORCA、Turbomole等程序里，由于DFT的密度拟合做得很出色，纯泛函比杂化泛函在耗时上对大体系能低一个数量级以上（见《大体系弱相互作用计算的解决之道》http://sobereva.com/214），因此如果你是这些程序的用户，推荐用纯泛函来充分节约高耗时步骤的计算时间。“xx 步骤使用纯泛函是为了节约时间”也是一个切换泛函的非常正当的理由，足以说服懂行的审稿人。对于ORCA用户，比如使用BLYP-D3(BJ)或B97-3c做优化和振动分析，而用ORCA支持的很好的wB97M-V算单点，这是完全可以接受的而且很恰当的做法。纯泛函有一些显著的软肋，比如计算(超)极化率会明显高估、计算激发能（特别是电荷转移态）会明显低估等等，算这些任务可千万别为了图便宜用纯泛函。
有些人，算不同问题的时候切换泛函不是出于计算量的考虑，也不是出于计算精度的考虑，而是为了凑实验数据，觉得计算结果与实验越接近就越好、越容易发表。其实这种做法是不好的，毕竟理论计算并不是实验的陪衬，没必要刻意去一味地迎合实验。盲目去让结果贴合实验，那做理论计算还有什么意义？每种泛函对于某类体系的某种问题的计算终究是有误差的，只要选用的泛函在原理上合适、计算过程各方面要素都已经考虑周全，存在一定误差应当坦然接受，只要文章中的理论计算能把要考察的问题分析、讨论清楚就行了。众所周知，泛函的HF成分越高，算出来的激发能就越高，于是有很多人利用这个趋势，偏要挑一个激发能和实验值最接近的泛函来算激发能，然后在文章里得意地说，本文计算结果和实验相符非常好。这样的做法其实真没什么意思，而且如果你为了达到这个目的，用了一个对当前问题很少使用的泛函，那么你凑数据的意图很可能会被有经验的审稿人识破。比如算个有机体系的UV-Vis光谱，之前优化的时候你用的B3LYP（HF成份20%），但可能你发现B3LYP略微高估了激发能，于是尝试使用HF成份更低（10%）的TPSSh泛函来算，并发现其结果误差相对于实验来说比B3LYP更小一点点，于是就在文章中用TPSSh算光谱了。这么做就给人感觉在动机上非常可疑，毕竟这个泛函极少用于有机体系的激发能计算，若是审稿人追究用这个泛函的原因，那么你很难给出正当理由，如果坦白交代就是单纯为了凑实验数据，这必然不会给审稿人好印象。因此，对当前这个例子，索性优化和激发能计算都用B3LYP，这就不会令审稿人觉得有什么疑点，就算误差大一点也比起惹质疑好。还有一种挺常见的情况是研究者比较了几种泛函，发现计算吸收光谱时候A泛函最接近实验，但计算发射光谱的时候B泛函最接近实验，如果因此就在文章中算吸收用A，算发射用B，那也是非常不合适的，这会让审稿人觉得动机不纯，毕竟不管是吸收还是发射都是激发能计算问题，用俩泛函是什么鬼？这种情况，就应该用一个对当前类型体系用得比较多的，而且实际发现算吸收和发射的误差都在可容忍范围内的泛函，别凑实验凑得太露骨。
前些天有人在思想家公社QQ群里问“请教用m062x优化，b3lyp-d3做单点，这个配置合理吗，能量的方法比优化的小了，主要觉得m062x算能量有点低估”（我不知道他算的什么体系，姑且视为是有机体系吧）。像这种做法，显然非常不合理。有经验的人都知道B3LYP-D3本身就比M06-2X便宜，通常精度还更低，这么做明显违背一般逻辑。若碰到懂行的审稿人，要么会说你缺乏常识而批你一通，要么会怀疑你在刻意凑数据而要求你给个说法，说不定还会觉得你在计算时由于忽略了什么重要因素，导致本应该精度较好的方法M06-2X算出来的结果差，逼得你把M06-2X的数据隐藏起来，于是要求做更多解释...假设对于当前问题，由于M06-2X自身存在的不足，可能算能量时表现得就是比B3LYP-D3要差，就是应当用B3LYP-D3才最合适，那你何不在优化时候也用B3LYP-D3，使得整个计算过程用的泛函统一？这样从逻辑角度上显得更自然，也免得因切换泛函被别人质疑。
关于加DFT-D3的问题，在《谈谈“计算时是否需要加DFT-D3色散校正？”》（http://sobereva.com/413）中有专门讨论，本文不再详述。诸如优化用M06-2X而单点用M06-2X-D3这种做法是纯属自找麻烦的，虽然通常来说对M06-2X加不加D3影响很小，但你似乎都知道加D3会带来好处于是在算单点时候加了，有什么正当理由在优化时候不加D3？何况是完全免费的校正，目前主流程序又不是不支持优化时候带D3。而诸如优化弱相互作用复合物用B3LYP，算相互作用能的时候才改用B3LYP-D3(BJ)，那文章被秒拒一点都不冤枉。众所周知B3LYP完全不能描述色散作用，而弱相互作用问题又和色散作用密切相关，B3LYP优化这种体系往往连结构的定性正确都没法保证。
如果你用的后校正方法会明显增加计算量，或者用于优化过程有显著障碍，而且在优化时候其必要性又没那么高，那么只在计算能量的时候加校正才是恰当的。比如解决BSSE问题的Counterpoise (CP)校正，用了它会导致失去解析梯度而令优化任务计算耗时极大地提升，而且从实际上看只要基组不是很小且理论方法合理，不用CP也不至于令优化结果误差明显，因此CP只适合在最后算能量的阶段才考虑。而至于Grimme提出的gCP方法解决BSSE问题，它类似于DFT-D3一样是免费的后校正，因此在用的时候也应该像DFT-D3一样，优化和单点要加就都加，要不加就都不加。
众所周知，几何优化、振动分析、产生IRC这三个任务，用的计算级别必须严格相同，显然更是绝对不能切换泛函，不管你有什么理由都不行。比如在ORCA 4.0里用M06-2X优化了分子，但是在ORCA 4.0里面M06-2X没有解析Hessian而只能用数值Hessian，因此算频率巨慢，故想改用B3LYP做振动分析，这显然是绝对不行的。
凡是需要横向对比的计算，在计算时也不要切换泛函。不同泛函给出的绝对能量没有丝毫可比性，没法互相求差，比如过渡态用M06-2X，极小点用B3LYP，它们的能量求差算出来的势垒毫无意义。如果比较的问题是基于相对能量算的，比如算一批体系的反应能垒来比较谁的反应活性强，有的体系用wB97XD，有的体系用M06-2X，这也是不行的，因为不同泛函都存在着不同程度的系统性误差，只有相同泛函的计算结果之间对比，从原理上讲系统误差抵消得才较好，从而能反映出研究的不同体系自身特点的不同。但为了遵从这个原则，有时候可能会遇到一些麻烦的情况。比如研究一个含过渡金属的复合反应，有的阶段的基元反应只涉及主族原子，而有的阶段则涉及过渡金属。众所周知，M06-2X适合主族问题，而反应直接牵扯到过渡金属时，相对于M06-2X而言，其同门的M06则更为适宜。那么计算这种体系的整个反应过程，能不能有的用M06-2X，有的用M06？这样做是不合适的，给审稿人解释起来也困难。像这种情况，索性用一个倾向性没那么强的泛函，比如后来提出的明尼苏达系列泛函MN15，它号称对主族和过渡金属考虑得比较均衡；或者用传统的且普适性较强的B3LYP、PBE0之类。用它们来算，有可能整个反应中涉及主族的基元反应精度不如用M06-2X，涉及过渡金属的基元反应精度不如用M06，这没关系，最后都用精度高普适性也好的双杂化算个单点就完了（ORCA里双杂化计算由于能利用密度拟合技术，并不很昂贵）。
再看一个情况，优化一个体系基态用B3LYP，而做电子激发计算发现它的某个激发态是电荷转移激发态，由于B3LYP算电荷转移激发态不好（而且有文章证明其优化电荷转移激发态结构也同样不好），于是在优化这个激发态的时候改用了长程校正泛函wB97XD。这种做法正当不？毕竟这是同一个体系，又都是优化任务，却用不同泛函来做，是不是这样不妥？实际上这样做是没有明显问题的，因为激发态和基态的研究首先就不属于同一类问题，而且你切换泛函的理由也比较明确。而如果你就是希望统一，完全避免潜在的麻烦，那你基态和激发态都用wB97XD也完全可以，本身wB97XD算基态（至少对于有机体系来说）也是很适用的。注意如果牵扯到基态和激发态电子结构特征直接比较，特别是要求差的时候，那么切换泛函是不行的。比如，要通过Multiwfn绘制垂直跃迁过程的密度差图，用wB97XD算激发态密度，而用B3LYP算基态密度，这样的结果没有任何意义，因为不同泛函下电子结构还是存在不小差异的，这么求的密度差就没法如实体现出电子跃迁本身特征了，而会引入显著的人为因素，甚至使得密度差图是完全误导性的。
上面虽然说的都是切换泛函，但上述讨论也可以推广到切换计算级别上，即DFT与波函数方法间的切换、基组的切换。比如优化时候用6-31G**，计算单点时候用def2-SVP，这就不合适了，本来都是同一档次的基组，切换干啥？意图何在？再比如，某有机体系，优化时候用M06-2X，算单点时候用MP2，这几乎是作死行为。在有机体系能量计算上，MP2的精度明显不及M06-2X，耗时还比M06-2X高很多，这和“算单点用的精度应该至少不低于几何优化”这个量化计算常识完全相违背！如果你在优化和单点计算时都用MP2，可能有的内行审稿人只是在心里吐槽一下：都什么年代了还用MP2，这作者有没有看过近10年的研究文章啊。但未必会写comment；如果优化和单点都用M06-2X，有经验的审稿人可能心想：其实单点用精度更高的双杂化更理想，也未必会写comment；而如果你优化用M06-2X算单点却用MP2，这位审稿人看到作者连基本常识都没有，绝对忍不住要在审稿意见里写comment批一顿了。
以上讨论还可以进一步推广到切换溶剂模型、DFT积分格点等方面，同样也是如果没有特别原因，建议在整个计算中保持这些设定的统一。对于溶剂环境中的问题的研究，我们多数情况用隐式溶剂模型，而隐式溶剂模型一般也就增加耗时的几分之一（用线性响应溶剂模型时），即曰只要你在真空中算得动，你在溶剂模型下也肯定算得动。虽然大多数情况，溶剂效应对分子结构的影响是很小的，在气相下优化就够了，但是毕竟优化时考虑溶剂模型对耗时增加才那么丁点，何不在优化时也用溶剂模型，令整个计算保持统一，省得被审稿人找碴？而且，有些体系（特别是局部带显著电荷的体系），溶剂效应对结构影响非常明显，甚至是定性的，你始终都带溶剂模型，也消除了因为某些情况没考虑溶剂效应导致结果不合理的风险。也有另一种情况，比如我要研究某个分子在8种溶剂下的垂直电离能，我根据经验可以判定溶剂效应对这个分子结构的影响是可以忽略不计的，于是我只在气相下优化这个分子，省得在各个溶剂下都优化一遍了，大大节约了计算量。这样做是完全可以的，理由是非常正当的。
显然，本文不可能把各种各样理论研究过程中切换泛函的情况全都包罗万象讨论一遍，对其它的情况大家应当根据自己积累的理论知识判断切换泛函的合理性、必要性，如果实在难以判断，应求助真正懂行的人。
如何向审稿人证明自己用的泛函（或其它类型的理论方法）是合适的？如何选择合适的泛函？这个问题和本文有一定联系，在下面笔者进行了归纳（前述的《简谈量子化学计算中DFT泛函的选择》里笔者对泛函的选用推荐也是根据以下方式来的）。
(1)看专门的对不同泛函计算某种问题的精度的对比测试（benchmark）文章。笔者还专门有相关博文：《谈谈量子化学研究中什么样的benchmark才有意义》（http://sobereva.com/554）。这类benchmark文章在主流的理论/计算化学刊物，如JCTC、JCC、PCCP、JCP等期刊上可谓多如牛毛，哪怕是那些很冷门的问题（比如转动常数、原子核处的电场梯度）都有人做过横测。从中你会知道对这类问题哪种泛函表现比较突出从而值得使用、哪些泛函精度烂或者用起来很危险。如果审稿人怀疑你用当前泛函的合理性，直接把相关的横测文章塞给他即可。但要注意，哪怕是类似体系的类似问题，由于不同测试文章选取的测试集不同，往往出现两篇文章结论大相径庭的情况，这就需要你多看测试文章，看一篇往往不够，起码得看至少三篇，从而逐渐形成自己的认识，不要盲目被某一篇测试文章的结论牵着走。
(2)看大部分类似体系和问题的计算文章都用的什么泛函。如果某个泛函在某类体系的某种问题上的计算被广为使用，成为了研究者们优先使用甚至默认使用的泛函，而且结果也确实不错，那么这个泛函大可以放心去使用。虽然未必这个泛函对这类情况真的是最最理想的泛函，但至少一般不会出现定性错误的情况，比起用其它乱七八糟的泛函放心得多。如果一个泛函对这种情况的计算表现得经常很差，它也注定不会流行起来。这些相关计算的文章看多了，当审稿人问你为什么用当前的泛函算当前这类问题，塞给审稿人一大把较高水准的用这个泛函做的相关研究文章往往就够了。
(3)进行实际对比测试。比如当前研究了一大批体系，如果其中某些（或者之前考察过的类似的体系）有实验值（或者自己用高精度级别计算了其可靠的理论值），你同时用了好几种泛函做了计算，发现其中某个泛函精度高且稳健性高，那么你用这个泛函去计算类似的体系一般也是靠谱的，也不怕被审稿人怀疑，因为论据很充分。
(4)根据研究经验和理论分析。如果上述三类论据你都找不到，那就只能告诉审稿人，根据我的大量研究经验，用xxx泛函是合理的（最好能给出已发表的文章）。或者，如果你有一定理论水平，你可以表示根据xxx泛函的特点（比如HF成份、拟合参数用的数据集的特点、对某些问题的计算表现出来的普遍特点等），此泛函用于此类计算在原理上是恰当的。但这种说辞主观性太强，论据不很充分，审稿人未必会信。最终能否说服审稿人，这在很大程度上看你的语言表达能力和理论功底了。
强烈不建议大家效仿极个别文章，或者道听途说，用一些冷门、不知名，而且缺乏充分论据以证明其适用性的泛函。比如前一阵，群里有个人问BVP86和BP86有什么区别，为什么用BP86而不用BVP86。这个问题很莫名奇妙，BVP86又没什么特别的好处，又非常冷门，干嘛想要用它？用google学术搜索去搜BVP86，结果才不到300条，而搜索BP86，结果多达23000多条。虽然流行程度和泛函的好坏没有必然关系，但你去用一个审稿人可能闻所未闻的极小众泛函，负责一些的审稿人不可能不让你解释原因，而你又无从解释，也没法利用已发表的大量文章从侧面体现合理性，显然是在给自己挖坑，被要求重算或者最终被拒也是自找的。
有些审稿人缺乏基本计算常识，对于合理的泛函切换瞎吐槽。比如我在网上答疑时看到有人问“老师请教您审稿人提出这样的要求，怎么回复？：Regarding the computational method used, the B3LYP/6-31G(d,p) method was used for geometry optimization and M062X/def2-TZVP was used for energy calculation. Using different functionals for geometry optimization and single-point energy calculation may lead to inconsistent energies and errors in the relative energy calculations. The authors should check the consistency of their results at the same level of DFT methods, e.g., B3LYP/6-31G(d,p) optimization and B3LYP/def2-TZVP single-point calculation.””
像这种情况，在回复审稿人意见时就按照本文讨论，明确说清楚为什么切换了泛函，当前显然是为了得到明显更好的结果而在单点用了更昂贵也更好的泛函（而这个泛函用于当前体系的优化则没有明显优势所以优化时候不用）。建议同时引用一些文献切换泛函的例子作为这种行为的支撑，比如可以引用《理论设计新颖的基于18碳环构成的双马达超分子体系》（http://sobereva.com/684）里介绍的笔者的研究文章，其中用wB97XD做几何优化，而用算弱相互作用精度更好的wB97X-V算单点来得到弱相互作用能。还可以再引用《理论设计由18碳环与锂原子构成的电场可控的光学开关》（http://sobereva.com/630）介绍的笔者的文章，其中用wB97XD做了几何优化和过渡态搜索，而用了明显更精确的DLPNO-CCSD(T)算单点得到锂离子迁移的势垒和能量变化。
PS：本文涉及的话题较杂，一些内容和讨论与以下文章有关，有兴趣可以看看
谈谈隐式溶剂模型下溶解自由能和体系自由能的计算
http://sobereva.com/327
浅谈为什么优化和振动分析不需要用大基组
http://sobereva.com/387
各种后HF方法精度简单横测
http://sobereva.com/378
大体系弱相互作用计算的解决之道
http://sobereva.com/214
乱谈DFT-D
http://sobereva.com/83
DFT-D色散校正的使用
http://sobereva.com/210
谈谈BSSE校正与Gaussian对它的处理
http://sobereva.com/46
乱谈激发态的计算方法
http://sobereva.com/265
不同DFT泛函的HF成份一览
http://sobereva.com/282