---
orphan: True
---

# 电子关联

```{card}
写于2025年3月4日

整理发布于2025年3月16日
^^^

之前整理了一下电子关联问题，也涉及到关联能，交换能等化学计算的重要问题。

在量子计算发展的早期，也就是当下时代，选择合适的化学问题来求解是重要的，强关联问题应当是这类问题的一个特征。
```

设非相对论电子哈密顿量为：

$$
H = \sum_{i} \left( -\frac{1}{2}\nabla_i^2 - \sum_A \frac{Z_A}{r_{iA}} \right) + \sum_{i<j} \frac{1}{r_{ij}}
$$

```{admonition} 相同表述形式：
---
class: seealso dropdown
---
可以参加电子哈密顿量从BO近似的来源[](/docs/量子化学/通用概念/au_boa.md#born-oppenheimer-approximation)
```

## 一、基本近似

任何量子化学方法均可分解为以下步骤：

### 基组截断 (Basis Set Truncation)

```{note}
:class: margin

电子薛定谔方程的维度过于大，不能严格求解，必须通过近似方法。这一系列近似方法就导致了计算精准度的缺失。如图所示，使用有限的基组只能在一定范围内逼近波函数，但总会有误差。越大的基组误差越小，但对算力消耗大。选择合适的基组，在 HF 和 DFT 方法中都很重要。

```

