电子结构计算简介

电子结构计算简介

求解薛定谔方程是定义在 ${\mathbb{R}}^{3N}$ 高维线性特征值问题. 除极少数体系外, 无法直接求解. 因此, 人们需要寻找与其等价或简化的可计算模型. 这些等价或简化的可计算模型中最典型的是 Hartree-Fock 方程和 Kohn-Sham 方程. 求解这两种方程，其对应了波函数方法和密度泛函理论方法。前者是量子化学领域的主流，后者在凝聚态物理、材料计算等领域广泛使用。

问题背景

基本概念

电子结构计算，也称作分子模拟或计算量子化学, 是通过对微观尺度物质的模拟计算来理解其性质和规律. 几乎所有的电子结构计算问题都是非凸的. 电子结构计算的理论依据是奠定于二十世纪初的量子力学理论. 非相对论量子力学最流行的表述形式是薛定谔的波动力学形式, 它的核心是波函数的概念及其运动方程薛定谔方程. 电子结构计算本质上就是通过求解薛定谔方程来预测物质的性质和规律.

薛定谔方程是描述微观物质世界的典型的第一原理数学模型. 薛定谔方程通常被认为是描述电子结构的 “第一原理 (first principles)”. 求解薛定谔方程，根据是否使用经验参数和实验数据主要分为三类: 第一原理方法、半经验方法和经验方法. 相比较而言, 第一原理方法是最精确的计算方法, 能对不同的体系都得到十分精确可信的结果, 但是它的计算量很大, 对大规模问题需要非常长的计算时间. 半经验和经验的方法虽然很好地降低了计算量, 而且容易实现. 但是这两类计算依赖于特定参数的选定, 往往只是对某种特定的电子体系有比较精确的刻画, 而不能对各种系统都得到同样好的结果.

基于第一原理的电子结构计算方法主要包含波函数方法, 密度泛函理论, 量子蒙特卡罗, 约化密度矩阵方法等, 其中使用最广泛的是波函数方法和密度泛函理论. 波函数方法和密度泛函理论通过近似在保证一定精度的同时很好地降低了计算复杂度, 但是也带来相应的代价: 非线性. 使用波函数方法和密度泛函理论近似后得到的问题本质上都是带约束的非线性极小化问题, 它们对应的 Euler-Lagrange 方程是一类非线性方程, 如 Hatree-Fock 方程和 Kohn-Sham 方程. 把原来的线性问题转化成非线性问题给我们带来新的问题. 一言以蔽之, 电子结构可通过包括 Hartree-Fock 方程和 Kohn-Sham 方程在内的非线性特征值模型(或其等价形式)来刻画. 求解这些非线性特征值问题既可通过反复求解线性化后的线性特征值问题来实现, 也可转化为最优化问题计算.

密度泛函理论方法

由于 (定态的) 薛定谔方程是一高维特征值方程, 难以直接求解, 人们须寻求其简化或等价模型. 电子密度泛函理论则是其中最成功的代表. 通过密度泛函理论, 人们将高维薛定谔方程转化为三维空间上的 Kohn-Sham 方程. 这样, 基于密度泛函理论的第一原理电子结构计算的核心就转化为求解 Kohn-Sham 方程.

电子结构计算方法研究涉及离散方法、非线性迭代法、代数特征值求解器以及相应的程序实现技术等方面.

现有的非线性特征值问题(或其等价形式)的离散格式大体上可分为三类: 倒空间方法、局部基集法以及实空间方法. 这些方法各有其优势和不足, 其中倒空间方法和局部基集法研究相对成熟, 而实空间方法尚处于发展阶段. 实空间方法主要包括有限差分方法、小波方法、有限体方法与有限元方法等. 实空间方法具有很好的局部性, 适应并行计算, 在高性能计算机上实现时能极小全局通信时间, 近年来越来越受到人们的重视. 有限元方法基函数构造简单、基底完备, 易实现自适应计算. 无论用什么方法来离散, 都得到有限维非线性算子的特征值问题或其等价形式. 其数值解收敛与否以及收敛快慢程度不仅取决于数值方法本身而且还依赖于分子体系构成.

波函数方法

