激发态计算中的溶剂效应

2020-07-27 15:48:59 浏览数 (1)

已经有不少小伙伴在催更,非常感谢大家的支持,也有了更新的动力。关于隐式溶剂模型的介绍,可参见《理论计算中的溶剂效应模型》一文。本文着重介绍在激发态计算中使用隐式溶剂模型时的相关要点,为后面介绍荧光和磷光的计算打个基础。

处理隐式溶剂时,我们使用自洽反应场(self-consistent reaction field, SCRF)模型。在SCRF中,溶剂形成的连续介质会对溶质的电子密度分布产生影响,而溶质的电子分布又会影响多数溶剂模型的参数,因此需要使两者达到自洽。溶质对溶剂会产生两方面的影响:(1)影响电子分布,这是一个快速的响应过程;(2)影响溶剂分子的朝向,这是一个慢响应过程。对于某个电子态的体系,如果溶剂有足够的时间在两方面都形成响应,则称为平衡溶剂(equilibrium solvation,简记为eq);若没有足够的时间响应,只有第(1)部分进行响应,则称为非平衡溶剂(non-equilibrium solvation,简记为neq)。

在激发态的计算中,垂直吸收和发射的过程非常迅速,因此要使用非平衡溶剂,而对几何结构优化过程,有足够的时间进行响应,使用平衡溶剂模型,所以,垂直吸收和发射以及绝热过程的准确计算方式为

隐式溶剂模型中,对响应过程的处理有如下三种方式:

(1)线性响应(linear-response, LR)模型,这是高斯中SCRF使用的默认方法,对TD计算有解析梯度,可用于激发态的几何结构优化。

(2)校正的线性响应(corrected linear-response, cLR)模型,是对LR模型的修正,高斯中使用SCRF(CorrectedLR)关键词实现,可用于结构优化。LR和cLR的反应场都是根据基态的电子密度定义的。

(3)态特定(state-specific, SS)方法,使激发态的电子密度与溶剂的反应场达到自洽,高斯中使用SCRF(ExternalIteration)关键词实现,精度高,但比较耗时,且对TD没有解析梯度,不能用于结构优化。

在高斯中,使用SCRF关键词时,默认使用LR方法,且对激发能的计算默认使用非平衡溶剂(LR, neq),而对激发态的几何结构优化,则使用平衡溶剂(LR, eq),这正好符合激发能和几何结构优化的计算要求,因此,在高斯中大部分计算任务使用默认的SCRF就可以了,无需额外的关键词。若使用SS方法,则默认使用平衡溶剂(SS, eq),因此要实现非平衡溶剂的计算则需要一定地处理。

以下我们尝试重复文献C. Adamo, D. Jacquemin, Chem. Soc. Rev., 2013, 42, 845中图6的结构IV和表1中的第一行数据,以加深对本文内容的理解,其余的数据读者可以自行尝试,其中对发射过程的(SS, neq)计算我们留到下一篇介绍荧光的计算时再做介绍。

1. Gas

1.1. 气相基态的结构优化:

代码语言:javascript复制
%chk=BODIPY.chk
#p opt pbe1pbe/6-31 G* freq

[No Title]

0  1
N  -1.163158   -0.537565    0.211682
C  -1.138190    0.737323    0.194738
C  -2.404453    1.178384    0.166602
C  -3.181710    0.085080    0.180786
C  -2.356907   -0.971556    0.205312
B   0.071316   -1.295105   -0.200072
N   1.171290   -0.504380    0.458008
C   1.143719    0.720606    0.115365
C   0.006866    1.435149    0.123070
C   2.358326   -0.929181    0.315459
C   3.157830    0.047513   -0.138145
C   2.365620    1.120109   -0.269686
F   0.214050   -1.325671   -1.737970
F   0.032138   -2.744157    0.334397
H  -2.741977    2.221496    0.125571
H  -4.278296    0.059146   0.163720
H  -2.601657   -2.039315    0.217958
H  -0.004981    2.519662   -0.066070
H   2.672037   -1.956405    0.547797
H   4.230507   -0.018013   -0.356233
H   2.659595    2.114865   -0.626617
 

