专栏首页量子化学用Gaussian寻找圆锥交叉点

用Gaussian寻找圆锥交叉点

一、圆锥交叉的含义

圆锥交叉(conical intersection, CI)简言之是指两个态的势能面交叉的地方,此时两个态的能量简并。在圆锥交叉区域,体系可以从激发态以无辐射形式回到基态。关于圆锥交叉的更深入的理论细节可参考Conical intersections: theory, computation and experiment一书。此外,点击文末“阅读原文”可打开一份非常不错的关于光化学计算的讲义。如果打开速度较慢,可在留言区的百度网盘链接获得。

在下图所示的势能面中,反应物经过垂直跃迁成为激发态,并弛豫到激发态势能面的能量极小点。在激发态势能面上还可能跨过能垒(若有的话),到达圆锥交叉点,随后便进入基态势能面,逐步转化为产物。

需要注意的是,S0与T1的势能面交叉称为MECP点更为合适,CI点最初只用于同一自旋多重度两个势能面之间的交叉,现在文献中把圆锥示意图泛泛地标在任两个自旋多重度间的势能面交叉,只表示交叉示意,不代表其是真的是圆锥形状的交叉,读者在读文献时务必注意。

如果仅是搜索S0与T1之间的MECP点,一般基态(U)DFT就能做,有不少程序(easyMECP, sobMECP,及二者的核心MECP程序)支持,计算量不大、操作较为方便。而对于同一自旋多重度的CI点,若用TDDFT,其激发态虽是多行列式的,但是基态仍是单行列式,既然基态、激发态简并,那么从原理上讲基态也应该考虑多个行列式的线性组合。而CASSCF是一种多组态方法,可以同时算多个根,每个根都是多行列式的,因此既可以搜索CI点,也可以搜索MECP点,若只是探究光化学、光物理过程的机理、对精度要求不十分苛刻时,是理想的选择。需要更高精度的话,后续仍需做(X)MS-CASPT2或MRCI。不过这些计算耗时都大,学习门槛也不低,本文暂只涉及用CASSCF寻找圆锥交叉点。

二、计算实例

本文以Exploring Chemistry with Electronic Structure Methods (Third Edition)中的例8.12为例,介绍相关过程的计算。这个例子计算的是苯到盆苯(benzvalene)的异构化过程中的光化学过程,盆苯是基态下能量较高的极小点。所有计算用Gaussian 16 C.01完成。

(1) 用CASSCF优化苯的基态结构。这一步建立在《用Gaussian做CASSCF计算》一文的计算基础之上。在挑选了合适的活性空间后,做结构优化及频率计算:

%oldchk=benzene_cas.chk
%chk=benzene_gs_opt.chk
#p opt freq cas(6,6)/6-31G(d) geom=allcheck guess=read
 

benzene_cas.chk是用DFT或其他方法优化结构,并做CASSCF单点计算生成的chk文件。

(2) 用CASSCF计算第一激发态,获得垂直激发能,同时观察轨道,确保活性空间的正确性。

%oldchk=benzene_gs_opt.chk
%chk=benzene_es.chk
#p cas(6,6,nroot=2)/6-31G(d) geom=allcheck guess=read
 

(3) 用CASSCF优化激发态的结构并做频率计算。

%oldchk=benzene_es.chk
%chk=benzene_es_opt.chk
#p opt freq cas(6,6,nroot=2)/6-31G(d) geom=allcheck guess=read
 

优化得到的激发态与基态相比,C−C键键长从1.396 Å增加到1.434 Å。

(4) 寻找圆锥交叉点。这是一项需要经验和技巧的任务。在exploring3中,使用的方法是利用柔性扫描的结果作为优化交叉点的初始结构。对C1−C2−C3−C6二面角进行扫描:

%oldchk=benzene_es_opt.chk
%chk=benzene_scan.chk
#p opt=modredundant cas(6,6,nroot=2)/6-31G(d) geom=allcheck guess=read nosymm
 
1 2 3 6 S 5 7.5
 

共扫描得6个结构,从输出文件中可以看到每个结构的垂直激发能在逐渐变小。这是符合预期的,因为在圆锥交叉点处基态和激发态的能量简并,垂直激发能应该趋近于0。