当处理小规模的体系, 人们主要关心计算精度而计算时间次之的时候, 波函数方法是非常成功有效地方法. 波函数方法的原型是 Hartree-Fock 近似, 它将多电子体系的电子函数用单个电子波函数的 Slater 行列式来表示, 从而导出 Hartree-Fock 方程进行计算. 在 Hartree-Fock 方法的基础上, 人们发展出了很多波函数方法, 使得计算精度进一步提高. 这些方法包括组态相互作用方法(configuration interaction, 简称 CI), 耦合簇方法 (coupled cluster, 简称 CC), MøllerPlesset 微扰方法 (MP2, MP3, MP4 等), 多组态自洽场方法 (multi-configuration self-consistent field, 简称 MCSCF), 多组态双激发组态相互作用方法 (MRDCI), 二次组态相互作用方法 (QCI) 等. 由于这些方法是在 Hartree-Fock 方法基础上进一 步改进变量空间和近似表达式, 因此也被称为后 Hartree-Fock 方法.

Hartree-Fock 方程形式上为本征方程, 但是福克算符中的库仑算符和交换算符都与分子轨道有关, 因此只能够通过自洽迭代的方法近似求解, 即哈特里-福克自洽场（HF-SCF）方法. HF-SCF 方法是组态相互作用方法、多体微扰理论、半经验量子化学计算等现代量子化学计算方法的基础.

最基础的 Hartree-Fock 算法

薛定谔方程

我们考虑由 $M$ 个原子核与 $N$ 个电子构成的体系, { R j : j = 1 , ⋯   , M } \left\{R_{j}: j=1, \cdots, M\right\} {Rj​:j=1,⋯,M} 是原子核的位置, { Z j : j = 1 , ⋯   , M } \left\{Z_{j}: j=1, \cdots, M\right\} {Zj​:j=1,⋯,M} 是原子核电荷数, { r i : i = 1 , ⋯   , N } \left\{r_{i}: i=1, \cdots, N\right\} {ri​:i=1,⋯,N} 是电子的位置. 描述这样一个多体系统的定态 (不依赖于时间变量的) 薛定谔方程为，
$$\mathcal{H}\psi =E\psi$$
方程中 $\psi \left(r_1,\cdots ,r_N\right)$ 为电子波函数, $E$ 称为能级, $\mathcal{H}$ 是描述电子的 Hamilton 量.
$$\mathcal{H}=-\sum _{i=1}^N\frac{{\hslash }^2}{2m_e}{\Delta }_{r_i}-\sum _{i=1}^N\sum _{j=1}^M\frac{Z_j{e}^2}{\left|R_j-r_i\right|}+\sum _{1\le i<j\le N}\frac{e^2}{\left|r_i-r_j\right|}$$

其中 $\hslash$ 是 Planck 常数除以 $2\pi ,e$ 是电子的电荷, $m_e$ 是电子的质量. 上式中的 三项分别代表电子的动能算符, 电子和原子核之间的相互作用, 以及电子之间的 相互排斥作用. 如果我们全采用原子单位, 即: $e=1;m_e=1;\hslash =1$, 长度用 bohr 作单位 $(1\mathrm{b}\mathrm{o}\mathrm{h}\mathrm{r}=0.529177249A)$, 能量单位为 hartree $(1$ hartree $=2$ rydbergs $=27.11\mathrm{e}\mathrm{V})$, 则上式变为
$$\mathcal{H}=-\sum _{i=1}^N\frac{1}{2}{\Delta }_{r_i}-\sum _{i=1}^N\sum _{j=1}^M\frac{Z_j}{\left|R_j-r_i\right|}+\sum _{1\le i<j\le N}\frac{1}{\left|r_i-r_j\right|}$$
在绝大部分情况下, 物理或化学系统处于最稳定的状态, 即能量最低的态, 称之为系统的基态, 即我们要求解薛定谔方程的最低能级. 所以要获得系统的基态等价于求解一个极小化问题
$$E_0=\inf \left\{⟨\psi ,\mathcal{H}\psi ⟩:\psi \in L^2\left({\mathbb{R}}^{3N}\right)\bigcap \underset{i=1}{\overset{N}{\bigwedge }}{H}^1\left({\mathbb{R}}^3\right),\Vert \psi {\Vert }_{L^2}=1\right\}$$
其中求得的最小值 $E_0$ 极为系统的基态能量, 符号 $\bigwedge$ 表示张量积. 定态薛定谔方程也可以看成极小化问题上述问题对应的 Euler-Lagrange 方程.

