1. LAMMPS分子动力学建模工程导论

分子动力学(Molecular Dynamics, MD)是一种以经典牛顿运动方程为核心、在原子尺度上追踪粒子随时间演化的计算方法。它的基本假设是:原子核的运动可用经典力学描述,电子效应被“折叠”进经验势函数中。给定初始位置和速度后,每个原子受到的力由势函数对坐标的负梯度给出,即 F_i = -∂U/∂r_i,再通过数值积分求解 m_i·a_i = F_i,逐步推进时间。因此,MD 的输出本质上是“轨迹”——每一时刻所有原子的坐标、速度与受力,而宏观量(温度、压力、能量、应力、扩散系数等)都从轨迹统计得到。

理解 MD 需要抓住三个关键量:势函数决定原子间如何相互作用,是精度的根源;系综决定模拟处于何种热力学约束(如 NVE、NVT、NPT);时间步长决定积分稳定性,通常取原子振动周期的十分之一量级,对含氢体系常在 0.5~1 fs,对金属体系可到 1~2 fs。三者相互制约:势函数越“硬”,步长必须越小;温度越高,步长也需相应减小。

LAMMPS(Large-scale Atomic/Molecular Massively Parallel Simulator)是当前材料与力学工程领域使用最广的开源 MD 引擎之一。它的定位可以概括为:面向大规模并行、以势函数和命令脚本驱动、可灵活扩展的原子尺度计算平台。在材料工程中,它用于研究金属/合金的力学响应、位错与孪晶演化、界面与晶界行为、辐照损伤、热传导与相变;在力学工程中,它常被用来提取弹性常数、屈服强度、层间剪切强度、摩擦系数等本构参数,为有限元等连续介质模型提供输入。LAMMPS 本身不提供图形界面,而是通过输入脚本(input script)逐条命令构建体系、施加载荷、输出结果,这种“脚本即模型”的方式正是建模工程化的基础。

一个完整的 LAMMPS 建模工程,通常不是“写一个脚本跑一次”,而是围绕可复现、可验证、可扩展展开的闭环。其核心流程可概括为:明确物理问题 → 选择势函数与单位制 → 构建初始构型 → 设置系综与边界条件 → 弛豫平衡 → 施加载荷或扰动 → 采样与输出 → 后处理与判据分析。其中“弛豫平衡”常被初学者忽略:若初始构型能量过高或密度不合理,直接加载会得到虚假的力学响应。工程上一般先做能量最小化(minimize),再在目标温度下 NVT/NPT 弛豫足够长时间,直到温度、压力、势能围绕均值稳定波动,才进入生产阶段。

输入输出文件体系是 LAMMPS 工程的骨架。输入侧主要包括:in.* 输入脚本(命令序列)、data.* 数据文件(原子坐标、类型、盒子尺寸)、势函数文件(如 EAM 的 .eam、Tersoff 的 .tersoff、ReaxFF 的 ffield.reax)。输出侧主要包括:log.lammps 运行日志(热力学量随步数变化)、dump.* 轨迹文件(原子坐标/速度/受力)、thermo 屏幕与日志中的温度-能量-压力序列,以及自定义的 fix ave/time 统计文件。后处理常用 OVITO 观察结构演化、用 Python 脚本读取 dump 计算 RDF、MSD、应力-应变曲线等。

LAMMPS 建模工程:输入 → 计算 → 输出 → 判据 输入文件体系 in.* 输入脚本 data.* 原子构型 势函数文件 单位制 / 系综设置 边界条件 p p p (脚本即模型) LAMMPS 计算引擎 力场求力 F = -∂U/∂r 数值积分 Verlet / 步长 Δt 系综控制 NVE/NVT/NPT 弛豫平衡:minimize → NVT/NPT 稳定 生产阶段:加载 / 扰动 / 采样 并行分解 · 邻域列表 · 时间推进 输出文件体系 log.lammps 日志 dump.* 轨迹 thermo 热力学量 fix ave/time 统计 restart 重启文件 (轨迹即数据) 后处理与结果判读 结构演化 OVITO / 位错 / 相变 统计量 RDF / MSD / 应力应变 判据与验证 能量守恒 / 收敛 / 对照 关键判据:温度/压力稳定波动 · 势能收敛 · 结果可复现 学习路线:命令语法 → 单晶拉伸 → 缺陷与界面 → 多物理场耦合 反馈调参

从工程视角看,LAMMPS 建模的难点往往不在“会不会写命令”,而在“参数是否合理、结果是否可信”。例如,用 EAM 势模拟面心立方金属铜时,若初始晶格常数取 3.615 Å、温度设为 300 K,在 NPT 系综下弛豫 50 ps 后,体系应回到接近零压的平衡密度;若压力长期偏离零值,说明势函数或初始密度需要检查。再如单轴拉伸模拟,通常沿 z 方向施加应变,应变率取 10⁸~10¹⁰ s⁻¹(远高于实验值,这是 MD 的尺度限制),通过 fix deform 逐步拉伸,同时用 fix ave/time 输出应力,最终得到应力-应变曲线,从斜率提取杨氏模量、从峰值判断屈服强度。这些数值属于教学示例量级,实际应以具体势函数和体系为准。

学习建议:先跑通一个最小可复现案例(如 500 个原子的铜单晶弛豫),再逐步增加复杂度。每改一个参数,只改一处,并保留 log 与 dump,便于对比。把“输入脚本 + 势函数 + 数据文件 + 后处理脚本”作为一个工程单元管理,而不是散落的文件。

本课程后续将围绕这条主线展开:从 LAMMPS 命令语法与单位制、势函数选择与验证,到构型构建、系综与边界条件、加载与采样、后处理与结果判读,最终形成可复现、可扩展的建模工程能力。掌握这条链路,才能在材料与力学问题中把原子尺度的计算真正转化为可信的工程结论。