辐照损伤模拟后处理：OVITO缺陷统计与位错环提取实战

模拟跑完那一刻，我始终觉得真正的工作才刚开始。dump文件里躺着几百万个原子坐标，能量曲线平滑得像什么都没发生，可材料内部早就天翻地覆。做辐照损伤、离位级联、位错运动的同行应该都有体会：模拟结束后的第一件事不是画能量曲线，而是把原子轨迹里的缺陷一个个找出来，看看晶体究竟伤成了什么样。

这篇文章围绕这个场景展开，核心是用OVITO从分子动力学轨迹里找出弗伦克尔缺陷，再用位错提取算法分析位错环。VESTA在前期建模和参考构型准备里的作用我也会一并说清楚——没有一份干净的参考结构，后面所有缺陷统计都是空中楼阁。文章里写到的参数、流程和判断方法，都来自我实际跑辐照损伤项目时的经验，尤其适合刚接触LAMMPS+OVITO这条链路、准备做点缺陷统计的同行参考。

1. 为什么辐照模拟结束，首先要交一张弗伦克尔缺陷清单

1.1 弗伦克尔缺陷为什么是辐照损伤的“基础货币”

在晶体材料里，一个原子离开原有格点位置、跑到间隙位置，原地就留下一个空位。空位加间隙原子，合称一个弗伦克尔缺陷（Frenkel pair）。别小看这样一个点缺陷，它几乎是所有辐照损伤演化的起点。

在分子动力学模拟中，PKA（初级碰撞原子）获得一定能量后，会在晶格中引发一连串碰撞级联。级联核心附近会同时产生大量空位和间隙原子，有的立刻复合，有的会扩散开，有的则稳定存在下来。统计一级碰撞后残余的弗伦克尔缺陷数量，是判断材料抗辐照性能最直接的量化指标。如果你连这个数都统计不准，后面关于缺陷演化、位错环尺寸分布、辐照硬化机理的讨论，全部缺乏根基。

实际项目里，我见过有人直接用“缺陷原子数目”来当结论，但模型里并不只有孤立的空位和间隙。一个级联区往往同时存在几十个空位和数量接近的间隙原子，它们位置随机、分布不均。只有通过Wigner-Seitz这类严格对比格点的方法，才能从海量坐标里把真正的缺陷数量捞出来。

1.2 位错环的两种出身：空位凝聚和间隙凝聚

点缺陷在足够高的温度和足够长的时间下会发生迁移，辐照损伤模拟里常见的位错环，主要对应两种物理过程，分析时必须分开对待。

一种来自间隙原子聚集，形成间隙型位错环。这类环在BCC金属里常以伯格斯矢量1/2<111>出现，也可能形成<100>矢量环，和滑移面、层错能都有关系。另一种来自空位聚集，级联中心往往空位浓度极高，退火演化时容易塌缩成空位型位错环，这类环的伯格斯矢量、稳定性和间隙型环有明显差异。

因为两类环的来源和后续演化规律不一样，做完DXA后只报告“体系里出现了位错环”是远远不够的。你需要把环的伯格斯矢量、环面法向、环内缺陷成分三个信息结合起来，才能判断环到底是空位型还是间隙型。这也是我在后面专门花篇幅讲伯格斯矢量读取的原因。

1.3 为什么不能靠肉眼在三维视图里找缺陷

刚接触分子动力学后处理时，很多人会想：把原子渲染出来，缺陷不就是“空了一块”吗？实际情况完全不是这样。密排晶体里，即便一个原子离开了格点，周围原子也会在热运动下立刻部分回填或畸变，从宏观渲染图上看不出某个位置“空了”。

更麻烦的是热振动。在几百K的模拟温度下，原子偏离平衡位置约0.1到0.2埃，而间隙原子产生的特征位移也就在亚埃到埃量级。靠鼠标旋转视角去肉眼识别缺陷，基本等于大海捞针，而且谁也没法复现你的“目测标准”。

