激发态计算中的溶剂效应

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

处理隐式溶剂时,我们使用自洽反应场(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。