量化编程（四）：周期体系的能带计算（理论篇）

nagami

量化编程(四)：周期体系的能带计算（理论篇）

文/Nagami   2016-1-27

       如何把求解隔离态分子的Hartree-Fock方程或者DFT方程推广到周期体相材料？这段时间查阅了一些资料，但收效甚微，没有得到一个令人满意的结果。现把了解到的内容做一番整理，一方面是希望熟悉这块的大咖能指点一二，另一方面是小编的凝聚态知识不怎么到家，如果下面的内容哪里有认识上错误，请尽情指出来，好改正

      下面，小编把分子的HF方程与固体的HF方程对应起来做一番简单的整理。

一）固体的周期性

晶格需要一个basis来表征其周期性，假设为a、b、c，其两两可以正交，也可不正交（会引发什么问题？），倒易空间的basis通过与a、b、c建立正交关系获得。

周期的平移特征通过集合R = {ua+vb+wc：u，v，w是整数}描述。周期代表单元记为Ω，体积为|Ω|，可以为单胞，超晶胞等，下面简称Cell

二）固体的能带

能带结构可从周期势的多体Schrodinger方程出发系统建立，详细的表述与Hilbert空间下的紧自伴算子的谱理论有关。见文献Reed和Simon[1975]。

三）固体HF方程的建立

3-1、晶体的Thermodynamic Limit

事实上小编查找了不少资料，但下面这本书的介绍和做法作为切入点让人比较满意：

【1】Dr. Mireille Defranceschi, Prof. Claude Le Bris (auth.)-Mathematical Models and Methods for Ab Initio Quantum Chemistry-Springer-Verlag Berlin Heidelberg (2000)

做法是：选取N个Cell的有限体系，计算一般有限分子体系的HF能量I_N，虽然无限体系的能量会发散，但I_N/N可能会趋近于一个有限值，可被认为是Single Cell的能量。

这引发了一系列问题，下面是书中的段落：

(i) Does there exist a limit for the energy per unit volume I_N/N when N goes to

infinity?

(ii) Does the minimizing electronic density rho_N approach a limit rho_inf(in a sense to

be made precise later) when N goes to infinity?

(iii) Does the limit density rho_inf have the same periodicity as the assumed periodicity

of the nuclei ?

这些问题表面上看过于基础化，但却很关键，而且相关讨论的数学文献也不少。

至少有两个问题已经引起了我们的注意：

1）无限周期体系的能量会发散（内在的），有意义的做法是考虑单个Cell的能量。

2）库伦势的无限叠加会发散。

第2个问题是来至于对离子晶体的静电能处理，即所谓的Ewald Summation技术，我们知道电中性的无限点阵电荷的静电能是条件收敛的，因此其绝对收敛级数会发散。关于Ewald求和的文献见：

【2】W. de Leeuw, J. W. Perram, and E. R. Smith (1980) Simulation of electrostatic systems in periodic boundary conditions. i. lattice sums and dielectric constants. Proc. R. Soc. Lond. A 373, pp. 27–56

3-2、Cell的能量计算

如何计算cell与其镜像相互作用的Hartree-Fock能，是接下去关键处理的问题。目前有两种方案

1）黄美纯老师的PPT课件提到的”电荷密度重整化“，见http://www.doc88.com/p-2866711033303.html，首先积分域虽然在Cell上，但受到的镜像Cell对其的库伦势是无限多叠加的，如果计算单个能量（核吸引，电子排斥能，核排斥能）不引入背景电荷去中和，单项能必然会发散（注意整体是电中性，是不会发散，条件收敛）。这种模型小编觉得还不够严密，物理上的做法。

2）数学上的做法是，利用Possion方程，建立Ewald周期势，将库伦势替换为Ewald势，下面详讲

3-3、Ewald势

经典HF能的公式如下，分别是动能，核吸引能、电子排斥能、HF交换能、核排斥能，积分区域是全域R3,

固体中的HF可做类比推广，

1）积分区域改为Unit Cell，只计算单胞能

2）库伦势用周期边界条件的Laplace算子的Green函数取代，注意到库伦势是全域Laplace算子的Green函数。

3）Bolch定理赋予固体新的量子数k，根据态叠加原理，需要计及。因为Bloch波矢k是连续周期变化，故只需对第一BZ区积分叠加。

4）周期性边界条件下，由于流守恒要求，需要满足相容性边界条件，Ewald函数需要引入背景电荷，做电中性化处理，规范化也很重要。

考虑以上条件，方程如下所示，关于Ewald势的更多的细节，见如下文献

【3】V. R. Saunders, C. Freyria Fava, R. Dovesi, L. Salasco and C. Roetti On the electrostatic potential in crystalline systems where the charge density is expanded in Gaussian Functions Molecular Physics, 77, 629-665, 1992

关于Ewald势的周期性，规范化等诸多性质见文献【3】，文献【3】还给出了Ewald函数的级数式，文献【1】也给出有对应的解析式，对于数值计算级数形式想必更好，如下

G的级数解形式做了收敛加速处理，见文献【2】。最终数值求解的方程形式如何，如下的Crystal小组的文献具有很好的参考：

【4】Dovesi, R., Orlando, R., Roetti, C., Pisani, C., Saunders, V.R.: periodic Hartree-Fock method and its implementation in the Crystal code. Phys. Stat. Sol. (b) 217, 63–88 (2000)

这篇文献给出了如下LCAO离散的Hartree-Fock矩阵元，但小编觉得，为了可做数值计算，应该把下式库伦势的叠加（求和M）用Ewald势取代，做类似的离散化。

这里的局域基可采用CGTO，但要改造成Bloch波形式，这样基函数就满足Bloch定理，求得的解自然也满足Bloch定理。积分域是Cell上的，周期性扩展使得可遍历到全域上。可自行推导验证。

密度泛函中的交换关联项由于是局域项，所以处理起来相对方便，积分域改为Cell上就可以，密度函数是周期的。但Becke分区的数值积分由于限制在Cell上，权函数需周期化，如此可扩展到R3的积分，也可采用合适的容积求积法，如下的文章是有用的

【5】Towler MD, Zupan A, Causa M (1996) Density functional theory in periodic systems using local Gaussian basis sets. Comput Phys Commun 98(1–2):181–205

关于K空间离散积分，MP方法和CC方法是众所周知的。

以上就是大致的内容了。基本的框架就是如此，小编理解浅陋，也想听听其他人的意见。

dindindindindin

好文章，学习了，

清微

老师您好，我有个困扰很久的问题：假如我们使用了周期性边界条件，建立一个立方盒子，包含100个水分子和1个质子，模拟质子转移过程，那么，1.这个晶胞中的H+和周围晶胞中的H+会产生排斥作用。如果盒子边长不太大的话，那么是否稳定？2.对于这样的体系，怎么计算其PH值，如果按照一个晶胞来算，其PH值为0.3，不合化学常理，该怎么考虑？