☰
LAMMPS fix surface/global 详解:为颗粒模拟搭建全局颗粒表面(GRANSURF 包)
2026/10/16 5:49:07 网站建设 项目流程
  • 科研
  • 科学计算
  • 高性能计算

【免费下载链接】lammps

Public development project of the LAMMPS MD software package

项目地址:https://gitcode.com/gh_mirrors/la/lammps
点击查看免费下载

本文以 LAMMPS 官方文档 fix_surface_global.rst 为主体,完整讲解fix surface/global命令的语法、input/model关键字与可选参数的全部含义,并结合 src/GRANSURF/fix_surface_global.cpp 的源码实现与 examples/gransurf 下的完整示例,说明如何在全局颗粒表面(三角形/线段)上施加胡克/赫兹/通用 granular 接触模型,实现颗粒-表面接触力计算、表面整体运动与dump image渲染。读完后你可以独立配置一个基于 STL 或分子模板文件的颗粒容器/螺旋喂料器等边界条件,并正确设置fix_modify move等进阶操作。

命令定位:全局表面与局部表面的区别

fix surface/global用于在颗粒(granular)模拟中把颗粒表面(granular surface)作为边界条件。颗粒表面由一组三角形(3d 体系)或线段(2d 体系)定义。LAMMPS 中颗粒表面分为两种风格,官方使用文档 Howto_granular_surfaces.rst 给出了总体使用指引:

  • global(全局):表面不存为粒子、不跨处理器分布,而是由本 fix 统一存储,每个处理器都保存全部表面的副本,并由该 fix 直接计算所有颗粒与表面之间的力。适合表面数量适中的场景;每个三角形/线段可以任意大,甚至可以和模拟盒的一个维度一样大。
  • local(局部):由 fix surface/local 命令处理,表面作为特殊粒子参与并行分布,适合表面数量极大的场景。

从源码结构看,fix surface/global注册在 GRANSURF 命名空间下(src/GRANSURF/fix_surface_global.h 中通过FixStyle(surface/global,FixSurfaceGlobal)宏注册),其构造函数强制要求原子属性中必须包含radius和omega(fix_surface_global.cpp),因为颗粒接触模型需要粒子半径与角速度。

命令语法总览

完整语法如下(ID 与 group-ID 的含义见 fix 命令文档):

fix ID group-ID surface/global input args ... model args ... keyword value ...
  • ID、group-ID:fix 的标准标识符与受作用的原子组;
  • surface/global:风格名;
  • input:可指定一个或多个,定义表面的几何来源;
  • model:可指定一个或多个,定义颗粒类型 × 表面类型的接触模型;
  • 末尾可追加零个或多个keyword value对。

input 关键字:表面几何的两种来源

input子句的形式为source source-args,source取mol或stl:

input args = source source-args source = mol 或 stl mol arg = template-ID template-ID = molecule 命令中定义的分子模板 ID, 该模板由一组三角形(3d)或线段(2d)构成 stl args = stype smol stlfile stype = 分配给 STL 文件中所有三角形的数值表面类型 smol = 分配给 STL 文件中所有三角形的数值分子 ID stlfile = 定义三角形集合的 STL 文件名

来源一:molecule 文件(mol)

若三角形/线段已经通过 molecule 命令读入,则source用mol,template-ID参数填写 molecule 命令使用的模板 ID。注意:

  • 一条molecule命令可以把多个分子文件挂到同一个 template-ID 上;每个分子文件必须定义三角形或线段,不能是原子。
  • 多个分子文件时,最终表面集合是所有文件中三角形/线段的并集。
  • 分子文件中的每条线段/每个三角形都会带有自己的 type 和 molecule ID,这是与 STL 的重要区别——mol 来源可以逐面指定类型。

来源二:STL 文件(stl)

STL(stereolithography)文件只定义一组三角形,不包含逐三角形的类型或标记,因此:

  • stype:赋给该文件全部三角形的表面类型(数值);
  • smol:赋给全部三角形的分子 ID(数值);
  • stlfile:STL 文件名,支持文本或二进制格式,命令会自动检测格式;
  • STL 文件不能用于 2d 模拟,因为它只能定义三角形(源码中 fix_surface.cpp 对dimension == 2直接报错)。