Hartree-Fock 方法

Hartree-Fock 方法是一种简单的第一原理电子结构计算方法. 它的本质思想是把问题
$$E_0=\inf \left\{⟨\psi ,\mathcal{H}\psi ⟩:\psi \in L^2\left({\mathbb{R}}^{3N}\right)\bigcap \underset{i=1}{\overset{N}{\bigwedge }}{H}^1\left({\mathbb{R}}^3\right),\Vert \psi {\Vert }_{L^2}=1\right\}$$
中的变量空间 $L^2\left({\mathbb{R}}^{3N}\right)\bigcap {\bigwedge }_{i=1}^N{H}^1\left({\mathbb{R}}^3\right)$ 限制到一个子空间 $S_N$ 上, 其中的波函数可以表示成由 $N$ 个单电子波函数构成的行列式, 即
$$S_N=\left\{\psi =\frac{1}{\sqrt{N!}}\det \left({\varphi }_i\right):{\varphi }_i\in H^1\left({\mathbb{R}}^3\right),{\int }_{{\mathbb{R}}^3}{\varphi }_i{\varphi }_j={\delta }_{ij},1\le i,j\le N\right\}$$
我们可以把 $S_N$ 中的波函数写成更具体的 Slater 行列式的形式
$$\psi =\frac{1}{\sqrt{N!}}\det \left({\varphi }_i\right)=\frac{1}{\sqrt{N!}}\left|\begin{array}{cccc}{\varphi }_1\left(r_1\right)& {\varphi }_1\left(r_2\right)& \cdots & {\varphi }_1\left(r_N\right)\\ {\varphi }_2\left(r_1\right)& {\varphi }_2\left(r_2\right)& \cdots & {\varphi }_2\left(r_N\right)\\ ⋮& ⋮& \ddots & ⋮\\ {\varphi }_N\left(r_1\right)& {\varphi }_N\left(r_2\right)& \cdots & {\varphi }_N\left(r_N\right)\end{array}\right|.$$
所以在 Hartree-Fock 近似下, 我们是求解如下极小化问题,
E 0 H F = inf ⁡ { ⟨ ψ , H ψ ⟩ : ψ ∈ S N } E_{0}^{H F}=\inf \left\{\langle\psi, \mathcal{H} \psi\rangle: \psi \in S_{N}\right\} E0HF​=inf{⟨ψ,Hψ⟩:ψ∈SN​}

将系统的 Hamilton 量与波函数的表达式代入, 我们可以得到 Hartree-Fock 方法更为具体的表示形式. 记

$$\rho (r)=\sum _{i=1}^N|\varphi (r){|}^2,\tau \left(r,r^{\prime }\right)=\sum _{i=1}^N\varphi (r)\varphi \left(r^{\prime }\right)$$

其中函数 $\rho$ 与 $\tau$ 分别称为波函数 $\psi$ 的密度和密度矩阵. 我们还简记
$$V_{ext}(r)=-\sum _{j=1}^M\frac{Z_j}{\left|R_j-r\right|}$$
于是，HF 极小化问题就等价于
$$E_0^{HF}=\inf \left\{E^{HF}\left({\varphi }_1,\cdots ,{\varphi }_N\right):{\varphi }_i\in H^1\left({\mathbb{R}}^3\right),{\int }_{{\mathbb{R}}^3}{\varphi }_i{\varphi }_j={\delta }_{ij},1\le i,j\le N\right\}$$
其中，能量的定义为,
$$\begin{array}{rl}{E}^{HF}\left({\varphi }_1,\cdots ,{\varphi }_N\right)=& \sum _{i=1}^N{\int }_{{\mathbb{R}}^3}{\left|\nabla {\varphi }_i\right|}^2+{\int }_{{\mathbb{R}}^3}{V}_{ext}\rho +\frac{1}{2}{\int }_{{\mathbb{R}}^3}{\int }_{{\mathbb{R}}^3}\frac{\rho (r)\rho \left(r^{\prime }\right)}{\left|r-r^{\prime }\right|}drd{r}^{\prime }\\ & -\frac{1}{2}{\int }_{{\mathbb{R}}^3}{\int }_{{\mathbb{R}}^3}\frac{{\left|\tau \left(r,r^{\prime }\right)\right|}^2}{\left|r-r^{\prime }\right|}drd{r}^{\prime }\end{array}$$
上述极小化问题对应的 Euler-Lagrange 方程即为 Hartree-Fock 方程,