所以专业工具选择用确定性几何判据，把“某个格点有没有原子”变成纯粹的计算问题。OVITO里的Wigner-Seitz缺陷分析就是干这件事的。

2. 从VESTA到OVITO：一条完整的缺陷分析工具链

2.1 VESTA在分子动力学流程里的真实价值

标题把VESTA放在最前面是有原因的。很多人把它只看成“画晶体结构的小工具”，但在分子动力学前后处理里，VESTA承担了三个关键任务：构建初始晶体模型、检查结构合理性、导出标准结构文件。

以最常见的BCC铁辐照模拟为例，我会先在VESTA里用空间群Im-3m生成一个单胞，再通过扩展功能设置2×2×2甚至10×10×10的超胞。超胞尺寸直接决定后面MD盒子的大小，而盒子大小又关系到级联缺陷数量统计是否可靠——盒子太小会出现级联与自身周期性镜像重叠，缺陷数被严重低估；盒子太大则计算成本失控。前期建模时多花几分钟，后期分析能少走很多弯路。

VESTA本身支持CIF、POSCAR、XYZ、LAMMPS data等格式导出，是衔接建模软件和MD模拟器最常见的桥梁。安装也简单，官方提供Windows、macOS和Linux版本，学术用途免费，解压即可用，不需要额外配置环境。

2.2 用VESTA导出参考结构时最容易踩的坑

VESTA导出的参考结构必须和MD模拟的初始结构完全一致，这句话我在项目里反复强调。具体涉及四个方面：原子种类和数量一致、晶格常数一致、盒子方向和晶向一致、坐标单位一致。VESTA默认使用埃，LAMMPS在metal单位下同样用埃，一般对得上，但换到别的软件或单位体系时一定要检查。

最隐蔽的坑是晶向排列。如果你在VESTA里手动旋转过模型再导出，坐标方向和MD设定就会错位。Wigner-Seitz缺陷分析依赖当前构型和参考构型之间的逐格点对应关系，参考构型的晶向一旦转了方向，每个格点的归属都会错乱，统计出的空位和间隙原子会散布在整个盒子里，完全没有参考价值。

我的习惯是：在VESTA里把超胞的a、b、c三轴方向显示出来，记录对应晶向，再在LAMMPS输入文件里严格按相同方向写lattice命令。这样OVITO做分析时，参考构型和模拟轨迹的盒子可以无缝对接。

2.3 OVITO能做什么，先用什么后用什么

OVITO是分子动力学后处理里最常用的可视化与分析工具，支持LAMMPS、VASP、Gromacs等格式，也能处理超大轨迹。基础版承担结构显示、粒子选择、渲染图生成；Wigner-Seitz缺陷分析、CNA（共邻分析）、DXA（位错提取）这类高级分析一般在OVITO Pro里完整支持，也可以通过Python接口调用对应算法模块。

实际分析时，我建议把能自动化的步骤全部挂进OVITO的修改器流水线。OVITO的核心逻辑是：加载轨迹后，依次添加多个修改器，每个修改器输出一项中间结果，后面的修改器消费前面的结果。比如先加Wigner-Seitz缺陷分析，再加CNA，最后挂DXA，全套下来就是一条完整的缺陷识别链路。这种流水线设计对参数调整和结果回放都非常友好，比在命令行里反复改参数重算省心得多。

3. Wigner-Seitz缺陷分析：把空位和间隙原子一个个捞出来

3.1 Wigner-Seitz胞到底在数什么

Wigner-Seitz缺陷分析的核心，是把“晶体里有没有缺陷”这个模糊问题，转化成“每个格点对应的空间区域里住了几个原子”这个精确的计数问题。

可以把完美晶格想象成一片由无数工位组成的办公室，每个工位对应一个原子的原始位置。Wigner-Seitz胞就是这些工位的严格几何定义：以某个格点为中心，取所有相邻格点连线的垂直平分面围成的区域。OVITO把这个定义推广到三维周期边界下的任意晶格。

