力场与势函数：LAMMPS模拟的灵魂

做分子动力学模拟，说白了就是在算力。算什么呢？算原子之间的相互作用力。而这个「怎么算」的规则，就是力场，也叫势函数。

我刚开始学LAMMPS那会儿，总觉得势函数就是个黑盒子。选一个，设几个参数，跑就完了。直到有一次，我模拟一个金属纳米颗粒的熔化过程，用了LJ势，结果熔点差了300多K。嗯，从那以后，我再也不敢随便选势函数了。

这一章，咱们就把势函数这件事彻底讲透。

常见势函数类型：LJ、EAM、ReaxFF

LAMMPS里支持的势函数有上百种，但常用的其实就那么几类。我按应用场景给你捋一捋。

LJ势：最经典的入门级选手

Lennard-Jones势，简称LJ势。它描述的是两个中性原子之间的范德华相互作用。公式长这样：

E = 4ε [(σ/r)^12 - (σ/r)^6]

其中ε是势阱深度，σ是原子直径。r^12项描述排斥，r^6项描述吸引。

适用场景：惰性气体（氩、氪）、简单分子液体、粗粒化模型。

我的经验：LJ势虽然简单，但千万别小看它。我做过一个聚合物熔体的粗粒化模拟，用LJ势描述单体间的非键相互作用，结果和全原子模拟吻合得很好。关键是参数调对了。

小技巧：LJ势的截断半径一般取2.5σ。太小了能量不连续，太大了计算量暴增。这是经验值，别问我为什么，试出来的。

EAM势：金属模拟的标配

Embedded Atom Method，嵌入原子法。它和LJ最大的区别是：EAM考虑了原子所处的「环境」。一个原子在金属块体内部和表面，能量是不一样的。

EAM势的能量由两部分组成：

  • 对势项：两两原子间的排斥
  • 嵌入能项：原子嵌入到电子云中的能量

适用场景：面心立方（FCC）金属（Cu、Al、Ni、Au等）、体心立方（BCC）金属（Fe、W等）、合金体系。

避坑指南：我曾经用EAM模拟Cu-Ni合金，直接从网上找了个Cu的势函数和Ni的势函数拼在一起用。结果算出来的混合焓完全不对。后来才知道，合金体系必须用专门拟合的交叉势函数，不能自己拼凑。

ReaxFF：化学反应模拟的利器

ReaxFF是反应力场。它最大的特点是：键可以断裂，也可以生成。这是传统力场做不到的。

ReaxFF的能量项非常复杂，包括键能、角能、二面角能、共轭能、非键相互作用等等。它通过键级（bond order）来判断原子间是否成键。

适用场景：燃烧、爆炸、催化、材料表面化学反应。

注意：ReaxFF的计算量非常大。我算过一个含10000个原子的燃烧体系，用ReaxFF跑了1纳秒，在128核的集群上跑了整整3天。所以，能用经典力场解决的问题，别轻易上ReaxFF。

pair_style与pair_coeff设置

在LAMMPS里，势函数通过两个命令来设置：pair_style 和 pair_coeff。

pair_style 定义势函数的类型和全局参数，pair_coeff 定义具体的原子类型参数。

基本语法

# LJ势
pair_style lj/cut 10.0
pair_coeff 1 1 0.0104 3.40

# EAM势
pair_style eam
pair_coeff * * Cu_u3.eam

# ReaxFF
pair_style reax/c NULL
pair_coeff * * ffield.reax C H O

参数说明：

  • lj/cut 10.0：LJ势，截断半径10.0 Å
  • 0.0104 3.40：ε=0.0104 eV，σ=3.40 Å
  • Cu_u3.eam：EAM势函数文件
  • ffield.reax：ReaxFF力场参数文件
重要：pair_coeff中的原子类型编号，必须和data文件中Atoms部分的原子类型一一对应。我见过太多人在这里栽跟头了。

势函数文件格式与获取途径

文件格式

不同的势函数，文件格式差别很大。

LJ势：参数直接写在in文件里，不需要单独的文件。

EAM势：文件格式有标准规定。第一行是注释，后面是元素符号、原子量、晶格常数。然后是电子密度函数和嵌入能函数的表格数据。

ReaxFF：参数文件包含几十种原子类型的键、角、二面角、非键等参数。格式非常复杂，一般不需要手动修改。

获取途径

势函数文件从哪里来？我一般用这几个渠道：

  1. LAMMPS官网：https://lammps.sandia.gov/potentials/ 有大量现成的势函数文件
  2. NIST Interatomic Potentials Repository：美国国家标准与技术研究院维护的势函数库
  3. OpenKIM：一个开源的势函数知识库，支持自动验证
  4. 文献：很多论文会提供势函数文件作为补充材料
我的习惯：下载势函数文件后，先用LAMMPS自带的测试例子跑一遍，确认能正常计算再用于自己的体系。这一步能省去后面大量的调试时间。

混合力场设置技巧

实际模拟中，经常遇到不同区域需要不同势函数的情况。比如：金属表面吸附有机分子，金属用EAM，有机分子用力场。

LAMMPS支持混合力场，通过 pair_style hybrid 来实现。

基本设置

# 定义混合力场
pair_style hybrid lj/cut 10.0 eam

# 金属-金属相互作用用EAM
pair_coeff 1 1 eam Cu_u3.eam

# 有机-有机相互作用用LJ
pair_coeff 2 2 lj/cut 0.0104 3.40

# 金属-有机交叉相互作用用LJ
pair_coeff 1 2 lj/cut 0.0052 3.40

关键点：

  • 交叉相互作用需要单独指定
  • 不同势函数之间的参数要匹配
  • 混合规则（如Lorentz-Berthelot规则）要明确
我曾经踩过的坑：混合力场中，不同势函数的截断半径不一致，导致能量不守恒。后来我统一把所有截断半径设为相同值，问题就解决了。

知识体系结构图

下面这张图，帮你理清本章的知识脉络：

力场与势函数 常见势函数类型 LJ势 EAM势 ReaxFF pair_style与pair_coeff 语法规则 参数匹配 原子类型对应 文件格式与混合力场 文件格式 获取途径 混合力场设置 核心：势函数决定模拟的准确性

这张图把本章内容分成了三大块：势函数类型、命令设置、文件与混合力场。你顺着这个结构往下学，思路会清晰很多。

好了，关于力场和势函数，核心内容就这些。记住一句话：选对势函数，模拟就成功了一半。另一半是什么？是参数的正确设置和充分的测试。别偷懒，每个新体系都要先做验证。

最后提醒：势函数文件不要随便从网上下载就拿来用。一定要确认文件格式、单位制、元素类型是否和你的体系匹配。我见过有人用Cu的EAM势去模拟Al，结果算出来的晶格常数差了10%。这种低级错误，完全可以避免。

公众号:蓝海资料掘金营，微信deep3321

← 上一章 📖 返回目录 下一章 →