![](https://pic3.zhimg.com/v2-691e624b5ed9a1bd0556ae8f03a3cc98_1440w.jpg)

- 数学操作：将单电子轨道空间限制为有限维子空间 $ \text{Span}\{\chi_\mu\}_{\mu=1}^M $
- 后果：引入**基组误差**（Basis Set Error），属于技术性近似。

**例子: Cusp条件与基组问题**

在非相对论量子力学中，当两个电子靠近时（$r_{12} \to 0$），波函数需满足**Kato cusp条件**：

$$
\left. \frac{\partial \Psi}{\partial r_{12}} \right|_{r_{12}=0} = \frac{1}{2} \Psi(r_{12}=0),
$$
这一条件确保电子间的库仑排斥发散（$\frac{1}{r_{12}}$）被动能项精确抵消。然而，常用高斯型基组（如cc-pVDZ）在$r_{12}=0$处过于平滑，无法自然满足cusp，导致动态关联能被低估。  

![img](https://figures.semanticscholar.org/339be0adeb55deabc94f2412d6618aee7b0903a8/10-Figure2-1.png)

**改进方法**：  

- 显式相关方法（如F12）：在基组中引入显式依赖$r_{12}$的函数项。  
- 基组外推：通过增大基组尺寸逼近完备基极限。

### 波函数反对称化(Antisymmetrization)

```{note}
:class: margin

对电子这种费米子来说，严格的解应当只有反对称解。但是电子薛定谔方程的解存在对称解和它简并。因此，如果在求解时对称解和反对称解耦合，就会引起“交换能”。HF 在波函数形式上保证了反对称，而 DFT 没有，xc 泛函是包含交换能的。
```

- **HF**，通过构造Slater行列式 $ |\Phi\rangle = \mathcal{A}(\phi_1 \phi_2 \cdots \phi_N) $, 这是**完备**的。
- 严格满足Pauli原理，无信息损失。
- **KS波函数**：形式上为Slater行列式，满足反对称性，但**仅用于构造电子密度**，不直接参与能量计算。
- **交换关联势**：通过泛函$E_{xc}[\rho]$隐式包含交换与关联效应，无显式反对称项。例如：
  - LDA泛函：$E_{xc}^{\text{LDA}}[\rho] = \int \rho(\mathbf{r}) \epsilon_{xc}^{\text{hom}}(\rho(\mathbf{r})) d\mathbf{r}$，其中$\epsilon_{xc}^{\text{hom}}$为均匀电子气模型下的交换关联能密度。
  - 杂化泛函（如B3LYP）：混合HF交换与DFT交换，如$E_x^{\text{B3LYP}} = a E_x^{\text{HF}} + (1-a) E_x^{\text{DFT}}$。

### 单参考近似(Single-Reference Approximation)

- 数学操作：假设真实波函数可展开为单个参考态 $ |\Phi_0\rangle $ 及其激发态的线性组合：

$$
|\Psi\rangle = c_0 |\Phi_0\rangle + \sum_{ia} c_i^a |\Phi_i^a\rangle + \sum_{ijab} c_{ij}^{ab} |\Phi_{ij}^{ab}\rangle + \cdots
$$

- 后果：引入**关联能缺失**（Correlation Energy Error）。

以下表述来自 Garnet Chan 的 lecture note.

```{note}
:class: margin

在 HF/DFT 类方法中，使用 SCF 的计算手段，我们首先得到的是一个基态。这个基态包含了单电子近似引起的问题。为了消除单电子近似的误差，HF 里使用激发态来补偿，而 DFT 使用 xc 泛函来补偿。Garnet Chan 用纠缠熵来表述，实际上是一种量子计算的视角。单电子近似只能表示直积态。
```

![image-20250228130816744](./0316/image-20250228130816744.png)

![image-20250228130837906](./0316/image-20250228130837906.png)

## 二、关联(correlation)

从 HF 的视角出发，可以更清晰地看出关联之间的区别，我们假设通过 HF 已经选定了一个态做参考态，以此出发，通过它的激发态耦合出更低能量的态。

在给定基组中，如果所有激发态的线性组合被遍历，应该来说，这种搜索在这个子空间是完备的。不过这在计算上是不现实的。一般来讲，有两类在实际体系中常见的情况：

###  动态关联（Dynamic Correlation）

```{note}
:class: margin

一般来说，动态关联指的是参考态有效，但高阶激发态较多，总体影响不小。这类情况常见于分子的平衡键长位置。可以通过 CC 方法快速收敛。
```

- **定义**：在单参考态 $ |\Phi_0\rangle $ 的邻域内，可通过微扰展开或低阶激发态线性组合有效逼近的关联效应。

- **数学特征**：
  - **能量尺度**：激发态与参考态能量差较大（$ \Delta E = \epsilon_a + \epsilon_b - \epsilon_i - \epsilon_j \gg |\langle \Phi_{ij}^{ab} | H | \Phi_0 \rangle| $）
  - **波函数权重**：主组态权重 $ |c_0|^2 \approx 1 $，激发态权重快速衰减（如 $ |c_i^a| \sim 10^{-2} $）
  - **微扰收敛性**：满足Brillouin定理（单激发态贡献为零），MP2、CCSD(T)等微扰方法收敛。

- **物理根源**：电子瞬时位置涨落导致的库仑排斥修正（短程效应）。

### 静态关联（Static Correlation）

```{note}
:class: margin

静态关联指的是参考态失效。也可以认为是 HF 的基态和单激发态能量非常接近，也可以说单电子轨道排布可能是不对的。真实基态往往是多个不同排布纠缠形成的，其能量显著低于 HF 基态。这种情况下，要考虑多参考态方法，例如 CAS.
```

- **定义**：需要多个参考态非微扰叠加才能描述的关联效应，无法通过单参考激发态展开有效逼近。

- **数学特征**：
  - **能量简并**：存在多个Slater行列式 $ \{|\Phi_I\rangle\} $ 满足 $ |\langle \Phi_I | H | \Phi_J \rangle| \sim |E_I - E_J| $
  - **波函数权重**：多个组态权重相当（如 $ |c_1|^2 \approx |c_2|^2 \approx 0.5 $）
  - **微扰发散**：微扰展开失效（如H₂解离时MP2能量发散）。

- **物理根源**：量子多体系统的非定域纠缠，如化学键断裂、过渡金属d电子强关联。（长程关联）

---

## 三、哈密顿量矩阵结构的可视化区分

考虑组态相互作用（CI）的哈密顿矩阵：

$$
H_{\text{CI}} = 
\begin{pmatrix}
\langle \Phi_0 | H | \Phi_0 \rangle & \langle \Phi_0 | H | \Phi_1 \rangle & \cdots \\
\langle \Phi_1 | H | \Phi_0 \rangle & \langle \Phi_1 | H | \Phi_1 \rangle & \cdots \\
\vdots & \vdots & \ddots
\end{pmatrix}
$$

- **动态关联主导**：矩阵呈对角优势（对角元远大于非对角元），微扰展开有效。

  $$
  |H_{ii} - H_{jj}| \gg |H_{ij}| \quad (i \neq j)
  $$

- **静态关联主导**：存在子矩阵块，其中对角元接近，非对角元不可忽略。

  $$
  \begin{pmatrix}
  H_{11} & H_{12} \\
  H_{21} & H_{22}
  \end{pmatrix}, \quad |H_{11} - H_{22}| \sim |H_{12}|
  $$

---

## 四、方法设计的数学指导

### 1. 动态关联处理方法

- **微扰理论（MP2, CCSD(T)）**：  
  利用Møller-Plesset拆分 $ H = H_0 + V $，其中 $ H_0 = \sum_i \epsilon_i a_i^\dagger a_i $。  
  动态关联能通过二阶微扰项捕获：

  $$
  E_{\text{corr}}^{(2)} = \sum_{i<j, a<b} \frac{|\langle \Phi_0 | V | \Phi_{ij}^{ab} \rangle|^2}{\epsilon_i + \epsilon_j - \epsilon_a - \epsilon_b}
  $$

- **耦合簇（CC）**：  
  指数参数化 $ |\Psi\rangle = e^{T} |\Phi_0\rangle $，其中 $ T = T_1 + T_2 + \cdots $。  
  动态关联通过双激发算符 $ T_2 $ 主导。

### 2. 静态关联处理方法

- **多参考态展开（CASSCF）**：  
  选择活性空间 $ \text{Span}\{\phi_p\}_{p=1}^m $，构建多参考波函数：

  $$
  |\Psi_{\text{CAS}}\rangle = \sum_{I} c_I |\Phi_I^{\text{active}}\rangle
  $$
  其中 $ |\Phi_I^{\text{active}}\rangle $ 为活性空间内的全组态展开。

- **密度矩阵重整化群（DMRG）**：  
  通过矩阵乘积态（MPS）参数化多体波函数：

- $$
  |\Psi_{\text{DMRG}}\rangle = \sum_{\{n_i\}} A^{n_1} A^{n_2} \cdots A^{n_N} |n_1 n_2 \cdots n_N\rangle
  $$

  适用于强关联一维体系。

---

## 五、统一分类框架

```{list-table} 关联的分类
---
width: 80%
name: aaa
header-rows: 1
stub-columns: 1
---

* - 特征
  - 动态关联
  - 静态关联
* - 能量尺度
  - 激发态与参考态能量差大
  - 多个参考态能量接近
* - 波函数结构
  - 单参考基态, 系数 $c_0$ 主导
  - 单激发态系数 $c_i$ 接近基态
* - 哈密顿矩阵
  - 对角优势
  - 分块对角化（强非对角耦合）
* - 处理方法
  - 微扰理论、CC
  - CASSCF、DMRG
* - 物理实例
  - H₂基态、小分子振动能级
  - O₂解离、过渡金属配合物

```

## 六、Python示例

使用以下例子，可以更加真切的看到因为各种效应导致的能量差别:

### 基组截断

```{code-block} python
:linenos:

import numpy as np
import matplotlib.pyplot as plt
from pyscf import gto, scf, dft

# 定义键长范围和基组
bond_lengths = np.linspace(0.5, 3.0, 20)
basis_sets = ['sto-3g', '6-31g', 'cc-pvdz', 'cc-pvtz']

plt.figure(figsize=(10,6),dpi=300)
for basis in basis_sets:
    hf_energies = []
    dft_energies = []
    
    for r in bond_lengths:
        # 构建分子
        mol = gto.M(
            atom=f'H 0 0 0; H 0 0 {r}',
            basis=basis,
            spin=0,
            verbose=0
        )
        
        # HF计算
        mf_hf = scf.RHF(mol).run()
        hf_energies.append(mf_hf.e_tot)
        
        # DFT计算 (LDA)
        mf_dft = dft.RKS(mol, xc='lda,vwn').run()
        dft_energies.append(mf_dft.e_tot)
    
    # 绘制曲线
    plt.plot(bond_lengths, hf_energies, '--', label=f'Hartree-Fock/{basis.upper()}')
    plt.plot(bond_lengths, dft_energies, '-', label=f'DFT-LDA/{basis.upper()}')

plt.title('several BASIS-SET in HF&DFT calc')
plt.xlabel('bond length (Å)')
plt.ylabel('energy (Ha)')
plt.legend(bbox_to_anchor=(1.05, 1), loc='upper left')
plt.grid(True)
plt.show()
```

```{figure} ./0316/p1.png
---
figclass: margin-caption
alt: orbital_cutoff
name: orbital_cutoff
---
$H_2$ 解离曲线，在不同基组下。
```

这里我们模拟了不同基组下，HF 和 DFT 的结果差异。虚线是 HF 的结果，实线是 DFT. 可以看到，基组越大能量越低，但是到达一定程度也就够了。

### 交换能

```{admonition} 一个有趣的事实
---
class: important margin
---
这里可以看到常见泛函"B3LYP"的定义式，来自于若干个基础的泛函。从这里也可以感受到 DFT 方法简洁而高效。
```

```{code-block} python
:linenos:

# 1. 仅Hartree项（无交换）
class NoExchangeHF(scf.hf.RHF):
    def get_veff(self, mol=None, dm=None, dm_last=0, vhf_last=0, hermi=1):
        vj, _ = self.get_jk(mol, dm)
        return vj  # 只保留库伦项

# Note B3LYP functional mf.xc = 'HF*0.2 + .08*LDA + .72*B88, .81*LYP + .19*VWN'
# 2. 仅LDA交换（无关联）
class LDA_X(dft.rks.RKS):
    def __init__(self, mol):
        super().__init__(mol)
        self.xc = '0.0*LDA + 0.0*HF, 0.0*VWN'  # 只使用交换项

# 计算不同方法
methods = {
    'HF-NoX': [],
    'DFT-NoXC': [],
    'HF': [],
    'DFT': []
}

# energies = {name: [] for name in methods}
for r in bond_lengths:
    mol = gto.M(
        atom=f'H 0 0 0; H 0 0 {r}',
        basis='6-31g',
        spin=0,
        verbose=0
    )
    
    mf = NoExchangeHF(mol).run()
    methods['HF-NoX'].append(mf.e_tot)

    mdft = LDA_X(mol).run()
    methods['DFT-NoXC'].append(mdft.e_tot)

    mf = scf.RHF(mol).run()
    methods['HF'].append(mf.e_tot)

    mdft = dft.RKS(mol).run()
    methods['DFT'].append(mdft.e_tot)

# 绘图部分已略去
```

```{figure} ./0316/p2.png
---
figclass: margin-caption
alt: exchange_e
name: exchange_e
---
$H_2$ 解离曲线，在交换能被关闭的情况下。
```

这里，HF 和 DFT 在**同一基组**下，都被关闭了交换项。DFT 额外关闭了关联项，因为它不能单独关闭交换项。可以看到实际上如果不考虑交换关联，DFT 和 HF 的结果是一致，也就是直积态(单电子近似态)的能量。这是一个神奇的结果。而 HF 保留了交换项，能量大幅降低，DFT 考虑了交换和关联项，能量得以进一步降低。

### 关联能

```{code-block} python
:linenos:

from pyscf import mp, mcscf, fci

plt.figure(figsize=(10,6),dpi=300)
basis = 'sto3g'

methods = {
    'HF': [],
    'MP2': [],
    'CASSCF(2,2)': [],
    'DFT-LDA': [],
    'DFT-B3LYP': []
}

for r in bond_lengths:
    mol = gto.M(
        atom=f'H 0 0 0; H 0 0 {r}',
        basis=basis,
        spin=0,
        verbose=0
    )
    
    # HF
    mf = scf.RHF(mol).run()
    methods['HF'].append(mf.e_tot)
    
    # MP2
    pt = mp.MP2(mf).run()
    methods['MP2'].append(pt.e_tot)
    
    # FCI
    mc = mcscf.CASSCF(mf, 2, 2).run()
    methods['CASSCF(2,2)'].append(mc.e_tot)
    
    # DFT-LDA
    mf_lda = dft.RKS(mol, xc='lda,vwn').run()
    methods['DFT-LDA'].append(mf_lda.e_tot)
    
    # DFT-B3LYP
    mf_b3lyp = dft.RKS(mol, xc='b3lyp').run()
    methods['DFT-B3LYP'].append(mf_b3lyp.e_tot)

# 以HF为基准做相对能量
ref_energy = min(methods['HF'])
for method in methods:
    plt.plot(bond_lengths, np.array(methods[method]) - ref_energy, 
             label=method, marker='o')

# 绘图部分已略去
```

```{figure} ./0316/p3.png
---
figclass: margin-caption
alt: corr
name: corr
---
$H_2$ 解离曲线，在不同的方法下，可以显示其关联效应。
```

这里，HF 和 DFT 以及其升级方法，在**同一基组**下进行了计算。其中平衡位置 ($ \sim 0.74 \mathring {\mathrm A} $) 主要是动态关联，HF 和 LDA 处理的都不好，但 MP2 一定程度上处理了这个问题。多参考，或者说静态关联区域在 $1.5\sim 3 \mathring {\mathrm A} $ 左右，这里就需要 CASSCF 来处理，B3LYP处理的也不是很好。另一方面要注意 DFT 方法不是严格变分(因此猜测出来的泛函不来自于电子薛定谔方程)的，因此在平衡位置能量低于 CASSCF，这是错误的。