扫描得到的第6个结构如下图所示:

再以此结构为初始,进行圆锥交叉点优化:

%oldchk=benzene_scan.chk
%chk=benzene_ci.chk
#p opt=conical casscf(6,6,nroot=2)/6-31g(d) guess=read nosymm
 
Benzene
 
0 1
 C         -0.00080300    1.41440000    0.65507800
 C          1.05274800    0.71961300   -0.07247700
 C          1.17424100   -0.73123000    0.00730800
 C          0.00097200   -1.47535000    0.15754600
 C         -1.17472000   -0.73263700    0.01547600
 C         -1.05143700    0.71680700   -0.07124800
 H         -0.00175900    2.49003400    0.65966300
 H          1.70776200    1.26435800   -0.73164000
 H          2.13962100   -1.19172400   -0.09740800
 H          0.00241700   -2.54058600    0.29322400
 H         -2.14003800   -1.19374700   -0.08641800
 H         -1.70900300    1.26006300   -0.72910300
 
0.500000000.50000000
 

优化圆锥交叉点的关键词为opt=conical,此时CASSCF默认为态平均CASSCF,因此需要在文件末尾加上权重。这一点在Exploring3中未提到,所给的例子中也没有写权重,实际程序运行时会报错。此外,这一步不能直接使用geom=allcheck读取scan步骤中的坐标信息,因为会将扫描的冗余内坐标沿用过来,影响圆锥交叉的优化。所以要重新保存结构,或使用geom=(allcheck,newdefinition)关键词。

优化得到的圆锥交叉点的垂直激发能为0.0012 eV,结构如下图所示:

(5) 观察圆锥交叉点的结构,可以看到,若将C2−C6键拉长,则接近苯的结构,而将C2−C6键缩短,则向盆苯的结构靠近,因此若想得到盆苯的结构,我们尝试手动将键长缩小,如至1.8 Å,以此结构为初始,进行结构优化。

%oldchk=benzene_ci.chk
%chk=benzene_pr.chk
#p opt freq casscf(6,6)/6-31g(d) guess=read
 
Benzene
 
0 1
 C     0.00046100    1.40464600    0.71242300
 C     0.89962637    0.71088527   -0.12565894
 C     1.16331800   -0.73594400    0.00659700
 C     0.00117700   -1.49102900    0.16332400
 C    -1.16241300   -0.73730500    0.01061800
 C    -0.90097037    0.70984973   -0.12239106
 H    -0.00007100    2.47952800    0.74928700
 H     1.53720837    1.27221027   -0.78782094
 H     2.14440600   -1.16412600   -0.07883700
 H     0.00209900   -2.55251600    0.32615000
 H    -2.14327100   -1.16667000   -0.07145300
 H    -1.54156837    1.27047273   -0.78223506
 

不幸的是,优化得到的结构有虚频,是一个过渡态的结构。观察其振动模式:

可知这可能是两个盆苯的异构过程的过渡态。因此,一个比较合适的做法就是借助IRC分析,得到靠近盆苯能量极小点的结构,再做最终的结构优化。可以先做一步粗糙的IRC分析:

%oldchk=benzene_pr.chk
%chk=benzene_irc.chk
#p irc=(rcfc,recorrect=never) casscf(6,6)/6-31g(d) guess=read geom=allcheck

结果如下图所示,IRC曲线很对称,观察两端的结构,确实印证了上述猜想。但是由于IRC所走的步数不够多,还不是很接近盆苯的能量极小点,如果直接用此时两端的结构进行优化,很有可能还是会回到过渡态的结构。

此时,可以重新做一次IRC分析,并只向一个方向走步,增大步数:

%oldchk=benzene_pr.chk
%chk=benzene_irc2.chk
#p irc=(rcfc,forward,maxpoints=50,recorrect=never) casscf(6,6)/6-31g(d) guess=read geom=allcheck
 

得如下结果:

最后一点的分子结构为

最后,以这一步的最终结构为初始结构,优化能量极小点:

%oldchk=benzene_irc2.chk
%chk=benzene_pr2.chk
#p opt freq casscf(6,6)/6-31g(d) guess=read
 
Benzene
 
0 1
 C    -0.30169200    0.96943500    0.72327700
 C     0.79642600    0.70088900   -0.26809600
 C     1.17355300   -0.76059400   -0.20523300
 C     0.04989100   -1.40825900    0.13313200
 C    -1.00131100   -0.33628900    0.28834500
 C    -0.69727600    0.80329900   -0.64511000
 H    -0.48556700    1.70567400    1.48406900
 H     1.49687100    1.48752000   -0.52796400
 H     2.14667400   -1.16482800   -0.42415000
 H    -0.10416600   -2.45960300    0.28151200
 H    -2.01412200   -0.53455000    0.59125500
 H    -1.27067600    1.38386100   -1.34125200
 

最终得到盆苯的结构如下:

在第一段提到的教程中正好有一个描述此过程的势能面,如下图所示:

总结一下便是,基态的苯在光照下得到激发态,在激发态势能面上行走到达能量极小点,并可能进一步跨越能垒(也可能没有能垒)到达圆锥交叉点,顺着圆锥交叉点进入基态势能面上的两个盆苯互变异构的过渡态,再得到盆苯的能量极小结构。同时,在图中也同时显示了在基态的势能面上也存在苯到盆苯过程的过渡态结构。但这个过渡态的能垒很高,不易发生。

本文分享自微信公众号 - 量子化学(quantumchemistry),作者:zhigang

原文出处及转载信息见文内详细说明,如有侵权,请联系 yunjia_community@tencent.com 删除。

原始发表时间:2020-09-02

本文参与腾讯云自媒体分享计划,欢迎正在阅读的你也加入,一起分享。

我来说两句

0 条评论
登录 后参与评论

相关文章

  • 综述:团簇结构全局优化的方法、应用和挑战

    本文用笔记的形式介绍团簇结构全局优化软件ABCluster(点击文末“阅读原文”可进入程序下载页面)的作者Jun Zhang在Int. J. Quantum C...

    用户7592569
  • 在Gaussian16中同时扫描两个反应坐标

    是依次扫描两个坐标,无法做到同时,因此得到的是一张二维势能面,总扫描点数是两个坐标扫描点数的乘积,计算量较大。然而有时候我们只想同时扫描两个反应坐标,即两个坐标...

    用户7592569
  • 用Gaussian做CASSCF计算

    完全活性空间自洽场(complete active space self-consistent field, CASSCF)是一种组态相互作用(configur...

    用户7592569
  • Xiuno如何修改首页和版块列表页默认排序为发帖时间排序

    这里举一个例子: 目前Xiuno首页和版块列表页排序是根据发帖时间+回复时间;按照以下方法修改首页和版块列表页默认排序为发帖时间排序; 找到:/model/th...

    奇梦
  • 通过预测API窃取机器学习模型

    由于机器学习可能涉及到训练数据的隐私敏感信息、机器学习模型的商业价值及其安全中的应用,所以机器学习模型在一定程度上是可以认为是机密的。但是越来越对机器学习服务提...

    FB客服
  • Entity Framework 迁移

    这一篇文章主要讲解EF的迁移,我们前面的文章一直是使用新增数据的方式生成数据库,但是在实际开发过程中,我们会使用代码迁移的方式生成数据库,下面我们来讲解一下代码...

    喵叔
  • 细致入微:Oracle中执行计划在Shared Pool中的存储位置探秘

    这两天我一直在想一个问题,那就是 Oracle 的执行计划到底存储在什么地儿?它会是一种什么样的格式? 这里我试图对这个问题做一点我自己认为的解释,这个解释可能...

    数据和云
  • 利用canvas绘制动态仪表盘

    该library的github地址:http://bernii.github.io/gauge.js/ 创建一个新的BSP application:

    Jerry Wang
  • phpMyadmin 服务简单安全加固

    老七Linux
  • 解决ClashR打开配置空白

    最近发现 ClashR 这款很好用的代理客户端,客户端可以在我的资源站下载:Frytea's Res.

    宋天伦

扫码关注云+社区

领取腾讯云代金券