分析时，OVITO先在参考构型上建立Wigner-Seitz胞的划分，然后把当前构型里的每个原子分配进对应胞内数数：

• 某个胞里没有原子：这个格点丢了原子，记为1个空位；
• 某个胞里有2个或多个原子：多出来的原子就是间隙原子；
• 某个胞里正好1个原子：这个格点被正常占据。

这正是Wigner-Seitz方法的最大优势：它直接比较的是“原子有没有待在正确的格点区域”，而不是看原子离某个位置多远。热振动造成的原子位移不会轻易被误判成缺陷——位移再大只要还在自己的胞内，就不算跑位；只有真正跑到别的格点区域里，才会被记为间隙原子。

3.2 实际操作：加载轨迹、绑定参考构型、运行分析

在OVITO里的操作流程大致如下，我按BCC铁的LAMMPS轨迹举例。

第一步，启动OVITO，载入LAMMPS输出的dump文件。多帧轨迹会被当作时间序列载入。如果dump文件过大，加载时间会明显变长，建议先用部分帧调试参数，确认无误再载入完整轨迹。

第二步，在右侧修改器列表里添加“Wigner-Seitz缺陷分析”。参数面板里最关键的设置是参考构型。推荐选择“使用文件”，载入提前准备好的完美晶格构型文件；如果轨迹第一帧本身就是完美晶格，也可以勾选“使用输入构型的第一帧”，前提是你能确认第一帧没有原子删除或重叠。

第三步，设置“只考虑子集”。比如体系里除了铁的基体原子，还有用于标记PKA的特殊原子或者杂质，需要把它们排除掉，否则PKA原子会被当成额外间隙原子统计进去。通常只保留晶格原子参与缺陷识别。

第四步，运行修改器，在输出表里看空位和间隙原子数量，再用颜色编码把缺陷原子高亮显示。

这里有个非常实用的技巧：在“颜色编码”里选“粒子类型”，把空位和间隙原子分别设成不同颜色，基体原子调成半透明。这样从渲染图上一眼就能看清缺陷分布在哪个区域，判断级联损伤范围非常方便。

3.3 结果怎么读：空位、间隙原子和替换原子

Wigner-Seitz分析跑完后，OVITO会给每个粒子新增属性，包括粒子类型、所在格点胞编号、位移量等。通过粒子类型筛选，可以轻松把空位和间隙原子变成两个独立选择集。

以BCC铁、10 keV PKA级联模拟为例，典型结果是：残余空位几十个、间隙原子几十个，两者数量大体接近。物理原因很简单——弗伦克尔缺陷成对产生，少数差异来自周期性边界和统计口径。如果空位数和间隙原子数差得离谱，第一反应应该是检查参考构型有没有问题，再看边界条件设置是否合理。

这里要特别解释一下替换原子。在纯元素体系里，W-S计数不会出现“替换原子”概念；但在合金或含杂质体系里，一个格点上的原子可能从A换成了B，胞内仍然只有1个原子，W-S结果不会显示任何缺陷。这种化学上的置换需要结合化学短程有序参数或邻居类型统计来分析。格点被占据，不等于没有发生化学变化，这是两个层面的问题。

3.4 缺陷数量随时间演化的曲线怎么出

单帧统计只是起点，辐照损伤模拟里更有价值的输出，是缺陷数量随时间的演化曲线：级联产生阶段缺陷数快速上升，之后由于空位-间隙原子复合慢慢下降，最终趋于平稳或形成平台。

想要一条连续的演化曲线，最省力的做法是载入完整轨迹，添加W-S修改器，在输出设置里启用“为每帧记录结果”，把每一帧的空位数和间隙原子数导出为CSV，再交给绘图工具画图。

如果轨迹有几千帧，先用“时间轴”设置里的隔帧载入试跑一遍，看趋势是否正常，确认无误后再载入全部帧做精确统计。这一步能省下大量等待时间，dump文件超过几十GB时效果尤其明显。

4. 位错环提取：DXA参数调对了，环才是环

