Skip to content

Commit f496615

Browse files
committed
feat(deltap): per-atom Wilson loop decomposition with Born charge tests
Code changes: - Fix M^I conjugation bug: swap conj from left U(k_j) to right U(k_{j+1}) Correct formula: M^I = U^I_D(k_j) * O * U^{I\dagger}_D(k_{j+1}) - Add branch persistence: load_branch()/save_branch() read/write deltap_branch.dat for cross-SCF phase smoothness - Implement per-atom SVD with singular-value truncation (k_eff < nocc) - Implement global SVD + per-atom weight decomposition (exact sum rule) - Add SVD Procrustes alignment across k-points for gauge continuity - Add hybrid Wilson-loop total + Berry-connection proportions Design documents: - DeltaSpin_vs_DeltaP_Essential_Difficulties.md: five-level difficulty analysis - Wilson_Loop_Per_Atom_Decomposition.md: per-atom Wilson loop math - SMO_Projected_Wannier_Evaluation.md: evaluation methodology - DeltaP_Incremental_Design.md: original incremental design Tests: - 18_LCAO_DELTAP_BTO: update KPT to 10x10x10, add deltap_method - 19_LCAO_DELTAP_SI: update KPT to 10x10x10, add deltap_method - 20_LCAO_BTO_BORN: Born effective charge test (equilibrium + 5 atoms x +/-disp) Report: docs/superpowers/specs/2026-06-25-deltap-wilson-per-atom-report.md Documents 4 implementation attempts, Z* results, and key findings: - nproj=16 > nocc=15 causes Wilson loop degeneracy - Berry connection (trace) != Berry phase (det) for finite dk - Single k-string insufficient (need multi-string averaging)
1 parent ae4b531 commit f496615

44 files changed

Lines changed: 2633 additions & 103 deletions

File tree

Some content is hidden

Large Commits have some content hidden by default. Use the searchbox below for content that may be hidden.

DeltaP_Incremental_Design.md

