2. 分子动力学基本原理

分子动力学（MD）把材料看成由原子核与电子共同构成的离散体系，在玻恩-奥本海默近似下，电子效应被“折叠”进原子间势函数，原子核则按经典牛顿第二定律运动。因此，MD 的核心可以概括为四件事：用什么方程描述原子运动、用什么势函数描述原子间相互作用、在什么统计系综下演化、以及如何用有限尺寸和离散时间步长近似真实材料。这四者共同决定了模拟能否稳定、结果是否可信。

一、牛顿运动方程：MD 的“引擎”

对第 i 个原子，位置为 r_i、质量为 m_i，其运动方程为：

m_i · d²r_i/dt² = F_i = -∇_i U(r_1, r_2, …, r_N)

其中 U 是体系总势能，F_i 是其他所有原子对原子 i 的合力。求解时通常采用“速度 Verlet”算法，它同时更新位置与速度，数值稳定且能量守恒性好：

v(t+Δt/2) = v(t) + F(t)·Δt/(2m)
r(t+Δt)   = r(t) + v(t+Δt/2)·Δt
计算 F(t+Δt)
v(t+Δt)   = v(t+Δt/2) + F(t+Δt)·Δt/(2m)

时间步长 Δt 必须小于最快振动周期的约 1/10。对含氢体系，O–H 伸缩振动周期约 10 fs，故 Δt 常取 0.5～1 fs；对金属体系可取 1～2 fs。Δt 过大将导致能量漂移甚至原子“飞散”，这是判断模拟是否失稳的第一信号。

二、原子间势函数：相互作用的“语言”

势函数决定 MD 能描述哪些材料与哪些性质。常见类型包括：

• Lennard-Jones（LJ）势：U(r)=4ε[(σ/r)¹²−(σ/r)⁶]，适合惰性气体与范德华体系，参数少、计算快，但无法描述金属键与共价键。
• 嵌入原子法（EAM）：把能量写成嵌入能加对势，能较好描述金属的弹性常数、层错能与空位形成能，是金属 MD 的主力势。
• Stillinger-Weber、Tersoff、AIREBO：分别面向硅、共价晶体与碳氢体系，能描述键角与配位变化。
• 反应力场（ReaxFF）：允许键断裂与生成，适合燃烧、催化与界面反应。

选择势函数时，应优先核对它是否复现了与目标性质相关的关键量，例如晶格常数、弹性模量、表面能或缺陷能。势函数不是“越复杂越好”，而是“越匹配目标性质越好”。

三、系综理论：控制温度与压力

真实材料处于特定热力学条件下，MD 通过系综来对应这些条件。常用系综包括：

控温的本质是让体系与“热浴”交换能量。Berendsen 控温温和但采样不严格；Nosé-Hoover 能产生正确的正则分布，但参数需调节，否则温度会振荡。实际流程常先 NVT 弛豫 50～100 ps，再切换 NPT 平衡 100～200 ps，最后进入 NVE 或 NVT 采样。

四、边界条件：把“小块”变成“大块”

MD 体系通常只有几万到几百万原子，远小于真实试样。为消除表面效应，最常用周期性边界条件（PBC）：原子从一侧离开时从对侧进入，体系像无限重复的晶格。PBC 下必须使用最小镜像约定，即只计算原子与最近镜像的相互作用，截断半径 r_c 必须小于盒子最短边长的一半。

对纳米线、薄膜、裂纹尖端等，则采用非周期/自由边界，让表面真实存在；对冲击、拉伸等，常采用固定边界或柔性边界，在边界层施加约束或恒力。边界选择直接决定能否观察到表面重构、位错形核等关键现象。

五、时间积分与稳定性判据

积分算法的误差主要来自截断误差与舍入误差。速度 Verlet 的局部误差为 O(Δt³)，全局误差为 O(Δt²)，且时间可逆、能量长期漂移小，是 MD 的标准选择。判断模拟是否健康，可观察：

• 总能量在 NVE 下应仅小幅振荡，若持续上升或下降，说明 Δt 过大或势函数截断处理不当。
• 温度在 NVT 下应围绕目标值波动，波动幅度约与 1/√N 成正比。
• 均方根位移（MSD）随时间线性增长，说明体系处于扩散或塑性流动状态；若 MSD 突然跳变，可能发生了原子重叠或“爆炸”。

把上述五个环节串起来，就形成一次完整的 MD 模拟：给定初始构型与速度，选定势函数与系综，施加边界条件，用速度 Verlet 逐步推进，最后从轨迹中提取能量、温度、均方根位移、径向分布函数与应力等观测量。若能量漂移、温度失控或原子重叠，应优先检查 Δt、截断半径与初始构型，而不是盲目延长模拟时间。理解这五个环节的物理含义与相互制约，是后续学习纳米材料力学行为模拟的基础。