21. 化简 11：SCF 自洽迭代¶
1. 化简位置¶
12 步化简链的第 11 步。到这一步我们已经完成了：化简 1-5 把多体问题压缩到 Kohn-Sham 单粒子方程 + XC 泛函；化简 6-10 把 KS 方程离散化成平面波 + k 点网格 + 对称性上的有限矩阵特征值问题。
但还剩最后一个循环依赖没解开——这就是本步要处理的。
2. 上一步的困难：\(V_{eff}\) 与 \(\rho\) 的鸡生蛋¶
Kohn-Sham 方程形式上是一个单粒子本征问题：
其中有效势
Hartree 势 \(V_H[\rho](\vec{r}) = \int \frac{\rho(\vec{r}')}{|\vec{r}-\vec{r}'|} d\vec{r}'\) 和 XC 势 \(V_{xc}[\rho]\) 都是密度 \(\rho\) 的泛函。而 \(\rho\) 又是通过 KS 轨道求出的：
这是一个非线性耦合：\(\rho \to V_{eff}[\rho] \to H \to \phi_i \to \rho'\)。求解不能一步到位。
ML 视角：这就像 EM 算法里的隐变量——你需要参数才能算出隐变量的后验，需要后验才能更新参数。解法也类似：固定点迭代。
3. 引入的假设 / 解法¶
假设：该固定点存在且可达（对几乎所有常规体系成立）。
解法：Self-Consistent Field (SCF) 迭代 + 混合算法。
一次 SCF 迭代的完整流程（calculation='scf' 内部循环）：
- 初猜 \(\rho_0\)
- 原子密度叠加（superposition of atomic densities），或上一次计算的结果
- 构造 \(V_{eff}[\rho_{in}]\)
- 算 \(V_H\)（FFT 解泊松方程）、\(V_{xc}\)（逐点评估泛函）
- 构造 H 算符在平面波基组里的作用
- 注意是"作用"而非"显式矩阵"——基组维度可达 \(10^5\)-\(10^6\)，显式 H 根本装不下
- 解广义特征值问题 \(H \vec{c} = \varepsilon S \vec{c}\)
- QE 默认用 Davidson 对角化（
diagonalization='david'） - Davidson 只需要 H·ψ 的作用，属于 Krylov 子空间方法
- 类比：GMRES / Lanczos / CG，矩阵-向量积 O(N²) 而不是显式分解 O(N³)
- 填充占据态
- 绝缘体：按 Aufbau 原则，能量最低 N 个轨道填 1，其他填 0
- 金属：费米面附近用 smearing（Gaussian / Marzari-Vanderbilt / Methfessel-Paxton），把阶跃函数平滑成温度类似的分布
- 算新密度 \(\rho_{new}(\vec{r}) = \sum_i f_i |\phi_i(\vec{r})|^2\)
- 混合 \(\rho_{in}^{(n+1)} \leftarrow \text{Mixer}(\rho_{in}^{(n)}, \rho_{new}^{(n)}, \text{历史})\)
- 检查收敛
|ρ_new - ρ_old| < conv_thr（QE 实际用的是能量意义下的估计收敛性：estimated scf accuracy）- 或总能量变化
|E_new - E_old| < conv_thr - 若未收敛，回到第 2 步
4. 引入的新概念¶
4.1 Davidson 对角化¶
广义特征值问题 \(H\vec{c}=\varepsilon S\vec{c}\) 在平面波表示下维度 \(\sim 10^5\)。但我们只要最低的 \(N_{band}\) 个（\(\sim\) 价电子数），用 Davidson 算法：
- 维护一个子空间 \(V\)（初始含 N_band 个试向量）
- 扩展方向：\(H\vec{v}_i - \varepsilon_i S\vec{v}_i\) 的预条件残差
- 在 \(V\) 里做小规模 Rayleigh-Ritz（几百维），得当前特征值估计
- 迭代扩展到收敛
关键：只需要 H·ψ 的作用（在平面波基下是 FFT + 局部乘 + 投影），不构造完整 H 矩阵。复杂度 \(O(N_{band}^2 \cdot N_{PW} \cdot \log N_{PW})\)。
ML 类比：就像做 PCA 时的 power iteration / Lanczos，不需要算完整协方差矩阵。
4.2 密度混合 (mixing)¶
为什么要混合？如果直接 \(\rho_{in}^{(n+1)} = \rho_{new}^{(n)}\)，绝大多数体系会振荡或发散——电子气的介电响应很强，密度扰动被放大。
线性混合 (plain): $\(\rho_{in}^{(n+1)} = \alpha \rho_{new}^{(n)} + (1-\alpha) \rho_{in}^{(n)}\)$
其中 \(\alpha \in (0, 1]\) 就是 mixing_beta。它起的作用是低通滤波——抑制高频振荡。
Anderson / Broyden 混合：用历史若干步 \((\rho_{in}^{(k)}, \rho_{new}^{(k)} - \rho_{in}^{(k)})\) 估计一个近似的逆 Jacobian，把迭代从一阶收敛加速到准牛顿的超线性收敛。
ML 类比非常直接：
- 线性混合 = SGD with momentum（\(\alpha\) 小 ≈ 学习率小）
- Broyden / Anderson = quasi-Newton（类似 L-BFGS 里维护历史差分）
- mixing_ndim = 保留的历史步数，类似 L-BFGS 的 m
QE 默认用 Broyden 的一个变体（mixing_mode='plain' 其实已是 Broyden-风格 mixing，只是命名历史遗留）。
4.3 Smearing（金属特有）¶
费米面上的占据函数是阶跃 \(f(\varepsilon) = \Theta(\varepsilon_F - \varepsilon)\)，但在数值上：k 点附近的本征值会穿过 \(\varepsilon_F\)，造成占据数剧烈变化（0 跳到 1）→ 对 k 点采样极敏感 → 能量不收敛。
解法：用一个宽度 degauss (Ry) 的光滑函数替代阶跃：
- smearing='gaussian'：高斯
- smearing='mv' / 'cold'：Marzari-Vanderbilt（推荐）
- smearing='mp'：Methfessel-Paxton（可消除一阶误差）
物理上类比有限电子温度 \(T \sim\) degauss，但只是数值 trick，不代表真实温度。
4.4 SCF 不收敛的"难体系"¶
| 体系 | 为什么难 | 应对 |
|---|---|---|
| 金属 | 费米面态密度高，占据变化剧烈 | smearing + 较小 mixing_beta |
| 磁性（Fe, Ni, NiO） | 存在非磁解 / 反铁磁解多个局部极小；易塌到非磁解 | 显式 starting_magnetization |
| 强关联 (NiO, 过渡金属氧化物) | DFT 本身定性错误（把 Mott 绝缘体预测成金属） | 不是 SCF 的问题，需要 DFT+U、DMFT |
| 大晶胞 / 低维 | 长程电荷转移被混合算法误处理 | mixing_mode='local-TF' 或 'TF' |
5. 对应 QE 字段¶
&ELECTRONS 名字空间全部关于 SCF：
| 字段 | 含义 | 典型值 |
|---|---|---|
calculation='scf' |
触发单点 SCF（无原子移动） | — |
conv_thr |
收敛阈值（Ry，能量意义） | 1e-6（粗）→ 1e-8（常规）→ 1e-10（超严） |
mixing_beta |
线性混合系数 α | 0.7（绝缘体），0.3（金属），0.1（难体系） |
mixing_mode |
混合算法 | 'plain'（默认，含历史）/ 'TF' / 'local-TF' |
mixing_ndim |
混合历史长度 | 8（默认） |
electron_maxstep |
最大迭代步数 | 100（默认足够），磁性 200+ |
diagonalization |
对角化器 | 'david'（默认）/ 'cg'（更稳但慢） |
&SYSTEM 里与 SCF 相关的：
| 字段 | 含义 |
|---|---|
occupations |
'smearing' / 'fixed' / 'tetrahedra' |
smearing |
'mv' / 'gaussian' / 'mp' 等 |
degauss |
smearing 宽度（Ry），典型 0.01-0.02 |
nspin=2 |
打开自旋极化（磁性体必需） |
starting_magnetization(i) |
第 i 种原子初始磁矩（Bohr magneton 数） |
输出中的 SCF 监控¶
每轮迭代 QE 输出：
iteration # 3 ecut= 40.00 Ry beta=0.70
Davidson diagonalization with overlap
ethr = 1.00E-08, avg # of iterations = 2.3
total cpu time spent up to now is 1.8 secs
total energy = -15.84528512 Ry
estimated scf accuracy < 0.00000041 Ry
estimated scf accuracy：当前与自洽解的能量差估计，要降到conv_thr以下avg # of iterations：内层 Davidson 步数，如果一直 >20 说明对角化预条件不好- 最后
convergence has been achieved in N iterations表示成功
6. 对应 benchmark 条目¶
DFTBench 中涉及 SCF 收敛的字段：
questions/*.json里的输入模板有conv_thr参数questions/vc_relax_scf.json要求conv_thr = 1.0e-8（高于常规 scf，因为结构弛豫对力的精度要求高）- 输出
total energy→ 评分字段total_energy_ev_per_fu（每 formula unit 的总能量，eV 单位） - SCF 不收敛时 benchmark 直接失败（计 0 分）
见 tritonDFT-src/benchmark/questions/ 下 9 个任务模板。
7. ML 类比速查¶
| DFT/SCF 概念 | ML 对应 |
|---|---|
| SCF 迭代 | EM 算法（E 步：给定 V_eff 解对角化得轨道；M 步：由轨道算 ρ 更新 V_eff） |
| 线性混合 | SGD with momentum，mixing_beta ≈ 学习率 |
| Broyden/Anderson 混合 | quasi-Newton / L-BFGS，mixing_ndim = 历史长度 m |
| Davidson 对角化 | Krylov 方法（Lanczos / GMRES / CG），只需矩阵-向量积 |
conv_thr |
early stopping 的阈值 |
| SCF 不收敛 | 训练发散（梯度爆炸 / loss 振荡） |
| smearing | 类似 softmax 温度，消除不可导阶跃 |
starting_magnetization |
网络权重初始化（初始化不好会掉进差的 local min） |
特别要注意：SCF 是内层循环。下一步 BFGS 是外层循环。整个 vc-relax 是一个双层优化结构（像 MAML / bilevel optimization），本笔记处理内层。
8. 典型取值与常见坑¶
典型参数¶
| 体系类型 | conv_thr | mixing_beta | degauss | electron_maxstep |
|---|---|---|---|---|
| 半导体/绝缘体（Si, NaCl） | 1e-8 | 0.7 | — (fixed) | 100 |
| 简单金属（Al, Cu） | 1e-8 | 0.3-0.5 | 0.01 Ry | 100 |
| 磁性（Fe, Ni） | 1e-8 | 0.3 | 0.01 Ry | 200+ |
| 强关联（NiO, MnO） | 1e-7 | 0.1-0.3 | 0.01 Ry | 300+，常配 DFT+U |
常见坑¶
- SCF 不收敛 / 振荡
- 降
mixing_beta（0.7 → 0.3 → 0.1） - 增
mixing_ndim（8 → 20） - 换
mixing_mode='local-TF'（对长程电荷转移有效） - 增
electron_maxstep -
检查初始结构是否合理（原子距离过近会导致密度奇异）
-
磁性体塌到非磁解
- 必须显式设
starting_magnetization(1)=0.5（每种原子） -
反铁磁需要不同 site 给相反初始磁矩 → 用
nspin=2+ 人工 site 区分 -
金属 k 点不够 / degauss 不对
- 金属需要比绝缘体密得多的 k 网格（Si 8×8×8 够，Fe 可能要 16×16×16）
degauss太大 → 能量有系统偏差；太小 → 需要更多 k 点-
做
degauss收敛测试：外推到degauss→0 -
Davidson 卡住 /
avg # of iterations很大 - 初猜密度差 →
startingwfc='random'换成'atomic' -
对角化阈值
ethr自动调节，一般不用碰 -
Si SCF 实战参考（你之前跑过的）
- Si 体系在
ecutwfc=40, k=6×6×6, conv_thr=1e-8下约 10 轮收敛 - Ge 稍难，15-20 轮
- Fe 磁性体可能 50+ 轮
-
NiO DFT+U 可能要 100+ 轮
-
conv_thr太松导致力不准 — 这是下一步 BFGS 会踩的坑，见22-bfgs-relax.md。
下一步阅读¶
22-bfgs-relax.md— 化简 12：BFGS 结构弛豫（外层优化，内层就是本篇的 SCF）30-observables.md— SCF 收敛后能直接读出哪些可观测量