资讯动态

LAMMPS compute temp/sphere 详解:球形粒子平动与转动温度的计算原理与实战用法

发布时间:2026/10/17 6:08:05 来源:尧图企业网站定制
科研科学计算高性能计算【免费下载链接】lammpsPublic development project of the LAMMPS MD software package项目地址https://gitcode.com/gh_mirrors/la/lammps点击查看免费下载导读本文围绕 LAMMPS 官方文档 compute_temp_sphere.rst 展开系统讲解compute temp/sphere命令的语法、物理公式、自由度DOF处理、bias/dof两个可选参数以及它与压力张量、热浴 fix、compute_modify的配合方式。读完本文你将理解有限尺寸球形粒子温度的计算公式与源码实现路径能够为atom_style sphere体系正确配置温度计算、去除速度偏置并接入恒温器还能依据单元测试样例核验自己的输入脚本。一、功能定位从点粒子温度到球形粒子温度LAMMPS 默认的 compute temp 把原子视为点粒子point particle只统计平动动能。而compute temp/sphere专门面向有限尺寸球形粒子同时计入平动动能与转动动能两部分贡献compute ID group-ID temp/sphere keyword value ...例如文档给出的三个典型用法compute 1 all temp/sphere compute myTemp mobile temp/sphere bias tempCOM compute myTemp mobile temp/sphere dof rotate从源码看该命令注册为ComputeStyle(temp/sphere, ComputeTempSphere)定义于 compute_temp_sphere.h构造函数位于 compute_temp_sphere.cpp。它设置了scalar_flag vector_flag 1、size_vector 6、tempflag 1表明这是一个温度类 compute既输出一个全局标量温度也输出一个长度 6 的全局向量对称张量并可以被fix nvt、fix npt、fix temp/rescale等所有需要温度输入的命令使用。何时必须使用 temp/sphere体系使用atom_style sphere粒子带半径、质量/密度、角速度 omega与力矩 torque你需要把粒子的自旋转动自由度纳入温度与热浴的统计你希望温度张量中的转动分量参与压力张量计算。与之对比compute temp假想所有粒子为点粒子即使体系中有旋转也不会统计转动动能。二、语法与关键字详解完整语法来自 compute_temp_sphere.rstcompute ID group-ID temp/sphere keyword value ...ID、group-ID与其他 compute 一致见 compute.rsttemp/sphere本命令的 style 名称可附加零个或多个keyword/value对。支持的选项关键字取值含义biasbias-ID指定一个去除速度偏置的温度 compute 的 IDdofall统计平动 转动自由度的温度默认值dofrotate只统计转动自由度的温度默认值无 biasdof all。源码中的参数解析逻辑很直白compute_temp_sphere.cpp枚举常量enum { ROTATE, ALL }对应两种模式mode默认取ALL遇到未知关键字直接报错Unknown compute temp/sphere keyworddof后接既不是rotate也不是all时同样报错。一个容易被忽略的细节是当mode ROTATE时构造函数会把extra_dof强制置 0compute_temp_sphere.cpp因为只统计转动温度时不应再扣除默认的平动自由度修正项。三、温度计算公式平动 转动动能与自由度compute temp/sphere的温度遵循经典关系$$T \frac{2 E_{\rm kin}}{N_{\rm DOF}, k_B}$$其中 $k_B$ 为玻尔兹曼常数$N_{\rm DOF}$ 为组内贡献温度的自由度数。平动部分与 compute temp 完全相同转动部分采用 $\frac12 I \omega^2$$I$ 为球体转动惯量$\omega$ 为粒子的角速度角速度矢量。文档明确给出了各情况下的自由度约定点粒子半径 0.0不转动3d 下只有 3 个平动自由度3d 球形粒子每个粒子 6 个自由度3 平动 3 转动2d 球形粒子每个粒子 3 个自由度2 平动 1 转动。源码中的 DOF 统计dof_compute()是自由度统计的核心函数compute_temp_sphere.cpp其逻辑完全对应上述约定3d半径 0 的点粒子计 3仅ALL模式有限尺寸球体在ALL模式下计 6、ROTATE模式下计 32d点粒子计 2球体在ALL模式下计 3、ROTATE模式下计 1各进程局部计数后经MPI_Allreduce汇总为全局dof若使用了 bias会再减去被偏置 compute 移除的自由度最终dof - extra_dof fix_dof其中fix_dof是fix shake、fix rigid等约束分子运动的 fix 扣掉的自由度换算因子tfactor force-mvv2e / (dof * force-boltz)compute_temp_sphere.cpp若dof 0则tfactor 0.0并在运行时检测dof 0报错。转动惯量的实现细节转动动能求和的关键代码位于 compute_temp_sphere.cppt (omega[i][0]*omega[i][0] omega[i][1]*omega[i][1] omega[i][2]*omega[i][2]) * INERTIA*rmass[i]*radius[i]*radius[i];其中常量INERTIA 0.4compute_temp_sphere.cpp是均匀实心球体转动惯量 $I \frac{2}{5} m r^2 0.4, m r^2$ 的系数。文档还特别提醒在 2d 模型中粒子仍按球sphere而非圆盘disk处理转动惯量与 3d 相同。⚠️自由度假设的注意点上述 DOF 约定建立在所有有限尺寸粒子都能自由转动、充分采样全部转动自由度的假设之上。如果通过某种势函数或 fix 施加了不产生力矩、或约束了部分/全部粒子转动实际自由度会变少此时应使用 compute_modify 的extra/dof选项手动修正。compute_modify默认的extra/dof为 2d 取 2、3d 取 3对应整体平动不变性的修正支持传负数来增加自由度compute temp/asphere文档中就有相应示例。四、六分量对称张量动能张量与压力张量除了标量温度compute temp/sphere还输出两个长度 6 的全局向量分量顺序均为 $xx, yy, zz, xy, xz, yz$。第一个是动能张量公式与温度公式一致只是把 $v^2$、$\omega^2$ 替换为对应分量的交叉项例如 $xy$ 分量用 $v_x v_y$、$\omega_x \omega_y$。第二个是供compute pressure使用的对称张量公式与 $E_{\rm kin}$ 相同但不包含 1/2 因子且同样以速度/角速度交叉项代替平方项。因此这些分量是传统动能张量分量的两倍——这正是压力张量定义所要求的。源码中compute_vector()compute_temp_sphere.cpp展示了具体实现每个粒子的平动贡献 $m, v_i v_j$、转动贡献 $I, \omega_i \omega_j$ 累加到 6 个桶中其中转动部分同样使用inertiaone INERTIA * rmass[i] * radius[i] * radius[i]ROTATE模式下只累加转动分量跨进程MPI_Allreduce后整体乘以force-mvv2e换算到能量单位。该张量可以直接喂给 compute pressure 参与压力/应力张量计算。五、bias 选项去除平动速度偏置bias关键字让compute temp/sphere在计算温度前先通过另一个温度 compute 将原子速度中的偏置如流动速度剖面、质心漂移移除只统计剩余的热速度。典型用法compute tempCOM all temp/com compute myTemp mobile temp/sphere bias tempCOM支持作为 bias 的典型 compute 包括compute temp/com、compute temp/region、compute temp/profile、compute temp/ramp、compute temp/deform等。这些速度偏置类 compute 的介绍见 Howto_temperature.rst 与 Howto_thermostat.rst——该文档明确指出temp/sphere与temp/asphere本身不直接计算偏置而是通过这个可选参数间接复用其他温度 compute 的偏置移除能力。源码init()compute_temp_sphere.cpp对 bias 做了严格校验按 ID 查找偏置 compute找不到则报错Could not find compute ... for temperature bias偏置 compute 必须满足tempflag ! 0能算温度且tempbias ! 0能算速度偏置偏置 compute 的组必须与本 compute 的组一致否则报Bias compute group does not match compute group还针对加速变体做了特殊处理若偏置 compute 是temp/region的加速版本如temp/region/kk统一识别为tempbias 2走逐粒子统计被移除自由度的路径避免dof_remove(-1)索引错误。remove_bias()/restore_bias()系列方法compute_temp_sphere.cpp均直接转发给偏置 compute 的对应接口保证热浴在更新速度前后正确地去掉偏置 → 施加热噪声 → 恢复偏置。六、输出信息类型、量纲与取值方式输出类型属性单位全局标量温度intensive强度量标量温度单位见 units.rst全局向量6 分量对称张量extensive广延量向量索引 1–6能量单位因此标量温度与原子数无关、可用于任意温度类命令向量分量则随粒子数累加适合作为compute pressure的输入或通过compute reduce等命令聚合。全局输出机制的总览可参考 Howto_output.rst。另外文档要求组内原子数在运行期间视为恒定若体系存在增删原子如fix pour、fix deposit、fix gcmc或原子逃逸出模拟盒应使用 compute_modify 的dynamic/dof yes让自由度数在每次计算时重新统计。七、使用前提与限制粒子必须存储力矩torque、角速度omega与半径radius即必须使用 atom_style sphere。文档还提到若要让粒子因力矩而转动需要使用atom_style hybrid sphere与带/sphere后缀的积分器如fix nve/sphere、fix nvt/sphere、fix npt/sphere见 atom_style.rst 中关于 sphere 风格的说明。组内所有粒子必须是有限尺寸球体或半径 0.0 的点粒子不能混入椭圆体等非球体粒子椭圆体请改用compute temp/asphere。源码在构造阶段即做硬性检查compute_temp_sphere.cppCompute temp/sphere requires atom attribute omega Compute temp/sphere requires atom attribute radius若输入脚本的atom_style不含 sphere 属性命令会直接报错终止。八、加速变体与相关命令加速变体该 compute 提供 KOKKOS 加速版本temp/sphere/kk文档首页即标注Accelerator Variants: temp/sphere/kk可在启用 KOKKOS 包时通过-sf kk自动启用。相关命令compute temp点粒子平动温度compute temp/asphere椭球非球形粒子的平动 转动温度温度计算与热浴的整体讨论见 Howto_thermostat.rst球形/椭球粒子的积分器、热浴与相关 compute 的完整清单见 Howto_spherical.rst。该文档特别指出默认thermo_style的热力学输出只统计平动自由度若体系全部由有限尺寸粒子组成、希望温度/压力输出包含转动需要用thermo_modify指定本 compute见 thermo_modify.rst。只对转动自由度施加随机热浴时可使用fix langevin的omega或angmom选项见 fix_langevin.rst。九、实战示例与单元测试佐证最小配置片段units lj atom_style sphere region box block 0 20 0 20 0 20 create_box 1 box create_atoms 1 single 5.0 5.0 5.0 set group all diameter 1.0 density 2.0 # 赋予半径与质量 velocity group all set 0.5 0.0 0.0 set group all omega 0.4 -0.3 0.7 # 赋予角速度 compute tT all temp/sphere # 平动转动温度 compute tR all temp/sphere dof rotate # 仅转动温度 thermo_modify temp tT # 热力学输出改用含转动的温度 fix 1 all nvt/sphere temp 1.0 1.0 0.5 # 球形粒子恒温器 run 100单元测试中的数值印证仓库的力风格回归测试直接给出了该 compute 的可复现数值compute-temp_sphere.yaml前置条件atom_style spherepair_style zero 4.0后置命令对全部 27 个原子set diameter 1.0 density 2.0、set omega 0.4 -0.3 0.7并对原子 1–8 单独改写 omegacompute test all temp/sphere的期望全局标量温度为0.90307859372518466 分量向量各分量也有精确期望值。同一目录下的两个变体测试进一步验证参数行为compute-temp_sphere_dof_rotate.yamldof rotate模式下只统计转动动能期望温度大幅下降至0.029554242000437315compute-temp_sphere_region_bias.yaml先定义region half block INF INF INF INF 6.75 INF与compute tbias all temp/region half再以compute test all temp/sphere bias tbias验证偏置移除路径期望温度0.9879111478197116。这些 YAML 既是unittest体系的回归基准也可以当作核对自己输入脚本数值的现成参考。十、常见问题排查现象原因与对策Compute temp/sphere requires atom attribute omega/radius未使用atom_style sphere或 hybrid 中未含 sphere 属性改用 sphere 风格Bias compute group does not match compute groupbias指向的 compute 分组与本 compute 分组不一致调整 group 定义Temperature compute degrees of freedom 0自由度被过度扣除extra_dof、fix 约束或 bias 扣除超过实际值检查compute_modify extra/dof设置温度明显偏高/偏低粒子实际并未自由转动被势函数或 fix 约束自由度数应小于默认假设用compute_modify extra/dof修正或确认是否需改用dof rotate只统计转动热力学输出温度不包含转动默认thermo_style用只含平动的温度需thermo_modify temp指向本 compute见 thermo_modify.rst 与 Howto_spherical.rst结语compute temp/sphere是 LAMMPS 中把有限尺寸球形粒子的自旋纳入温度统计的标准入口通过统一的 $\frac12 m v^2 \frac12 I\omega^2$ 动能表达、按维度与半径精细统计的自由度、以及可选的bias/dof模式它能够无缝接入fix nvt/sphere、fix npt/sphere等球形粒子积分器与各类温度类输出。配合 compute_modify 的动态自由度修正与 unit test 数值基准开发者可以快速定位温度计算中的自由度假设偏差写出物理上自洽的球形粒子模拟脚本。赞分享科研科学计算高性能计算【免费下载链接】lammpsPublic development project of the LAMMPS MD software package项目地址https://gitcode.com/gh_mirrors/la/lammps点击查看免费下载相关推荐LAMMPS compute temp/asphere 完全指南非球形粒子的平动与转动温度计算LAMMPS compute temp/asphere 完全指南非球形粒子的平动与转动温度计算 compute temp/asphere 是 LAMMPS 为科研科学计算高性能计算LAMMPS compute temp/body 命令详解体粒子body particle平动与转动自由度的温度计算LAMMPS compute temp/body 命令详解体粒子body particle平动与转动自由度的温度计算 compute temp/body科研科学计算高性能计算LAMMPS compute erotate/sphere球型粒子旋转动能的计算原理、实战用法与源码剖析LAMMPS compute erotate/sphere球型粒子旋转动能的计算原理、实战用法与源码剖析 本文围绕 LAMMPS 的 compute erot科研科学计算高性能计算上一篇Newsboat RSS阅读器配置详解从基础到高级技巧下一篇自定义数据集训练OneFormer指南无全景标注也能玩转实例分割创作声明:本文部分内容由AI辅助生成(AIGC),仅供参考

读完文章,也想定制专属网站?

尧图设计师 24 小时内与您沟通定制方案

免费获取报价 →
↑