$$\{\begin{array}{l}-\frac{1}{2}\Delta {\varphi }_i+V_{ext}{\varphi }_i+\left({\sum }_{j=1}^N{\left|{\varphi }_j\right|}^2\ast \frac{1}{|r|}\right){\varphi }_i-{\sum }_{j=1}^N\left({\varphi }_i{\varphi }_j\ast \frac{1}{|r|}\right){\varphi }_j={\lambda }_i{\varphi }_{i}i=1,\cdots ,N\\ {\int }_{{\mathbb{R}}^3}{\varphi }_i{\varphi }_j={\delta }_{ij}\end{array}$$

求解该方程有自洽场 HF-SCF 方法等，本质就是简单的接待。我们注意到 $L^2\left({\mathbb{R}}^{3N}\right)\cap {\bigwedge }_{i=1}^N{H}^1\left({\mathbb{R}}^3\right)$ 中的任何波函数都可以表示成 $S_N$ 中 (有 限或无限个) 波函数的线性组合. 扩大变量空间 (比如做 Slater 行列式的线性组 合), 使之尽可能接近 $L^2\left({\mathbb{R}}^{3N}\right)\cap {\bigwedge }_{i=1}^N{H}^1\left({\mathbb{R}}^3\right)$, 可以很大地提高第一原理计算的精 度, 这就是所谓的后 Hartree-Fock 方法. 在 Hartree-Fock 方法的基础上发展起来 的这类方法包括组态相互作用方法 (configuration interaction, 简称 CI), 耦合簇方 法 (coupled cluster, 简称 CC), Møller-Plesset 微扰方法 (MP2, MP3, MP4 等), 多组 态自洽场方法 (multi-configuration self-consistent field, 简称 MCSCF), 多组态双激发组态相互作用方法 (MRDCI), 二次组态相互作用方法 (QCI) 等 [108]. 这些后 Hartree-Fock 方法和 Hartree-Fock 方法统称波函数方法. 波函数方法虽然计算精 度很高, 但是计算复杂度也很大 (在 $\mathcal{O}\left(N^4\right)$ 以上). 在处理小规模体系, 主要关注 计算精度而计算时间次之的时候, 人们更多地使用波函数方法.

数学的人和化学的人用的符号可能不太一样，但是大体表示的是相同的意思。化学的写法，参考量子化学入门读物：SZABO 的 modern quantum chemistry。

软件包

CP2K是一款著名的从头算分子动力学软件.它是由马克斯-普朗克研究中心早在2000年发起的一项用于固体物理研究的项目, 代码主要使用Fortran 编写.现在它已转由苏黎世ETH和苏黎世大学维护, 成为一个开源的项目, 遵从GPL协议, 用户可以从其官方网站下载到源代码.

CP2K基于密度泛函理论(DFT), 它使用混合的高斯平面波近似(GPW)以及多粒子势, 可以计算更大的体系.支持的理论层次包括DFTB, LDA, GGA, MP2, RPA, 半经验方法(AM1, PM3, PM6, RM1, MNDO)和经典力场(AMBER, CHARMM)等.主要用来研究计算固态、液体、分子和生物体系的性质, 如量子点的结构计算、材料表面结构驰豫、相变过程(例如冰的融解、液态水的模拟)等等, 所涉及的研究领域均是与人类的生产、生活密切相关的.

CP2K官方GitHub仓库（https://github.com/cp2k/cp2k）.

中科院戴小英课题组基于科学与工程计算国家重点实验室开发的并行自适应网格平台 PHG (parallel hierarchical grid)5), 研制了一套第一原理实空间并行自适应计算程序 RealSPACES (real space parallel adaptive calculation of electronicstructure). 经过十几年的发展和完善, 他们的程序功能在逐步健全, 计算效率也在逐步提高. 特别要指出的是, 现有的第一原理计算商业或开源软件由于受基函数影响, 通常要么只能进行赝势计算 (如基于平面波离散的软件), 要么只能进行全势计算 (如基于原子轨道基函数方法离散的软件). 与之不同的是, RealSPACES 没有这方面的限制. 他们的程序既能进行赝势计算, 也能方便地进行全势计算. 没开源！

很多其他不列举…