用 OVITO-Python 做 LAMMPS 后处理：统计 Type=1 原子沿 x 方向的密度分布（附可扩展到速度/势能/应力）

引言： 在分子动力学（MD）模拟的后处理中， 沿某个方向的空间分布（profile） 是非常常见的分析需求，例如数密度 / 质量密度、速度分布、应力（压力）分布、电荷分布、势能分布等。理想情况下，我们可以在 LAMMPS 的 in 文件里用 compute chunk/atom + fix ave/chunk 在运行过程中直接对空间分箱（binning）并做时间平均；而且 fix ave/chunk 还能对多种 per-atom 量 （如速度、力、电荷、势能、应力等）做按 chunk 的汇总 / 平均输出。 但很多时候模拟已经跑完，脚本里没写这些 profile 计算，手头只剩下轨迹文件（dump / xyz / lammpstrj）。这时再回去重跑既费时又不必要——更合适的做法是：使用 后处理工具 从轨迹中重新统计分布。 OVITO 的 Python 模块提供了一套非常清晰的“数据流水线（pipeline）”机制： 导入数据（ import file ）→ 依次应用 modifiers → compute() 得到结果 → export file 导出 ，非常适合把后处理流程脚本化、自动化。 4.0.png --- 一、为什么要做“后处理密度分布” 在 LAMMPS 里， 密度 / 电荷 / 速度 / 应力 等空间分布，本来可以在 in 文件里通过 chunk + 平均 直接算出来（典型做法是 compute chunk/atom + fix ave/chunk ），并输出为 profile 文件。LAMMPS 官方文档也明确：chunks 可以按空间 bin 切分，并对每个 chunk 统计包括密度、速度、应力、势能等多种量。 但现实经常是： 计算已经跑完了； 当时没有在 in 文件里写 profile 的计算； 现在只剩 dump / xyz 轨迹文件。 这时就需要 后处理 。 OVITO 的 Python 模块非常适合做这一类“按空间分箱统计”的任务： 读取轨迹 → 选择原子子集 → 做直方图 / 分箱 → 导出表格。 --- 二、核心原理：直方图（Histogram）= 沿 x 方向分箱计数 我们把盒子在 x 方向切成 bin count 个小区间（bin），对每个 bin： 计数 ： N bin （这个 bin 里有多少个目标原子） 若要“数密度”（number density）而不仅是计数： 4.1.png 其中 Ly 、 Lz 是盒子在 y、z 方向长度（假设截面积不随时间变化，或你逐帧计算）。 OVITO 的 HistogramModifier 支持多种归一化方式： AbsoluteCount （每个 bin 的样本数） RelativeFrequency （频率 = bin count / 总样本数 ） ProbabilityDensity （概率密度 = bin count / 总样本数 / bin width ） 这些模式在 OVITO 3.13.0 起提供。 4.2.png 注意： ProbabilityDensity 是“概率密度”，不是物理意义上的“数密度（ /体积）”。要得到真正数密度，仍建议用 V bin 做体积归一化。 --- 三、环境准备：安装 OVITO Python 模块 两种常见方式： 方式 A：pip（普通 Python 环境） PyPI 页面列出了支持的平台与 Python 版本范围（例如 Windows / Linux x86 64、Python 3.10–3.13 等）。 方式 B：conda（推荐给 conda 用户） OVITO 官方给了 conda 安装命令（含 ovito 模块；并提示 conda 包也可能带桌面端组件）。 --- 四、数据准备：你能统计什么，取决于 dump 里有什么 4.1 统计密度（位置分布）最低要求 dump / xyz 至少要有： Particle Type （或 type ） Position.X/Y/Z 如果你用的是 legacy XYZ （列含义不固定），需要用 import file(..., columns=[...]) 显式告诉 OVITO 每一列是什么含义，例如： python pipeline = import file( "dump.xyz", columns=["Particle Type", "Position.X", "Position.Y", "Position.Z"] ) OVITO 的 import file() 会自动识别格式并创建 pipeline。 4.2 想统计速度 / 势能 / 应力：必须提前写进 dump 速度：dump 里要有 vx vy vz 势能：需要 LAMMPS compute pe/atom 并用 dump custom 输出该列 应力：需要 compute stress/atom （或相关命令）并输出到 dump；LAMMPS 文档也给过“用 per-atom stress 做 1D 压力 profile”的示例思路 关键提醒： 后处理不能凭空“算出”你没输出的 per-atom 势能 / 应力。 OVITO 能做的是：对“文件里已有的每原子属性”进行统计、分箱、平均。 --- 五、最小可用：Type=1 原子沿 x 的直方图 你原始代码非常接近正确，但建议补上两点： 显式开启 fix xrange=True ，否则 xrange start/end 可能不生效 可指定 normalization （默认是绝对计数） python from ovito.modifiers import HistogramModifier, SelectTypeModifier from ovito.io import import file, export file pipeline = import file("dump.xyz") 或 dump.lammpstrj 等 选出 Type=1 pipeline.modifiers.append( SelectTypeModifier(operate on="particles", property="Particle Type", types={1}) ) 沿 x 分箱做直方图 hist = HistogramModifier(