Lines changed: 623 additions & 0 deletions
Large diffs are not rendered by default.
Lines changed: 205 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,205 @@
1+
# 从DeltaSpin到DeltaP:为什么极化约束比磁矩约束困难得多
2+
3+
> **问题**: 对比磁矩约束(DeltaSpin)与极化约束(DeltaP),分析复刻DeltaSpin方法到DeltaP的**本质困难**
4+
> **参考**: DeltaSpin论文 (Zheng et al., 2025), DeltaSpin实空间投影Note
5+
6+
---
7+
8+
## 1. DeltaSpin为什么能精确计算自旋相互作用
9+
10+
### 1.1 原子磁矩的局域定义
11+
12+
在DeltaSpin中,原子磁矩通过SMO投影从实空间密度矩阵直接计算:
13+
14+
$$\boxed{M^I_p = \sum_{\sigma\sigma'} \sum_{lmm'} \sigma^p_{\sigma\sigma'} \, n^{\sigma\sigma'}_{Ilmm'} \, \delta_{mm'}}$$
15+
16+
$$n^{\sigma\sigma'}_{Ilmm'} = \sum_{\mathbf{R}} \sum_{\mu\nu} \rho^{\sigma\sigma'}_{\mu\nu}(\mathbf{R}) \, \langle\phi^0_\mu|\alpha^I_{lm}\rangle\langle\alpha^I_{lm'}|\phi^{\mathbf{R}}_\nu\rangle$$
17+
18+
**这是DeltaSpin一切优势的根源**:$M^I$是一个**实空间的、局域的、规范不变的**物理量。
19+
20+
### 1.2 为什么它是精确的
21+
22+
逐条分析:
23+
24+
**(a) 局域性**:$M^I$仅依赖原子I附近的SMO投影。SMO本身的局域半径使每个原子的磁矩计算只涉及有限范围的实空间密度矩阵——**没有k空间积分**
25+
26+
**(b) 规范不变性**:$n^{\sigma\sigma'}_{Ilmm'}$ 是密度矩阵的二次型 $\rho \cdot |\alpha\rangle\langle\alpha|$ 的迹。在波函数规范变换 $\psi \to e^{i\varphi}\psi$ 下,$\rho_{\mu\nu} = \sum_n f_n C^*_{n\mu}C_{n\nu}$ 是不变的——相位 $e^{i\varphi}$ 在乘积 $C^*C$ 中消去。因此 $M^I$ **天然规范不变**
27+
28+
**(c) Sum Rule的精确性**:$\sum_I M^I = M^{\text{total}}$ ——当SMO集在自旋子空间上完备时严格成立。这里的"完备"只需要SMO能展开占据态的自旋密度——一个远比展开整个Hilbert空间宽松的条件。
29+
30+
**(d) 无模歧义**:$M^I$ 是**实空间向量的模**。它不需要取arg,不对2π取模。不存在类似极化量子的歧义。
31+
32+
### 1.3 H^λ的简洁结构
33+
34+
DeltaSpin的H^λ拥有最简洁的数学形式:
35+
36+
$$H^{\lambda,\sigma\sigma'}_{\mu\nu}(\mathbf{R}) = \sum_I f(I,\sigma\sigma') \sum_{lm} \langle\phi^0_\mu|\alpha^I_{lm}\rangle\langle\alpha^I_{lm}|\phi^{\mathbf{R}}_\nu\rangle$$
37+
38+
预存储 $H^{\text{pre},I}_{\mu\nu}(\mathbf{R}) = \sum_{lm} \langle\phi^0_\mu|\alpha^I_{lm}\rangle\langle\alpha^I_{lm}|\phi^{\mathbf{R}}_\nu\rangle$ 后:
39+
40+
$$\boxed{H^{\lambda} = \sum_I f(I) \cdot H^{\text{pre},I}}$$
41+
42+
**这是DeltaSpin H^λ的全部**——一个实数矩阵的线性组合,与k无关,与波函数无关。
43+
44+
而且 $H^{\text{pre},I}$ 的结构**恰好与DFT+U完全相同**
45+
46+
$$H^{+U,\sigma\sigma'}_{\mu\nu}(\mathbf{R}) = \sum_{I,mm'} \langle\phi^0_\mu|\alpha^I_m\rangle \, U\left(\frac{1}{2}\delta^{\sigma\sigma'}_{mm'} - n^{\sigma\sigma'}_{I,mm'}\right) \, \langle\alpha^I_{m'}|\phi^{\mathbf{R}}_\nu\rangle$$
47+
48+
DeltaSpin和DFT+U共享同一个SMO投影框架,共享同一个HContainer预存储,共享同一套力/应力修正公式。这是**零额外工程成本的约束方法**
49+
50+
---
51+
52+
## 2. 极化约束:五个层面的本质困难
53+
54+
当我们试图将DeltaSpin的框架复刻到极化约束时,在五个递进的层面上遭遇了根本性障碍。
55+
56+
---
57+
58+
### 困难层级一:物理量定义的全局性(物理层面)
59+
60+
**磁矩**
61+
$$M^I = \int_{\text{原子I附近}} \mathbf{m}(\mathbf{r}) \, d^3r$$
62+
积分区域是**实空间的一个局域球**。$M^I$ 是一个局域可加的物理量——体系的总磁矩等于各原子磁矩之和。
63+
64+
**极化**
65+
$$P^{\text{el}} = -\frac{ie}{(2\pi)^3} \sum_n \int_{\text{BZ}} d\mathbf{k} \, \langle u_{n\mathbf{k}}|\nabla_{\mathbf{k}}|u_{n\mathbf{k}}\rangle$$
66+
这是**k空间全BZ的积分**。$P$ 是一个全局量——它不由任何实空间局域积分给出。
67+
68+
**关键推论**:不存在"极化密度"$p(\mathbf{r})$使得 $P = \int p(\mathbf{r}) d^3r$ 且 $p(\mathbf{r})$ 可以按原子分区。现代极化理论(Resta, 1994; King-Smith & Vanderbilt, 1993)的核心结论就是:**极化是Berry phase,不是电荷密度的偶极矩**
69+
70+
这意味着DeltaSpin的路径——用实空间密度矩阵的原子分区来定义 $M^I$ ——在极化问题中**没有对应的操作**。P^I不能写成 $\int_{I} (\text{某局域量}) \, d^3r$ 的形式。
71+
72+
---
73+
74+
### 困难层级二:规范依赖性的本质差异(数学层面)
75+
76+
这是决定性的对比:
77+
78+
| | DeltaSpin 磁矩 M^I | DeltaP 极化 P^I |
79+
|---|---|---|
80+
| 对波函数的依赖 | $\rho = CC^\dagger$(相位消去) | $A = \langle\psi\|\partial_k\psi\rangle$(相位不消去) |
81+
| 规范变换 $C \to Ce^{i\varphi}$ | 不变 ✓ | $\partial_k \to \partial_k + i\partial_k\varphi$ ✗ |
82+
| 天然规范不变 |||
83+
84+
**数学根源**:$M^I$ 是密度矩阵 $\rho$ 的泛函,而 $\rho$ 是规范不变的。$P^I$ 涉及波函数的k导数 $\partial_k C$,这是规范依赖的。
85+
86+
**直接后果**:DeltaSpin不需要任何规范固定。独立SCF计算的波函数可以带有任意随机相位,$M^I$ 的计算结果不变。DeltaP必须解决规范固定问题——这直接引出了我们之前分析的所有连续性困难。
87+
88+
---
89+
90+
### 困难层级三:算符结构的代际差异(算法层面)
91+
92+
这是H^λ实现代价的根本差异:
93+
94+
**DeltaSpin H^λ**:$\propto |\alpha\rangle\langle\alpha|$ ——实空间投影算符
95+
96+
**DeltaP H^λ**:$\propto |\partial_k\alpha\rangle\langle\alpha| + |\alpha\rangle\partial_k\langle\alpha|$ ——k空间微分算符
97+
98+
两种算符结构的对比:
99+
100+
| 属性 | DeltaSpin H^λ | DeltaP H^λ |
101+
|------|:---:|:---:|
102+
| 实空间局域性 | ✅ 完全局域 | ⚠ 需Fourier变换到k空间 |
103+
| k点独立性 | ✅ 各k点独立 | ⚠ 涉及∂_k,连接相邻k点 |
104+
| 预存储 | ✅ HContainer(R) | ⚠ HContainer_grad(k) |
105+
| 与现有代码同构 | ✅ ≡ DFT+U | ❌ 无现有代码可复用 |
106+
| 最小实现代价 | 已有 | ~50行(仅Part A) |
107+
108+
**核心差异**:DeltaSpin的 $|\alpha\rangle\langle\alpha|$ 是**零阶算符**——它只"数"投影。DeltaP的 $|\partial_k\alpha\rangle\langle\alpha|$ 是**一阶算符**——它需要波函数的k导数信息。这是两个不同代数阶次的算符,计算代价天然不同。
109+
110+
---
111+
112+
### 困难层级四:响应矩阵的维度差异(数值层面)
113+
114+
**DeltaSpin的Jacobian** $\partial M^I / \partial \lambda^J$:
115+
- 对角占优(最近邻交叉项~30%)
116+
- 可用子空间微扰解析估计
117+
- CG/Broyden收敛稳定(DeltaSpin已验证)
118+
119+
**DeltaP的Jacobian** $\partial P^I / \partial \lambda^J$:
120+
- 同样有交叉项(~20-40%)
121+
-**不能**用子空间微扰简单估计——因为$P^I$需要k空间积分,而子空间微扰只能估计单个k点的响应。k空间积分将每个k点的误差累积。
122+
- Berry connection方案中,$\partial P^I/\partial\lambda^J$ 受到与 $\partial P^I/\partial\tau^J$ 相同的1/Δk放大效应
123+
- Wilson loop方案中,$\partial P^I/\partial\lambda^J$ 需要通过SVD导数链(5个环节),计算代价极高
124+
125+
---
126+
127+
### 困难层级五:2π歧义对约束框架的根本性冲击(理论层面)
128+
129+
这是DeltaP面临的最独特的困难——DeltaSpin完全没有对应的问题。
130+
131+
**Berry phase的模2π歧义**:$P^I$ 只能确定到模 $eR/\Omega$。
132+
133+
**这对约束框架意味着什么?**
134+
135+
在λ内循环中,算法试图找到λ使得 $P^I(\lambda) = P^{I,\text{target}}$。但如果 $P^I$ 在λ空间是多值的(不同分支差一个极化量子),那么:
136+
137+
1. **解的非唯一性**:对同一个 $P^{\text{target}}$,存在多个λ(位于不同极化分支上)
138+
2. **内循环可能收敛到错误分支**:λ优化器可能在两个分支之间振荡
139+
3. **跨结构的不连续性**:结构微小变化可能导致 $P^I$ 的分支跳变,使得收敛后的λ发生突变
140+
141+
这种歧义在最坏情况下是**灾难性的**——λ内循环可能收敛到一个λ值,该值在物理上对应的极化与目标极化差了一个极化量子。
142+
143+
**DeltaSpin完全不受此影响**:$M^I$ 是唯一确定的矢量,对其模没有任何模歧义。
144+
145+
---
146+
147+
## 3. 五个困难的纠缠关系
148+
149+
这五个困难不是独立的——它们构成一个**因果链**
150+
151+
```
152+
物理定义全局性(层级一)
153+
154+
├──→ 必须用Berry phase而非实空间积分
155+
│ │
156+
│ └──→ 引入 ∂_k 依赖(层级三)
157+
│ │
158+
│ └──→ 引入规范依赖性(层级二)
159+
│ │
160+
│ └──→ 响应矩阵的1/Δk放大(层级四)
161+
162+
└──→ 引入模2π歧义(层级五)
163+
164+
└──→ 约束框架的多解问题(与层级四耦合)
165+
```
166+
167+
**根本矛盾**:极化是Berry phase → 必须做k空间积分 → 必然涉及∂_k → 必然规范依赖 → 必然有模2π歧义。这不是任何一个实现细节的问题——它是极化物理本质的必然结果。
168+
169+
---
170+
171+
## 4. 什么可以复用,什么必须重新构建
172+
173+
| DeltaSpin组件 | 可复用到DeltaP? | 原因 |
174+
|:---|:---:|------|
175+
| SMO构造 (`build_orb_onsite`) | ✅ 完全复用 | SMO是几何定义,与约束量无关 |
176+
| 二中心积分 (`snap`) | ✅ 完全复用 | ⟨φ\|α⟩重叠积分不依赖约束类型 |
177+
| 双层循环框架 | ✅ 框架复用 | Lagrange乘子法的通用结构 |
178+
| λ CG/Broyden优化 | ✅ 策略复用 | 从M→P,Jacobian需调整 |
179+
| **H^λ构造** | ❌ 本质不同 | M^I的H^λ是零阶\|α⟩⟨α\|,P^I的是一阶\|_kα⟩⟨α\| |
180+
| **HContainer预存储** | ⚠ 部分复用 | H^pre(R)可做∂_k Fourier得到Part A,但Part B是全新的 |
181+
| **规范固定** | ❌ DeltaSpin不需要 | DeltaP的核心新增复杂度 |
182+
| **力/应力修正** | ⚠ 框架复用 | Pulay项的结构相同,但算符内容不同 |
183+
| **子空间对角化加速** | ✅ 部分复用 | 内循环加速策略可迁移 |
184+
185+
### 定量评估
186+
187+
DeltaSpin的核心代码(`dspin_lcao.cpp`, `lambda_loop.cpp`等)中,可直接迁移到DeltaP的部分约占**40%**(SMO构造、二中心积分、CG框架、双层循环逻辑、力修正框架)。需要全新实现的部分占**60%**(k空间Berry connection计算、规范固定、H^λ新结构、∂_k处理)。
188+
189+
但如果采用Wilson loop方案(后处理P^I而非内嵌H^λ),**H^λ构造部分可以完全跳过**——这是Wilson loop方案最吸引人的地方。代价是内循环无法通过修改Hamiltonian来直接驱动电子态,只能通过外循环间接约束。
190+
191+
---
192+
193+
## 5. 结论
194+
195+
### 5.1 DeltaSpin成功的核心原因
196+
197+
DeltaSpin之所以能精确计算自旋相互作用并高效实现,因为**磁矩天然是局域的、规范不变的、无模歧义的实空间量**。SMO投影提供了这个局域量的自然定义,而实空间密度矩阵的原子分区使得H^λ恰好与DFT+U同构——从而几乎零额外成本地复用了已有基础设施。
198+
199+
### 5.2 DeltaP的本质困难
200+
201+
极化是**全局的、规范依赖的、有模2π歧义的k空间量**。复刻DeltaSpin的路径在五个递进层面遭遇障碍:物理定义(全局vs局域)→规范依赖(∂_k引入了相位敏感性)→算符结构(一阶\|_kα⟩⟨α\| vs 零阶\|α⟩⟨α\|)→响应矩阵的放大效应→模2π歧义对约束框架的冲击。这些障碍不是实现细节,而是极化物理本质的直接体现。
202+
203+
### 5.3 最简可行路径
204+
205+
绕过而非对抗这些困难:**Wilson loop方案**跳过了H^λ构造和规范固定的全部复杂度,将P^I的计算留在后处理中。这种设计更接近Diéguez-Vanderbilt (2006)的原始哲学——在SCF外循环中通过调整电场来匹配目标极化,而不是在SCF内部修改Hamiltonian。

0 commit comments

Comments
 (0)