4.1 为什么点缺陷分析还不够，还要看位错环

弗伦克尔缺陷是点缺陷，位错环是线缺陷的闭合形态。辐照模拟中，点缺陷聚集到一定程度会以位错环的形式体现。评估辐照后材料硬化时，位错环的类型、尺寸和密度，往往比点缺陷总数更有解释力。

所以确认缺陷数量后，下一步通常是问：这些缺陷有没有形成环？环是间隙型还是空位型？尺寸多大？OVITO里处理这个问题的标准工具是DXA，全称Dislocation Extraction Algorithm。它的原理可以粗略理解为：在晶体内部划分三角形网格，通过比较网格结点的邻居关系，找出晶格不匹配的路径，从而追踪位错线走向，并计算出每段位错线的伯格斯矢量。

4.2 DXA和CNA的分工：先判断结构，再抽取位错

实际运行时，DXA强烈依赖它前面的CNA结果来区分原子属于哪种晶体结构。很多人遇到“DXA提取完什么都没显示”，问题根源常常不在DXA，而在CNA参数。

CNA通过分析每个原子最近邻居之间的连接关系，判断局部结构是FCC、BCC、HCP还是“其他”。对BCC铁来说，绝大多数原子应显示为BCC，缺陷和位错核心区域才会显示为“其他”或错误堆垛类型。

在OVITO里添加DXA修改器后，DXA会在CNA分类基础上搜索位错线。如果CNA把BCC误判成FCC，整个晶格背景就会七零八落，DXA自然找不到连续的位错线。我的习惯是：先单独跑CNA，看渲染颜色是否大面积符合基体结构，确认无误后再挂DXA。这一步检查只花一分钟，却能避免后面所有无效操作。

4.3 DXA核心参数和我的常用设置

DXA参数面板里几项关键参数，我来逐一说清它们的逻辑：

网格尺寸是最需要微调的一项。它的物理含义是搜索位错线时空间离散化的粗细。网格太小，位错线会被切得支离破碎，环的闭合性很差；网格太大，微小位错环可能被吞掉，只留下一个光秃秃的节点。

那怎么判断参数是否合适？直接看输出的位错线是否连续、是否闭合、伯格斯矢量是否稳定。最专业的做法是做一次小范围参数扫描，把网格尺寸按0.05 nm步长从0.2 nm试到0.6 nm，记录环数和总位错长度的变化，然后选取结果稳定区间的参数。

4.4 伯格斯矢量与位错环类型的判断方向

位错线提取完成后，OVITO会把它作为独立对象输出。选中任意一段位错线，右侧信息面板会显示方向和伯格斯矢量。伯格斯矢量是位错最核心的特征量。

常见对应关系有必要列一下：

• BCC结构：1/2<111>和<100>伯格斯矢量都比较常见；
• FCC结构：全位错通常是1/2<110>，不全位错可能是1/6<112>；
• HCP结构：基面位错对应最常见的伯格斯矢量类型。

判断位错环是空位型还是间隙型，光看伯格斯矢量并不够，因为同一伯格斯矢量可能对应两种环。更可靠的方法是把W-S分析和DXA结果叠加起来看：如果环附近存在明显的间隙原子富集，大概率是间隙型环；如果环附近以空位富集为主，则偏向空位型。有些项目还会通过退火演化时环的尺寸变化方向来辅助判断，两种环在吸收点缺陷时的长大规律不一样。

4.5 环长统计与时间演化分析

OVITO里提供了“提取位错环”功能，可以自动识别闭合位错线，并计算周长、长度等参数。把每条环的信息导出，就能得到一组环尺寸分布数据，这对辐照损伤模型验证非常关键。

做时间序列分析时，周期性边界条件是个容易忽略的坑。如果位错环跨越了盒子边界，从一侧出去又从另一侧进来，DXA在周期性图像下通常能识别出来，但直接导出坐标时可能出现截断。如果发现某条环长度突然异常大，先检查它是否跨过了盒子边界，打开周期性显示确认环是完整连续的，而不是被边界切成了两段。