1.2. 在第1.1得到的结构上,做TD计算,得到气相的垂直激发能,输入文件如下:

代码语言:javascript复制
%chk=BODIPY.chk
#p pbe1pbe/6-31 G* td geom=allcheck guess=read
 

输出如下:

代码语言:javascript复制
 Excited State  1:      Singlet-A      3.1634 eV 391.93 nm  f=0.4333  <S**2>=0.000
      48 -> 50         0.23678
      49 -> 50         0.66969
      49 <- 50        -0.11067
 This state for optimization and/or second-order correction.
 Total Energy, E(TD-HF/TD-DFT) =  -680.542822576
 Copying the excited state density for this state as the 1-particle RhoCI density.

激发能3.16 eV,与文献完全一致。

2.(LR, neq)

2.1. 以气相结构为初始,在溶液中优化基态结构:

代码语言:javascript复制
%chk=BODIPY.chk
#p opt pbe1pbe/6-31 G* scrf=pcm geom=allcheck guess=read

2.2. 在2.1的结构上做TD计算,得激发能,输入文件:

代码语言:javascript复制
%chk=BODIPY.chk
#p pbe1pbe/6-31 G* scrf=pcm geom=allcheck guess=read td

输出结果为

代码语言:javascript复制
 Excited State  1:      Singlet-A      3.0452 eV 407.15 nm  f=0.5396  <S**2>=0.000
     48 -> 50        -0.17334
     49 -> 50         0.68777
 This state for optimization and/or second-order correction.
 Total Energy, E(TD-HF/TD-DFT) =  -680.560533859
 Copying the excited state density for this state as the 1-particle RhoCI density.

激发能3.05 eV,与气相的差值为−0.11 eV。

3.(LR, eq)

这是不太常见的计算,在2.1的结构上,做td=eqsolv的计算:

代码语言:javascript复制
%chk=BODIPY.chk
#p pbe1pbe/6-31 G* scrf=pcm geom=allcheck guess=read td=eqsolv

输出结果为

代码语言:javascript复制
 Excited State  1:      Singlet-A      2.7720 eV 447.27 nm  f=0.7785  <S**2>=0.000
     49 -> 50         0.70180
 This state for optimization and/or second-order correction.
 Total Energy, E(TD-HF/TD-DFT) =  -680.570572987
 Copying the excited state density for this state as the 1-particle RhoCI density.

激发能2.77 eV,与气相的差值为−0.39 eV。

4. (SS, eq)

在2.1的结构上做TD计算,在SCRF中加入externaliteration选项:

代码语言:javascript复制
%chk=BODIPY.chk
#p pbe1pbe/6-31 G* scrf=(pcm,externaliteration) geom=allcheck guess=read td

使用SS计算时,应当在After PCM corrections, the energy is处读取激发态的能量,并与基态能量相减。本例的结果为

代码语言:javascript复制
After PCM corrections, the energy is  -680.556350794     a.u.

计算得激发能为(−680.556350794 680.672442961) * 27.21 = 3.16 eV,与气相相差0.00 eV。

5. (SS, neq)

该计算需要分两步进行:

5.1. 对2.1所得结构做基态计算,同时将基态溶剂信息写入check文件,使用关键词scrf(noneq=write):

代码语言:javascript复制
%chk=BODIPY.chk
#p pbe1pbe/6-31 G* scrf=(pcm,noneq=write) geom=allcheck guess=read

5.2. 做SS计算,并读入基态溶剂信息,使得关键词scrf(noneq=read):

代码语言:javascript复制
%chk=BODIPY.chk
#p pbe1pbe/6-31 G* scrf=(pcm,externaliteration,noneq=read) geom=allcheck guess=read td

注意,scrf(noneq=write)的写法为G16所支持,若使用G09,需要写成scrf(read),并在输入文件的分子结构后面写上noneq=write,read也是如此。

得到结果:

代码语言:javascript复制
After PCM corrections, the energy is  -680.555726445     a.u.

与基态能量相减,得激发能(−680.555726445 680.672442961) * 27.21 = 3.18 eV,与气相相差0.02 eV。

0 人点赞