AtomPPGen 教程 07：完整端到端示例（以 Al 为例）

本 notebook 以 Al（Z=13） 为例，把赝势生成流程从头到尾串联成一个“可一键 Run All”的最小闭环：

  1. 全电子原子解（AE）

  2. TM 伪化（构造平滑伪轨道）

  3. 半局域势反演（得到 \(V_l(r)\)）

  4. KB 可分离形式（得到 \(V_{\mathrm{loc}}\), \(\beta_l\), \(D_l\)）

  5. 可转移性验证（范数守恒 / 对数导数 / 幽灵态）

  6. 导出（JSON+NPZ）

如果你已经理解 00-05 的每个环节，这个 notebook 可以作为“把知识变成产物”的工作模板。

[1]:

from pathlib import Path
import platform
import numpy as np
import matplotlib.pyplot as plt

from atomppgen import (
    solve_ae_atom,
    tm_pseudize,
    invert_semilocal_potential,
    kb_transform,
    export_pseudopotential,
)
from atomppgen.validate import run_full_validation

# 中文字体配置（兼容多平台）
if platform.system() == 'Darwin':
    plt.rcParams['font.sans-serif'] = ['Arial Unicode MS', 'Heiti TC', 'STHeiti']
elif platform.system() == 'Windows':
    plt.rcParams['font.sans-serif'] = ['Microsoft YaHei', 'SimHei']
else:
    plt.rcParams['font.sans-serif'] = ['DejaVu Sans', 'WenQuanYi Micro Hei']
plt.rcParams['axes.unicode_minus'] = False

print('Ready')

Ready

Step 1：全电子原子解（AE）

AE 解提供参考能级与波函数，它决定了后续伪化的目标：

  • 外区（\(r>r_c\)）必须与 AE 波函数一致

  • 内区（\(r\le r_c\)）可以重写，但要满足模守恒与光滑性

[2]:

Z = 13
xc = 'PZ81'
lmax = 2

ae = solve_ae_atom(
    Z=Z,
    spin_mode='LDA',
    xc=xc,
    lmax=lmax,
    grid_type='exp_transformed',
    grid_params={'n': 800, 'rmax': 100.0},
    scf_params={'tol': 1e-6, 'maxiter': 150},
)

print('AE eigenvalues (Ha):')
for l in range(lmax + 1):
    print('l=', l, 'eps=', ae.eps_by_l[l][-1])

AE eigenvalues (Ha):
l= 0 eps= -0.2469907526930879
l= 1 eps= 0.0009146350242940937
l= 2 eps= 0.00718614024367731

Step 2：TM 伪化（构造平滑伪轨道）

TM 伪化的关键是“在内区抹平节点/尖峰，同时保持外区完全一致”，并用模守恒锁定散射性质。

[3]:

rc_by_l = {0: 2.1, 1: 2.2, 2: 2.4}
tm_dict = {}
for l in range(lmax + 1):
    tm_dict[l] = tm_pseudize(
        ae.r, ae.w,
        ae.u_by_l[l][-1],
        ae.eps_by_l[l][-1],
        l=l,
        rc=rc_by_l[l],
    )

for l, tm in tm_dict.items():
    print('l=', l, 'rc=', tm.rc, 'norm_error=', tm.norm_error)

l= 0 rc= 2.1 norm_error= 4.721090412432319e-16
l= 1 rc= 2.2 norm_error= 5.935084483710011e-16
l= 2 rc= 2.4 norm_error= 2.9757430941508053e-16

Step 3：势反演（得到半局域势 \(V_l(r)\)）

直觉上：当我们“换了波函数”，就必须同步“换一个势”，使它仍是同一能量下的解； 反演公式把这个要求变成逐点代数计算。

[4]:

inv_dict = {l: invert_semilocal_potential(tm_dict[l], ae.r) for l in tm_dict}
print('Inversion done:', list(inv_dict.keys()))

Inversion done: [0, 1, 2]

Step 4：KB 可分离形式（得到 \(V_{\mathrm{loc}}\), \(\beta_l\), \(D_l\)）

半局域势在平面波里昂贵，KB 通过秩-1 投影把非局域部分写成少量标量积，从而显著加速。

[5]:

kb = kb_transform(inv_dict, {l: ae.u_by_l[l][-1] for l in tm_dict}, ae.r, ae.w, loc_channel=2)
print('KB loc_channel =', kb.loc_channel)
print('KB nonlocal channels =', sorted(kb.beta_l.keys()))
print('KB D_l =', kb.D_l)

KB loc_channel = 2
KB nonlocal channels = [0, 1]
KB D_l = {0: np.float64(4.425529213473979), 1: np.float64(2.6664472703540274)}

Step 5：可转移性验证

验证不是“锦上添花”，而是赝势可用性的底线：

  • 范数守恒：锁定核内电荷与散射信息

  • 对数导数：在能量窗口内比较散射性质（相移）

  • 幽灵态：排除非物理的深束缚态

[6]:

report = run_full_validation(
    ae,
    tm_dict=tm_dict,
    inv_dict=inv_dict,
    r_test=3.0,
    E_range_Ry=(-0.5, 0.5),
    E_step_Ry=0.05,
)

print('overall_passed =', report.overall_passed)
print('norm pass =', all(r.passed for r in report.norm_results.values()))
print('log-derivative pass =', all(r.passed for r in report.log_deriv_results.values()))
print('ghost pass =', report.ghost_result.passed)

overall_passed = True
norm pass = True
log-derivative pass = True
ghost pass = True

Step 6：导出（JSON+NPZ）

导出后你将得到两类文件：

  • *.json：面向阅读/审计（参数、单位、验证报告、以及 NPZ 文件索引）

  • *.npz：面向数值计算（网格、半局域势、KB 投影子等数组）

[7]:

cwd = Path.cwd().resolve()
repo_root = next(p for p in [cwd, *cwd.parents] if (p / 'AtomPPGen' / 'pyproject.toml').is_file())
out_dir = repo_root / 'AtomPPGen' / 'outputs' / 'tutorial_07'
out_dir.mkdir(parents=True, exist_ok=True)
output_prefix = out_dir / 'al_lda_complete'

files = export_pseudopotential(
    ae_result=ae,
    tm_dict=tm_dict,
    inv_dict=inv_dict,
    validation_report=report,
    kb_result=kb,
    output_prefix=str(output_prefix),
    formats=['json', 'npz'],
)
files

[7]:

[PosixPath('/Users/gilbert/Downloads/同步空间/Work/Courses/大四上/电子结构理论与计算/代码作业/2-3/AtomPPGen/outputs/tutorial_07/al_lda_complete.json'),
 PosixPath('/Users/gilbert/Downloads/同步空间/Work/Courses/大四上/电子结构理论与计算/代码作业/2-3/AtomPPGen/outputs/tutorial_07/al_lda_complete.npz')]