5. 实战中容易翻车的几个环节与我的排查路径

5.1 参考构型漂移和盒子尺寸漂移

使用NPT系综或者其他允许盒子变化的系综时，盒子尺寸会随温度和应力改变。这会给W-S分析带来直接麻烦：参考构型晶格常数和当前帧不一致，W-S胞的划分整体失真，统计结果里全是假缺陷。

解决办法有两个方向。一是尽量用模拟早期、体系尺寸还没有明显热膨胀的构型做参考；二是如果盒子尺寸确实变了，在OVITO里加一个“仿射变换”修改器，把当前帧盒子对齐到参考构型，再做W-S分析。更稳妥的方案是，在分析目标的时间段刻意选用NVE或NVT系综，避免盒子尺寸变化带来的额外干扰。

如果发现W-S统计出的缺陷数量在某一帧后突然整体抬高，先别急着上升到物理机制，优先排查是否盒子尺寸或参考构型发生突变。

5.2 自由表面和晶界带来的假缺陷

表面和晶界的原子配位数与完美晶格不同。在W-S分析里，表面悬挂原子可能被计成间隙原子，表面下方的空缺格点会被计成空位。做纳米线或带空洞模型时，这些表面假缺陷完全可能淹没真实缺陷信号。

我的处理习惯有三个：只关心体内缺陷时，用OVITO的裁剪工具把表面层切掉；关心界面效应时，明确分开统计界面缺陷和体内缺陷；在论文方法部分注明统计区域和统计口径，保证别人能复现结果。缺陷总数本身没有意义，有统计边界的缺陷总数才有意义。

5.3 超大数据量下的性能优化

一次几百万到上千万原子的级联模拟，dump文件可能达到几十GB。直接在OVITO里全量加载再逐帧跑DXA，计算资源消耗非常大。几个亲测有效的优化方式：

• 降低轨迹输出频率。在MD输入文件里把dump间隔从0.1 ps调整为0.5 ps甚至1 ps，只要不丢失关键的级联演化信息即可；
• 使用OVITO的载入时裁剪功能，先用包围盒选中缺陷可能出现的区域，再加载；
• 批量分析时使用OVITO Python脚本的无界面模式，把结果输出成CSV或筛选后的目标文件。

DXA本身计算量较大，建议先用较粗网格快速跑一遍，确认缺陷大概在哪些位置、环数数量级是否合理，再用精细网格出最终结果。粗网格会漏细节，但能暴露出明显的参数错误，节省大量时间。

5.4 多修改器互相验证，避免单一结果误导

我在实际项目里很少只依赖一个修改器的结果。最常用的一套组合是W-S缺陷分析、CNA、DXA三者相互对照：

• W-S告诉你哪些格点少了原子、哪些格点多了原子；
• CNA告诉你哪些区域原子的局部堆垛次序已经改变；
• DXA告诉你位错线的走向和伯格斯矢量。

这三者强相关但不等价。如果在某个区域看到大量间隙原子，但CNA显示该区域仍是完美BCC结构，那很可能是参考配置错误或者表面假象。反过来，如果CNA显示一片区域已经非晶化，位错线在这里根本无法定义，而W-S却统计出大量缺陷，那就要判断这片区域到底出了什么物理过程。

把三个修改器叠加渲染，是我在辐照模拟论文出图前最后一道自检工序。有一次我差点把位错环分析中一个跨边界环的投影图画出去，就是因为同时打开了周期性显示才发现环其实跨过了盒子两侧，才避免了返工。

最后分享一个工作习惯：每次出图时，我会把OVITO的状态文件（.ovito）一起保存。分析参数、参考构型路径、颜色设置全部记录在状态文件里，过几个月回来看图，还能清楚知道当时的网格尺寸、平滑等级和具体分析流程。别人问起这张图怎么生成的，也能逐条复现，而不是靠记忆猜测。课题周期长、项目跨度大的时候，这套参数留痕的习惯能省掉大量解释成本。