一个值得注意的实现细节:工程设计软件导出的 STL 常含退化三角形(顶点重复或共线、面积为零)。自 2Sep2026 版本起,读取时会跳过这类三角形,并在 0 号处理器打印带跳过数量的警告。源码位于 fix_surface.cpp:对每个新三角形计算两条边向量的叉积,若叉积全零则nskip++并跳过,循环结束后若nskip非零则输出Skipped {} degenerate triangle(s) in STL file {}警告。若需要手工编辑二进制 STL,可用仓库中 tools/stl_bin2txt.cpp 工具把二进制转为文本;文本 STL 示例可参考 examples/gransurf/surf.tri.stl(一个由 4 个 facet 拼成的平面示例)。

命令允许多次使用input关键字,每次的来源可以是mol或stl,最终使用的表面是所有 input 输入三角形与线段的并集(源码构造函数中以哈希表按坐标去重合并点集,见 fix_surface_global.cpp)。

连接性(connectivity):表面如何被"拼"在一起

所有 global 三角形/线段定义完成后,该命令会计算它们的连接关系:

  • 3d:两个三角形"相连"当且仅当它们属于同一分子,且共享一个角点或一条边(两个角点);
  • 2d:两条线段"相连"当且仅当属于同一分子且共享一个端点。

一对相连的表面若同时与某个颗粒接触,将按 Howto_granular_surfaces.rst 文档中定义的规则相互作用,以产生物理上合理的合力。表面之间不要求全部连成一个整体,可以表示多个独立物体的表面。颗粒只要与某个(或几个)三角形/线段接触(重叠)即参与相互作用;三角形的两面都会与颗粒作用,因此表面可以无限薄,例如搅拌器的桨叶。

源码层面,连接关系分别存入Connect2d/Connect3d结构及neigh_p1/p2(2d 端点)与neigh_e1/e2/e3(3d 边)、neigh_c1/c2/c3(3d 角点)等锯齿二维数组(见 fix_surface_global.h),构造函数在connectivity2d_global()/connectivity3d_global()之后调用check_molecules():若一组互连的三角形/线段的 molecule ID 不一致,LAMMPS 会打印警告(但不报错)。

model 关键字:颗粒类型 × 表面类型的接触模型

每个model子句作用于"一个或多个颗粒类型 × 一个或多个表面类型"的笛卡尔积:

model args = ptype stype fstyle fstyle-args ptype = 数值颗粒类型(也支持通配符形式) stype = 数值表面类型(也支持通配符形式) fstyle = hooke 或 hooke/history 或 hertz/history 或 granular hooke / hooke/history / hertz/history 的 args = Kn Kt gamma_n gamma_t xmu dampflag (limit_damping) Kn = 法向颗粒排斥弹性常数(力/距离或压力单位,见 pair gran/* 文档) Kt = 切向接触弹性常数(同上) gamma_n = 法向碰撞阻尼系数(1/时间 或 1/时间·距离) gamma_t = 切向碰撞阻尼系数(同上) xmu = 静摩擦屈服判据(0.0–1.0e4 之间的无量纲值) dampflag = 0 或 1,是否排除/包含切向阻尼力 limit_damping = 可选参数,防止出现吸引相互作用 granular 的 args = 与 pair_style granular 的 pair_coeff 命令相同语法

要点(均与 pair_style gran/* / pair_style granular 文档保持一致):

  • 通配符:ptype/stype可以是单个数值,或*、*n、n*、m*n形式。若无数字,*表示 1 到 N 的所有类型;*n表示 1..n(含);n*表示 n..N(含);m*n表示 m..n(含)。
  • 全覆盖检查:model 关键字必须为"每个颗粒类型 × 每个表面类型"的组合都指定相互作用,否则报错;若某个组合被重复指定,则采用最后一次指定。从源码看,表面类型的上界取smaxtype或输入表面中的最大 type(上限 1024,见 fix_surface_global.cpp 中MAXSURFTYPE/MAXSURFMOL常量),随后在types2model二维表中逐格检查,发现nullptr即报Fix surface/global type pair is missing a granular model错误(fix_surface_global.cpp)。
  • 力的物理含义:三角形/线段与颗粒间作用力与pair_style gran/*中的公式相同,只是取其中一个粒子半径与质量趋于无穷大(即平面)的极限:delta = radius - r为颗粒与表面的重叠量,m_eff取颗粒质量,有效接触半径取颗粒自身半径。重叠量与法向量的几何计算细节见 Howto_granular_surfaces.rst。
  • NULL 缺省:Kt写 NULL 时取缺省Kt = 2/7 * Kn;gamma_t写 NULL 时取缺省gamma_t = 1/2 * gamma_n(含义与单位与pair_style gran/*页相同)。

从源码结构看,model解析区分两条路径(fix_surface_global.cpp):

  1. classic 路径:fstyle为hooke/hooke/history/hertz/history时调用model->define_classic_model(...),可再跟一个limit_damping标记;
  2. granular 路径:fstyle为granular时按pair_coeff语法调用model->add_sub_model(...),可追加tangential、damping、rolling、twisting、heat等子模型,以及limit_damping标记。若指定了heat子模型,则heat_flag置 1,从而必须提供temperature关键字,否则构造函数直接报错(fix_surface_global.cpp)。

文档中给出的最简示例(对应上面两种 input 来源):

molecule lines surf.line fix 1 all surface/global input mol linesurf model 1 1 hooke 4000.0 NULL 100.0 NULL 0.5 1 fix 1 all surface/global input stl 1 1 object.stl model * 1 hooke/history 4000.0 NULL 100.0 NULL 0.5 1

该命令可以多次使用(定义多个 fix surface/global),但通常没必要;多次使用时,不同命令定义的表面之间不会按上述规则建立连接关系。

可选关键字(keyword/value 对)

关键字速查表

关键字取值含义默认值
smaxtype数值允许的最大表面类型input 中出现的最大表面类型
smaxmol数值允许的最大表面分子 IDinput 中出现的最大分子 ID
flatmaxangle(度)相连三角形/线段之间被视为"平面连接"的最大夹角1 度
temperatureTsurf(开尔文)表面温度;使用含 heat 模型的 granular model 时必需无
neighbornsq或bin表面-原子邻居搜索算法bin

各关键字的行为与源码印证:

  • smaxtype/smaxmol:默认取 input 表面中的最大类型/分子 ID。若之后要用fix_modify type/region(见下文)把表面类型改成比默认更大的值,需要先用这两个关键字把上界调大。源码中smaxtype不能小于 input 已出现的最大类型,且不得超过内部上限 1024(fix_surface_global.cpp)。
  • flat:设置阈值maxangle(角度制)。两条夹角不超过阈值的相连表面在接触模型中被当作"flat",即使颗粒同时接触两面也只施加一个力。默认 1 度;源码把角度阈值换算为flatthresh = 1 - cos(π·flat/180)存储([fix_surface_global.cpp](https://link.gitcode.com/i/3f9eedea2644f836adbfd1402970e65f#L66-L67, L277-L283)),并对取值做了 0–90 度的合法性检查。
  • temperature:值Tsurf为表面温度(开尔文)。若任一 granular model 含依赖表面温度的 heat 模型则必填,否则被忽略。源码中该值还允许v_前缀的变量形式(fix_surface_global.cpp)。
  • neighbor:决定表面与原子之间邻居列表的构建算法,与 neighbor 命令类似:bin(默认)把原子与表面按空间分箱,仅检查相邻 bin 内的组合;nsq则对所有原子-表面对做 O(N²) 式遍历。文档建议两种都测试:一般而言表面数量越多bin相对越划算;当表面数量非常大时,fix surface/local 可能性能更好。

完整实战示例:向三角形容器倾倒颗粒并旋转喂料器

仓库 examples/gransurf 目录提供了 global/local/surf 三种风格的对照示例。下面以 examples/gransurf/in.tri.gran.global 为例,逐步说明一次典型配置:

# pouring spheres into box made of triangles units lj atom_style sphere boundary f f f dimension 3 comm_modify vel yes region box block -7 7 -7 7 0 9 units box create_box 2 box # simply for dump_modify to color lines = green pair_style granular pair_coeff * * hooke 4000.0 100.0 & tangential linear_history 1000.0 1.0 0.5 & damping mass_velocity neighbor 0.3 bin neigh_modify delay 0 every 1 check yes molecule tris surf.tri fix 20 all surface/global input mol tris & model * * granular hooke 4000.0 100.0 & tangential linear_history 1000.0 1.0 0.5 & damping mass_velocity fix 2 all nve/sphere fix 3 all gravity 1.0 vector 0 0 -1

关键步骤解读:

  1. atom_style sphere提供 radius/omega 属性,满足fix surface/global的构造要求;
  2. molecule tris surf.tri先把三角形文件 examples/gransurf/surf.tri 读成模板;
  3. fix 20 all surface/global input mol tris model * * granular ...:用 mol 来源引用模板tris,模型部分对"所有颗粒类型 × 所有表面类型"使用granular风格,主接触模型hooke 4000.0 100.0,切向子模型linear_history 1000.0 1.0 0.5(对应 Kt、gamma_t、xmu),阻尼子模型mass_velocity——这套参数与颗粒-颗粒的pair_coeff保持一致,使表面接触行为与颗粒间接触物理上协调;
  4. 颗粒动力学由fix nve/sphere+fix gravity提供,fix ins pour在盒顶上方持续补充颗粒;
  5. 示例结尾展示了本文后续详述的fix_modify move用法:
run 5000 fix_modify 20 move * rotate 0 0 0 0 0 1 25 run 5000 fix_modify 20 move * none run 5000

即先纯倾倒 5000 步,再让全部表面绕 z 轴以周期 25 旋转(相当于旋转的锥形喂料器),最后用move * none停转。

同一目录下的 examples/gransurf/in.tri.gran.global.stl 是 STL 版对照,差异只在表面来源一行:

fix 20 all surface/global input stl 2 1 surf.tri.stl & model * * granular hooke 4000.0 100.0 & tangential linear_history 1000.0 1.0 0.5 & damping mass_velocity

注意 STL 来源时三角形类型统一设为 2、分子 ID 统一设为 1,而 mol 版由 examples/gransurf/surf.tri 文件逐行给出。两份输入在 examples/gransurf 中均附有 g++.1/g++.4 的完整运行日志(如log.17Jul26.tri.gran.global.g++.4),可直接用于核对输出。

dump image 渲染支持

该 fix 支持 dump image 的fix关键字,会把表面几何信息传给dump image一并渲染:

  • 2d 体系:表面是一组线段,渲染为圆柱;fflag1设定圆柱半径;
  • 3d 体系:表面是一组三角形,可渲染为 wireframe、实心三角形或两者兼有;fflag1选择渲染风格(1 = triangles,2 = wireframe,3 = both),fflag2设定 wireframe 线径。

颜色规则:默认按三角形所关联原子类型的配色(color style 为 "type" 或 "element" 时);color style "const" 时默认白色,可用dump_modify fcolor修改;透明度默认不透明,可用dump_modify ftrans修改。示例输入中的注释块给出了实际用法:

#dump viz all image 100 tri.gran.global-*.ppm type diameter & # size 800 600 zoom 1.5 shiny 0.4 fsaa yes view 70 10 box yes 0.025 & # region zcylinder white transparent 0.5 fix 20 const 1 0 #dump_modify viz pad 6 acolor 1 orange boxcolor midnightblue backcolor black & # backcolor2 white lights 0.2 0.9 0.45 0.9 fcolor 20 silver ftrans 20 1

从源码看,几何导出由FixSurfaceGlobal::image()完成(fix_surface_global.cpp):2d 时每个条目打Graphics::LINE标志并输出两个端点坐标,3d 时打Graphics::TRI标志并输出三个顶点坐标,同时携带表面 type 供配色。

fix_modify:move、type/region、mol/region

本 fix 为 fix_modify 命令定义了三个专属关键字(只在本文档页描述,不在 fix_modify 总页中):move、type/region、mol/region,均可多次使用。通用格式:

fix_modify fix-ID keyword values ...

move:让表面按预设方式运动

move的值为smol mstyle args:

move values = smol mstyle args smol = 表面分子 ID(支持逗号分隔与通配符形式) mstyle = none / linear / wiggle / rotate / transrot / variable none args = (无) linear args = Vx Vy Vz wiggle args = Ax Ay Az period rotate args = Px Py Pz Rx Ry Rz period transrot args = Vx Vy Vz Px Py Pz Rx Ry Rz period variable args = v_dx v_dy v_dz v_vx v_vy v_vz
  • smol可以是逗号分隔的多个表面分子 ID:sm,sm,...,sm,每个子项可为数值或*/*n/n*/m*n通配形式(语义同前:N为分子总数)。
  • 必须包含相连物体的全部分子 ID:若一个物体由 2 个以上分子的表面构成,而move没有列全它们的分子 ID,会破坏连接关系。文档明确说明LAMMPS 不检查这一要求;同理,多次move造成的运动若违反表面几何限制(见 Howto 文档),LAMMPS 也不会检查。
  • none风格用于关闭之前启用的运动(例如停止一个物体的旋转)。
  • 除none外的各风格语义及其参数与 fix move 相同。
  • variable风格比 fix move 受限:只允许 equal 风格变量(variable 命令定义),不允许 atom 风格变量;且若某方向的位移变量与速度变量均为 NULL,则该方向的位置/速度分量不变——这与 fix move 对单个粒子的说明不同。

源码印证:modify_param中解析move子句,variable风格的每个变量名都经variable->exists()与equal_flag检查(fix_surface_global.cpp);运动状态存于Motion结构列表,none会取消该组分子的运动并在全部停转时释放pointmove内存(fix_surface_global.cpp)。对 local 表面,同样的运动操作可以直接用fix move命令加一个包含三角形/线段粒子类型的 group 实现。

type/region 与 mol/region:事后重设表面类型/分子 ID

type/region values = stype region-ID stype = 数值表面类型 region-ID = 之前用 region 命令定义的几何区域 ID mol/region values = smol region-ID smol = 数值表面分子 ID region-ID = 之前用 region 命令定义的几何区域 ID

这两个关键字都用于在input完成初始化之后重设表面属性,对 STL 三角形特别有用——因为 STL 文件无法逐三角形指定 type 或 molecule ID。判定规则:

  • mol/region中smol必须在 1 到 maxmol 之间(maxmol 为 input 读入的最大分子 ID 或smaxmol关键字的设置值);
  • region 可以定义为几何物体的内部或外部(如球、block 的 inside/outside);
  • 用三角形/线段的几何中心点判断表面是否位于区域内,若在内则把其分子 ID/type 重置为给定值(源码中即region->match(xsurf[i][0],xsurf[i][1],xsurf[i][2]),见 fix_surface_global.cpp)。

文档给出的示例:

fix_modify 1 move 2 rotate 0 0 0 0 0 1 25 fix_modify 1 move 1,3*5,8* rotate 0 0 0 0 0 1 25 fix_modify 1 type/region 3 myBlock

重启、输出与 run 行为

  • 不写二进制 restart:本 fix 的信息不会写入 restart 文件,重启时需重新执行该命令;
  • 不产生全局量或 per-atom 量:没有可供 输出命令 访问的 compute 式数据;
  • 无 start/stop 参数:本 fix 的任何参数都不能用于 run 命令的start/stop关键字;
  • 能量极小化期间不激活:minimize过程中不调用本 fix。

限制与默认值汇总

限制:

  1. 属于GRANSURF package,只有在 LAMMPS 编译时启用了该包才可用,见 Build_package 页;
  2. 分子 ID 目前尚未真正参与颗粒表面相互作用(未来可能用到),其设计意图是唯一指派给每个互连的三角形/线段集合(相当于把每个物体看作一个"分子"),但这不是强制要求,LAMMPS 不检查;若一组互连表面的分子 ID 不一致会收到警告;
  3. 多个fix surface/global命令之间定义的表面上互不连接(见前文);
  4. STL 不可用于 2d;move的完整性与几何合法性不被自动校验。

默认值:smaxtype(smaxmol)= input 中所有表面的最大类型(分子 ID);flat= 1 度;neighbor=bin;temperature= 无。

相关命令:fix surface/local(局部表面)、pair_style surf/granular(以 pair 风格计算表面接触力)、fix smd/wall_surface(MACHDYN 包中的表面墙)。

  • 科研
  • 科学计算
  • 高性能计算

【免费下载链接】lammps

Public development project of the LAMMPS MD software package

项目地址:https://gitcode.com/gh_mirrors/la/lammps
点击查看免费下载

相关推荐

上一篇:BOOTSTRA.386 实战:Bootstrap 5.3 Link 链接工具类完全指南
下一篇:Kronos金融时序预测模型:重新定义金融市场的"语言"理解

创作声明:本文部分内容由AI辅助生成(AIGC),仅供参考

需要专业的网站建设服务?

联系我们获取免费的网站建设咨询和方案报价,让我们帮助您实现业务目标

立即咨询