简谈Gaussian里找过渡态的关键词opt=TS和QST2、3

简谈Gaussian里找过渡态的关键词opt=TS和QST2、3

A brief introduction to the keywords opt=TS and QST2/3 for finding transition states in Gaussian

文/Sobereva@北京科音  2019-Feb-3


Gaussian里找过渡态的关键词有三:opt=TS、opt=QST2、opt=QST3。经常在网上看到有很多初学者居然把QST2当做找过渡态默认用的关键词,明显是被一些坑爹言教程毒害的后果。这里就抽一点时间非常简单地说说常用的opt=TS关键词,以及为什么我很不建议在一般情况下用QST2、QST3。本文并不举实际计算例子,我认为学习找过渡态必须通过大量例子才能培养出感觉,在北京科音初级量子化学培训班(http://www.keinsci.com/workshop/KEQC_content.html)我会给出涵盖各种类型反应的过渡态搜索例子,在相当于中级班的北京科音中级量子化学培训班(http://www.keinsci.com/workshop/KBQC_content.html)里我会讲得明显更深入、给出更复杂更高难的例子,并介绍诸多高阶的经验技巧。推荐参加这两个培训。关于过渡态搜索的具体算法,本文也不做介绍,感兴趣者参看各种方法的汇总《过渡态、反应路径的计算方法及相关问题》(http://sobereva.com/44)。


1 关于opt=TS

一般情况下我都建议用opt=TS找过渡态,其算法和优化极小点结构非常类似,只不过它要求优化的是势能面的一阶鞍点。具体来说一般是用opt(TS,calcfc,noeigen)关键词,calcfc代表精确计算初始结构的Hessian矩阵,noeigen代表优化过程中不对每一步做Hessian矩阵本征值数目的检测(这个要求太严而且没太大必要)。用opt=TS只需要提供一个初猜的过渡态结构就够了,很方便。只要初猜结构和实际过渡态比较接近,而且计算级别和模型合理,一般都能收敛到实际过渡态。一般建议再同时写上freq关键词使得优化结束后自动做个振动分析。

怎么验证过渡态找没找对?如果优化出来的过渡态结构在期望的反应物和产物之间,同时有且只有一个虚频,一般能有6成把握认为过渡态找对了。如果通过gview查看,发现仅有的那个虚频的振动方向对应反应坐标方向,那8成可以认为过渡态找对了;如果走IRC还能连通自己期望的反应物和产物结构,那么说明100%找对了。走IRC的相关信息看《在Gaussian中计算IRC的方法和常见问题》(http://sobereva.com/400)。

opt=TS常见问题如下:
·如果优化出现震荡、难收敛,但是用gview观看结构变化感觉收敛趋势对,可确认优化过程是在往过渡态模样的结构走,最后一帧结构有过渡态的大致模样,按照《量子化学计算中帮助几何优化收敛的常用方法》(http://sobereva.com/164)适当处理即可。
·如果是找出来的过渡态不仅存在一个对应反应坐标方向的虚频,还存在额外的虚频,按照《Gaussian中几何优化收敛后Freq时出现NO或虚频的原因和解决方法》(http://sobereva.com/278)消除多余虚频即可。
·如果优化不收敛或者报错,且最后一帧结构完全不像实际过渡态,或者虽然任务正常结束了但是根本都没有虚频,或者虽然只有一个虚频但是并不是对应当前反应方向的,那就是说明自己提供的初猜结构不合理(也有可能是计算级别不合理,导致当前级别势能面上根本没有对应实际过渡态的一阶鞍点),那就再反复尝试其它自己觉得靠谱的初猜结构,或尝试其它自己觉得合理的计算级别。

找过渡态和优化极小点一样,不需要用到太昂贵的计算级别,见《浅谈为什么优化和振动分析不需要用大基组》(http://sobereva.com/387)。找过渡态和优化极小点用的计算级别应当一致,因为并没有用不同级别的合理理由。对于溶剂下的反应,特别是极性溶剂下,优化过渡态的时候最好就带着隐式溶剂模型,因为有的时候可能有定性影响。对于某些反应,溶剂可能直接影响反应机理,此时用隐式溶剂模型时还得同时用显式溶剂模型。

哪怕是有丰富经验的研究者,也没有那个敢打包票说看见一个反应就能保证一次就找到过渡态。找过渡态的难度远大于优化极小点,因为对初猜结构要求高。显然初学者不应当因为对于一个反应尝试了几个过渡态初猜都失败就简简单单放弃,或者投奔QST2。如果体系比较大,每次试初猜都比较耗时,可以先用更便宜但对于当前体系也能给出定性合理结果的计算级别,等收敛到了期望的过渡态,再把这个结构当做过渡态初猜,结合用原本打算用的更好的级别去进一步准确优化过渡态。其实对于大部分有机反应(也有很多例外),用PM6、PM7这样的极便宜的半经验方法就可以给出定性合理的过渡态结构,由于很便宜,因此建议同时加上calcall增加收敛成功的几率。另外也可以考虑尝试用现在很火的Grimme的GFN-xTB,普适性比半经验方法更好而且速度飞快,见《将Gaussian与Grimme的xtb程序联用搜索过渡态、产生IRC、做振动分析》(http://sobereva.com/421)。

有时候不好确定过渡态初猜结构的时候可以借助柔性扫描,参看此文柔性扫描部分的讲解和例子:《详谈使用Gaussian做势能面扫描》(http://sobereva.com/474)。


2 关于QST2、3

网上有一堆乌七八糟的Gaussian学习资料,还有官方的exploring那本书的第3版,给初学净拿QST2举例子,对opt=TS居然反倒不怎么提,导致初学者以为Gaussian里找过渡态本来就应该用QST2,有经验的人都知道这真是赤果果地在坑初学者。相对于最常用的opt=TS,用QST2有很大弊端:
(1)使用QST2必须同时提供要研究的基元过程的反应物和产物的结构(不需要很精确也不需要非得事先优化,只要模样差不多对就可以)。这比起opt=TS只需要提供一个过渡态初猜结构显然在准备工作上更麻烦。而且有些情况你根本都不好估计其中一方是什么结构。比如双分子反应(或单分子裂解反应),你往往根本就不知道反应复合物(或产物复合物)长什么样子,得找出过渡态之后走完IRC才知道,如果你试图去猜的话,这比起你用opt=TS且直接给个合理的过渡态初猜的难度通常更大。
(2)用户必须确保提供的反应物和产物结构间原子顺序完全对应,否则任务一定失败,而保证这点往往很麻烦。虽然gview里提供了connection工具可以辅助你去把两帧结构里的原子顺序进行自动匹配,但是也往往匹配不成功,或者虽然表面上告诉你匹配成功了但其实还是有原子根本没有正确匹配,你还得手动根据两帧的结构图和原子序号通过反复交换原子使顺序一致,相当费事!
(3)QST2本身不是不需要过渡态初猜结构,只不过QST2计算过程中过渡态初猜结构是程序根据你提供的反应物和产物坐标根据插值产生出来的(优化过程中QST2算法还会根据反应物和产物结构对优化方向进行引导)。这样简单粗暴、机械化产生的过渡态结构,怎么可能比起有基本化学直觉的人摆出来的过渡态初猜结构更合理?QST2用多了你就会发现,QST2自动产生的初猜往往很差劲,虽然有时候很烂的初猜结构最终还是收敛到了实际过渡态,但是达到收敛所需的步数明显高于你用opt=TS并直接提供一个基本靠谱的过渡态初猜,因此浪费了很多计算时间。还有的时候,哪怕你提供的反应物和产物结构基本和实际相符,QST2也会给出一个一看就离谱的初猜结构,导致根本不可能收敛到实际过渡态去(即便QST2还会利用反应物和产物结构来引导优化也于事无补)。
(4)Gaussian的QST2代码本身也往往有问题,有时候明明输入文件完全正确,但是算起来却会莫名其妙报错。

老有一些初学者净以为QST2是万灵药,好像给了反应物和产物结构就一定能找出对应的过渡态,明显这是极度荒唐的,QST2可能哪有这么聪明?甚至于,实验上发现反应物A最终转化到产物B,而且一看就不是一个基元过程,居然有的初学者也拿A和B的结构直接输到Gaussian里跑QST2。这么做,整个复合反应中哪个过渡态的结构都优化不出来!PS:目前倒是也有一些自动找过渡态的算法,比如GSM、GRRM等等,思想、流程都远比QST2要复杂得多,可参考《寻找势能面上全部极小点、过渡态的算法GRRM简介》(http://sobereva.com/319)、WIREs Comput. Mol. Sci., e1354 (2017)、J. Chem. Theory Comput., 13, 5780 (2017)。

还有初学者用QST2的理由是“opt=TS算不出来”。我特别讨厌“算不出来”这种含糊不清的描述,因为等于什么都没说。要知道,只要恰当地使用opt=TS,什么反应的过渡态都能算出来!原理上,opt=TS找不出来的过渡态,用QST2也绝对不可能找到。如果你已经基于基本合理的过渡态初猜结构、在合理的设定下正确地做opt=TS任务,最终也没找出来过渡态,那么指望用QST2就能找出来你想要的过渡态,概率等于连续三次单抽都出UR。

QST2也不是一无是处,在恰当的时候使用也能帮助你更方便地找过渡态,但前提必须是你已经能够熟练使用opt=TS找过渡态,而且已经充分了解QST2的思想和优缺点。可以用QST2的情况有下:
(1)虽然你的头脑里知道要找的过渡态大概长什么样,自己在gview里慢慢摆也能摆出来其合适的初猜,但是涉及的几何变量太多,摆起来太花时间或者确实就是不好操作,与此同时反应物和产物的大致结构比较好确定、也容易确保原子顺序对应,那么可以先尝试用QST2,没准直接能优化到过渡态,从而节约了你反复尝试摆合理过渡态初猜的时间。
(2)对于一个基元反应,限于你的化学直觉,实在不知道过渡态应该长什么样,你可以用QST2碰运气,但是绝对不要抱太大期望。

值得一提的是,QST2计算的时候不需要像opt=TS那样需要由用户提供初始的Hessian,即哪怕calcfc都可以不写,因此对于Gaussian里一些没二阶解析导数的方法如CCSD也可以用(当然了,在CCSD级别下找过渡态一般是愚蠢的做法,贼贵而且没必要),也可以避免一些后HF方法在计算解析Hessian时的高昂代价。不过,QST2关键词里也写上calcfc(或者更贵也更好的recalc=x、calcall)有助于增加收敛到过渡态的几率。
PS:其实opt=TS也不是必须得明确提供初始的Hessian,比如用opt(TS,modRedundant,noeigen)或者opt(TS,gediis,noeigen)的话也可以不用calcfc,但当然收敛成功的几率要明显低于opt(TS,calcfc,noeigen)。

至于QST3,不仅像QST2一样需要提供反应物和产物结构,还得再提供过渡态初猜结构,麻烦得要死,而收敛到过渡态的成功几率并不比直接用opt=TS强多少,基本没有实际价值,因此无需考虑。

opt=TS、QST2、QST3找出来的过渡态结构在原理上是没有差别的。因此,不管用的哪种方法去找,只要最终的结构通过前面提到的方式验证可以确定找对了,那就绝对没有必要再用另外的方法再去计算一次。