lammps粘度计算

流体输运性质里，粘度是工业计算和材料设计常用的参数，但实验测高温或特殊体系的粘度往往很贵。我这个项目客户想算一种高温润滑流体的粘度随剪切速率的变化，我用了 LAMMPS 做非平衡分子动力学（NEMD）粘度计算，并把速度分布剖面和”粘度-剪切速率”的牛顿平台画了出来。

方法上有两条路：平衡态的 Green-Kubo（EK）和非平衡态的 NEMD。我这个项目要研究剪切稀化（非牛顿行为），必须用 NEMD。NEMD 的核心是给体系施加一个宏观剪切流，用 reverse peristaltic（反向多面体）或 SLLOD 算法驱动，让上下两半流体以相反方向运动，中间形成速度梯度。我建了一个长方体盒子，x 方向加周期边界、y 方向（或 z）做剪切，体系先 NVT 平衡 2 ns 再 NEMD 生产 4 ns 取平均。

第一步是看速度分布剖面。这是 NEMD 最直观的验证：在稳态剪切下，体系中间应该出现一条线性速度梯度（velocity profile），斜率就是剪切速率 γ̇。我画了沿剪切方向的速度剖面，确认它是直的——如果剖面弯曲或出现平台断层，说明体系还没到稳态或边界条件有问题。我这个项目里 2 ns 后就稳定成线性，说明平衡充分。

第二步是算粘度。粘度 η = τ_xy / γ̇，其中 τ_xy 是从应力张量算出的剪切应力分量（用 virial 公式），γ̇ 是施加的剪切速率。我扫了 5 个不同剪切速率（从 10⁸ 到 10¹⁰ s⁻¹，分子模拟的时间尺度决定了只能到这个范围），算出对应的 η，画成”粘度-剪切速率”曲线。

判读是重点：在低剪切速率区，η 基本是常数（牛顿平台），说明流体在该尺度表现为牛顿流体；当 γ̇ 超过某个阈值，η 开始下降，出现剪切稀化——这正是客户关心的高温润滑行为。我用这个结果解释了客户实验中”高速剪切下流体变稀、润滑膜变薄”的现象。

踩坑记录：第一，边界效应。NEMD 盒子太小，剪切流在边界附近会畸变，我用了足够大的盒子（> 40 Å 剪切方向）并保证中心区域线性。第二，热化。施加剪切后体系会生热，必须用 thermostat 控温，但 thermostat 不能干扰剪切流，我用的是沿流动方向的 NVT 而非全域。第三，应力涨落。NEMD 的 τ_xy 涨落大，生产跑短了误差能到 20%，我延长到 4 ns 并把数据每隔 100 ps 块平均。第四，与 Green-Kubo 交叉验证。NEMD 在零剪切极限应回到 EK 值，我用 EK 独立算了一个点，确认两者在平台区吻合，证明 NEMD 设置无误。

从我的工程经验看，lammps粘度计算 最有用的不是”给一个粘度数字”，而是画出完整的剪切速率依赖——告诉客户这个流体在哪些工况下还是牛顿的、哪里开始稀化。我给客户的建议是：用牛顿平台的粘度做常规 CFD 输入，但涉及高速工况必须用剪切稀化模型。更多流体输运计算，见 [lammps粘度计算](https://www.keyanxueshu.com/category/md/lammps/)；需要定制，[lammps粘度计算](https://www.keyanxueshu.com/) 上有入口。