前言
在上一期(「龙讯module小课堂」“光”怪陆离:PWmat计算光学性质(一))中,我们介绍了在DFT的基础上利用RPA求解介电函数的方法。虽然DFT不能真正的描述体系的对外场的动态响应,但是RPA的结果依然能在一定条件下满足需求。这是因为DFT计算低估的fundamentalgap有时恰好可以在数值上对应体系的opticalgap。作为代价,RPA可能会丢失一些介电函数的细节,例如峰值的相对大小。图1中展示了实验(红色点)和RPA计算(蓝色线)所得的吸收谱(介电函数的虚部),即使通过平移使计算与实验的吸收边对齐,RPA的结果也不能给出正确的第一吸收峰的相对大小,这就是“丢失细节”。
图1RPA计算的Si的吸收谱与实验结果的对比1由此引出了一个问题:是否存在一种不那么昂贵却能正确描述动态响应的方法?本文将从响应函数开始,介绍线性响应(linearresponse,LR)理论;随后介绍能求解和时间相关的电荷密度的方法——含时密度泛函理论(time-dependentdensityfunctionaltheory,TDDFT),以及在线性响应理论下的LR-TDDFT,并给出该方法的适用范围;最后介绍PWmat中基于LR-TDDFT计算光吸收的module。
介电函数与响应函数
响应函数χ是用来描述物体对外部激励δA的响应δB的函数
它是物体本身的属性,需要外部激励才能显现,但是它的存在并不依赖于外部激励。大家在中学时期学过欧姆定律,知道将纯电阻接在导通的电路后,电阻中会有电流通过。此处电阻两端的电压就是激励,电阻中的电流就是响应。中学老师最爱考的问题:是不是电压越大电阻越大?答案是否定的,电阻大小是纯电阻本身的属性,只是不放进电路的时候,我们不知道电阻有多大。按照定义,我们可以将电阻的倒数视为响应函数
常见的响应函数还有两种:
1.磁化率:磁矩对外加磁场的响应
2.电极化率:电荷密度对外加电势的响应
细心的读者可能已经在上一期发现
因此,介电函数的虚部和电极化率的虚部在形式上完全一致。在得到虚部后,就可以通过Kramers-Kronig(KK)变换得到介电函数的实部。
于是,我们可以将求解介电函数的问题转化为求解响应函数(电极化率)的问题。外加势可以分解为不含时的部分和含时的部分
而响应函数和电荷密度的变化存在多阶的展开
其中二阶及以后的项为非线性项,它们只会在扰动较强时出现。实际上,正是得益于激光这种高强度的扰动,人们才能观察到非线性光学的效应。当外场扰动较小时,可以近似认为电荷密度对外场的响应是线性的,即只考虑一阶项。真正的响应过程涉及在t’时刻施加在r’处的外场vi引起的t(tt’)时刻位于r处的电荷密度n的变化。于是在线性响应下可以有
Kohn-Sham(KS)方程可以视作有效势vKS下的单电子方程,其中的vKS是一个静态的势场,但是我们可以使用一种简单粗暴的方法在它和外场之间建立联系
其中
它的形式是已知的,可以直接用KS方程的计算结果表示
另一方面
于是
其中vH和vxc本来就是电荷密度n的泛函
全部经过Fourier变换到频域空间后
这样我们就在已知χKS的情况下,利用fxc这个与交换关联泛函有关的核心,构造了χ的自洽关系,可以由此迭代求解得到χ。这个自洽关系可以由KS本征态和占据数展开为一个矩阵,该矩阵的对角化问题又被称为Casida’s方程组,它的解正是激发态的能量以及χ的展开系数。
看起来万事俱备,只欠东风,但是这个“东风”来得并不容易,我们仍然需要一种理论来描述电荷密度和交换关联泛函随时间的变化——TDDFT。
Runge-Gross定理与TDDFT
上一部分已经包含了TDDFT的核心思想,但是电荷密度随时间变化的形式仍不明朗。DFT中有Hohenberg-Kohn定理能保证基态的势函数的泛函与基态的电荷密度存在严格的对应关系。问题在于,这种对应关系能否推广到含时的势函数中?
Runge,Erich和Gross在其年的工作(Physicalreviewletters,52,12)中严格证明了:将相互作用电子系统置于单电子含时有效外场veff(r,t)中,假定veff可以对t做Taylor展开,且已知体系的基态。那么,不能同时存在与体系的电荷密度n(r,t)对应的两个不同的veff(r,t)和veff’(r,t)(除非两者相差一个t的函数)。这就是Runge-Gross定理,它在TDDFT中的地位正如Hohenberg-Kohn定理在DFT中的地位。该理论说明,在含时体系中,我们仍然只需要通过计算体系的电荷密度来得到体系的一切信息!有了这一层保障,我们就可以大胆的把基态的KS方程向含时的情况拓展。
如果说DFT的核心是交换关联泛函,那么TDDFT的核心就是交换关联泛函和fxc核心。2获得准确的交换关联势已属不易,以此构造fxc核心的误差无法估计。为了使问题简化,人们再次引入了绝热近似(区别于Born-Oppenheimer近似)。在绝热近似下,体系的交换关联势对电荷密度变化的反应都是瞬时的,且不具备对上个时刻的记忆
在这个近似下,所有与时间有关的变化都是相对于初始时刻的变化,构成Casida’s方程仅仅需要初始时刻的静态交换关联势
fxc核心也与频率无关
而电荷密度的一阶(线性)含时修正则可以表示为
相比于DFT中的电荷密度和哈密顿量的自洽迭代,TDDFT中看重的则是电荷密度和响应函数的自洽迭代,在求得某个时刻的电荷密度后,自然可以构造该时刻的哈密顿量。至此,有的读者可能会问:这不就是耍了个滑头,把DFT的结果加工了一下然后编出个激发态吗?
TDDFT确实需要使用DFT计算得到的本征态和占据数,在绝热近似下,可以直接通过基态的交换关联能构造Casida’s方程组。但是,与早期的在DFT框架下实施的激发态方法3-5不同,LR-TDDFT直接求解了激发谱。在有限的系统中,能级是离散的
这就是电荷密度响应函数的Lehmann表示,区别于由KS本征值构建的χKS,此时响应函数的极点对应的能量Ej-E0是有着严格物理定义的激发态能量。有趣的是,根据Hohenberg-Kohn定理,这种线性响应的密度响应函数,也属于基态的电荷密度的泛函。在固体中,激发谱变成准连续的,响应函数中的求和需要相应的变成积分,或者是选定k-mesh加权再求和。此外,固体中还会存在更复杂的集体激发(元激发),例如等离激元(plasmon)和激子(exciton),从而使问题更复杂。
LR-TDDFT已经在计算小分子和团簇等孤立体系的低能激发谱中取得了成功。图2给出了使用独立粒子近似(IPA,紫色虚线),RPA(蓝色实线)和LR-TDDFT(绿色点-虚线)计算的Na8团簇的吸收谱与实验数据(红色散点)的对比。可见,无论是形状细节还是吸收边,LR-TDDFT的计算结果都明显优于IPA和RPA。我们可以用一个类Dyson屏蔽方程来联系IPA,RPA和LR-TDDFT
IPA即
它是最简单的情况,不仅没有交换关联核心fxc,也没有考虑库伦的核心,在IPA下介电函数在倒空间的矩阵元也只剩下ε00的对角元。它相当于不考虑局域场效应的RPA。真正的RPA是在长波近似下计算完整的εKK’-1,它只忽略了交换关联核心fxc。为了平衡精准度和计算量,中的module38使用的RPA实际上是介于IPA和RPA之间的一种方法,它计算了包含非对角元的ε00-1,而真正完整的RPA计算则包含在module59中,我们仅用其修正关联能,以得到更精确的总能。LR-TDDFT同时包含了交换关联核心和库伦核心,因此它具有最多的细节。不过,美中不足的是LR-TDDFT的精度依然受制于交换关联泛函的精度,尤其是对长程相互作用的描述。因此,LR-TDDFT在计算固体时并不理想。
图2分别使用IPA,RPA和LR-TDDFT计算Na8团簇的结果1此外,当所加的微扰较强时,非线性响应开始显现,LR-TDDFT也会失去作用。值得注意的是,“响应”本身就是个微扰概念,当所加的外场强到一定程度时,响应函数将再也不是之前所讲的形式,最典型的例子就是电容被足够强的电压击穿。
PWmat中利用LR-TDDFT计算光吸收的module
PWmat提供的module68可以使用LR-TDDFT计算光吸收。当计算对象为小分子或团簇时,module68比module38效果更好。
与PWmat内置的TDDFT计算不同,module68并没有使用JOB=TDDFT。相对的,我们只是求得了KS本征值和本征态,然后构造并求解了Casida’s方程组。用户只需要以下两个步骤即可实现:
1.做一个外场下的JOB=SCF得到KS本征态和本征值,需要用户自己使用可执行文件gen_confine_Vext.x得到IN.VEXT文件。加外场的目的是为了获得足够高的能量的本征态。为了保证后续TDDFT的计算结果,用户需要提供足够多的空态。用户可以在SCF时就设置一个很大的NUM_BAND,或者利用SCF的OUT.VR做一个设置了很大的NUM_BAND的JOB=NONSCF。具体采用哪种做法则取决于计算量和收敛效果。
2.使用并行计算程序calc_absorp_TDDFT_LR.x构造并求解Casida’s方程组。
计算结束后会得到absorp.TDDFT,它共有6列,如图3所示。该module支持HF,PBE0,HSE06和B3LYPS等计算的KS本征值和本征态。
图3输出文件absorp.TDDFT的内容以及每一列的意义以苯环为例,使用HSE06计算的吸收谱如图4所示。PWmat计算杂化泛函的速度非常快,相比其它同类型软件,用户可以用最少的机时获得准确的吸收边。
图4LR-TDDFT计算得到苯环分子的吸收谱相比于LR-TDDFT,PWmat最常用的TDDFT方法是realtimeTDDFT(rt-TDDFT)。这种TDDFT的应用范围很广。我们将在下一期介绍这种方法在计算光吸收中的应用。
1.Botti,S.;Schindlmayr,A.;Sole,R.D.;Reining,L.ReportsonProgressinPhysics,70,(3),-.
2.Kaur,J.;Ospadov,E.;Staroverov,V.N.JournalofChemicalTheoryandComputation,15,(9),-.
3.Gunnarsson,O.;Lundqvist,B.I.PhysicalReviewB,13,(10),-.
4.Gunnarsson,O.;Lundqvist,B.I.PhysicalReviewB,15,(12),-.
5.Perdew,J.P.;Levy,M.PhysicalReviewB,31,(10),-.
北京龙讯旷腾科技有限公司是成立于年的国家高新技术企业,是国内材料计算模拟工具软件研发创新的领导者,致力于开发满足“工业4.0”所需的原子精度材料研发Q-CAD(quantum-