势函数理论计算：从头算势函数开发与机器学习势的拟合策略

势函数理论计算的核心矛盾是精度与效率的权衡。项目组在铜纳米颗粒熔化模拟中面临选择：EAM势函数能跑百万原子、纳秒级模拟，但熔点预测偏差可达数十K（Cu实验熔点1358 K，EAM典型给出1270-1340 K，偏差20-90K视参数化而定）；DFT-MD精度够但只能跑几百原子、几十皮秒。最终用了一条中间路线：先开发针对性的EAM势，再用机器学习势（NEP）做精度修正——在10万原子体系上实现了近DFT精度（熔点偏差<10K）的MD模拟。

Composition of chemistry related design elements to serve as a supporting backdrop for projects on chemistry, scientific knowledge and education

势函数类型谱系

经典力场不是铁板一块，不同类型适合不同体系：

1. 对势（Pair Potential）

• Lennard-Jones、Buckingham、Morse
• 仅描述两体相互作用，不含多体效应
• 适合：惰性气体、简单离子晶体（NaCl）
• 局限：无法描述金属的结合能-配位关系

2. 嵌入原子方法（EAM/MEAM）

• 把每个原子的能量分为嵌入能（电子云密度的函数）和对势能
• EAM：各向同性，适合FCC金属（Cu、Al、Ni）
• MEAM：含角度依赖，适合BCC金属和共价体系
• 优势：正确描述金属的配位-能量关系，计算效率高

3. 经典价键力场（Tersoff、SW、ReaxFF）

• Tersoff/Stillinger-Weber：共价键体系（Si、Ge、C）
• ReaxFF：反应力场，可描述键的断裂和形成
• 优势：化学反应模拟，适合燃烧、催化
• 局限：ReaxFF参数多（>100个），拟合难度大

4. 机器学习势（MLP）

• NEP、GAP、NNP、MTP等
• 用DFT数据训练，精度接近DFT，效率接近经典力场
• 优势：可拟合任意体系的势能面
• 局限：训练数据需求大，外推性差

EAM势函数开发流程

项目组开发Cu-Ar界面势的流程：

Step 1: DFT计算生成训练数据

构建训练数据库，包含多种构型：

• Cu体相：FCC结构，扫描晶格常数（3.4-3.8 Å，步长0.05 Å）
• Ar体相：FCC结构，扫描晶格常数（5.0-5.6 Å）
• Cu-Ar界面：Cu(111)/Ar(111)界面，间距扫描（2.5-4.0 Å）
• Cu表面：不同取向的表面能计算
• Ar在Cu表面的吸附构型

每个构型在VASP中计算能量和力：PBE泛函，ENCUT=500 eV，EDIFF=1E-6，k-mesh根据体系大小调整。总共生成约200个训练构型。

Step 2: EAM势函数参数化

EAM势的总能量：

E_total = Σ F(ρ_i) + 1/2 Σ φ(r_ij)

其中F是嵌入函数，ρ_i是原子i处的电子云密度（邻居贡献之和），φ是对势函数。

项目组用原子嵌入式电子密度（通过DFT计算Bader电荷密度获取）拟合F(ρ)，用Lennard-Jones形式拟合Cu-Ar对势φ(r)。Cu-Cu和Ar-Ar对势直接使用文献中已有的EAM参数（Mishin的Cu势和LJ的Ar势）。

Step 3: 参数优化

用最小二乘法最小化训练数据上的能量和力偏差：

χ² = Σ w_E (E_DFT - E_EAM)² + Σ w_F |F_DFT - F_EAM|²

权重w_E和w_F需要平衡——能量拟合太重会导致力精度差，力拟合太重会导致能量曲面变形。项目组用w_F=10×w_E（力的数据点远多于能量，但每个数据点的权重应相当）。

拟合工具：使用fitEAM软件包，支持多种EAM势函数形式的拟合。优化算法用Levenberg-Marquardt。

Step 4: 交叉验证

200个训练构型中，180个用于拟合，20个留作验证集。验证集上的能量MAE（平均绝对误差）为2.3 meV/atom，力MAE为0.08 eV/Å——对于EAM势来说精度合格。

机器学习势：NEP的尝试

EAM势的根本局限是函数形式固定——无论参数怎么调，都无法描述超出”嵌入能+对势”框架的物理。项目组尝试用神经网络上势（NEP，Neuroevolution Potential）来突破这个限制。

NEP的优势：

• 用径向和角描述符作为输入，神经网络作为拟合函数
• 可以描述多体效应（三体、四体）
• 在GPUMD软件中高度优化，10万原子体系每步<0.1秒

训练数据：项目组用DFT计算了2000个Cu的构型（含体相、表面、缺陷、液态、纳米颗粒），包括能量、力和 virial张量。NEP的训练使用了神经进化算法（不需要梯度），训练时间约6小时（NVIDIA RTX 3090）。

验证结果（典型值范围，具体取决于训练数据和参数化）：

NEP在熔点和空位形成能上明显优于EAM——因为这两个性质涉及配位数的剧烈变化，EAM的固定函数形式无法准确描述。

势函数迁移性评估

一个好的势函数不仅在训练集上表现好，还要在未见过的构型上有合理表现。项目组用以下方法评估迁移性：

1. 温度外推：训练数据只含0K和300K的构型，测试在2000K液态下的表现
2. 应变外推：训练数据含±5%应变，测试在±10%应变下的表现
3. 缺陷外推：训练数据含空位，测试在间隙原子和位错下的表现

EAM在温度外推上表现尚可（液态结构RDF偏差<5%），但在位错核心结构的描述上偏差大（位错宽度偏差30%）。NEP在所有外推测试中偏差<10%——但这是在训练数据覆盖较好的前提下。如果遇到训练数据完全未覆盖的构型（如高压相变），NEP也会给出不可靠的结果。

反思：势函数不是万能的

势函数开发的本质是用有限的函数形式（或有限的训练数据）去近似无限维的势能面。EAM的近似太粗（忽略多体效应），NEP的近似依赖数据覆盖度。没有”最好”的势函数，只有”适合特定问题”的势函数。

选择原则：如果研究体相性质且体系大→EAM；如果涉及化学反应→ReaxFF；如果需要DFT精度但体系太大→MLP；如果体系小→直接AIMD。更多力场选择的讨论可以参考分子动力学栏目，或返回科研学术网首页。