LS-DYNA SPH超高速碰撞仿真：铝球穿靶K文件全流程解析

做超高速碰撞仿真LS-DYNA是绕不开的名字。尤其是弹体击穿铝靶板这类问题网格法经常算到一半网格就畸变、失效甚至负体积而SPH方法却能比较好地把碎片飞溅、穿孔演化模拟出来。我第一次用SPH做弹丸穿靶时最深的感受是理论不难难在把K文件里的每一个卡片配置明白。这篇文章从一个最简单的铝球撞击铝靶板模型讲起把SPH算法配置、K文件结构、材料参数、粒子间距选择、后处理判读全流程过一遍并给出一份可以直接改着用的K文件模板。适合刚接触LS-DYNA、想快速跑通超高速碰撞仿真的人。1. 超高速碰撞仿真为什么选SPH算法1.1 超高速碰撞的物理特征超高速碰撞一般指碰撞速度在2公里每秒以上的撞击过程弹体和靶板在微秒量级内经历强烈压缩、碎裂、甚至局部熔化。这个尺度下材料应变率可以到10的6次方每秒压力能到几十GPa远超出常规结构冲击的适用范围。常见铝合金在高应变率下表现出明显的强度提升和热软化失效模式也从静载下的韧性断裂变成绝热剪切带、层裂、崩落和碎片化。用传统Lagrange网格算这类问题非常吃力。网格单元在大变形下严重畸变时间步长会因为单元最小尺寸变小而不断下降最后出现负体积直接算崩。Erosion删除单元本质上是在人为“挖”材料删除之后动量守恒和能量守恒都会受影响碎片形态也不真实。Euler方法倒是能处理大变形但需要在空间中划分背景网格材料界面追踪复杂建模成本高。SPH被大量用于超高速碰撞正因为它的离散方式避开了这条死路。它不需要固定连接关系材料可以在粒子之间断开碎片飞散、穿孔开坑都是自然演化出来的结果不需要人为设定删除阈值。计算核心仍然是经典连续介质力学方程只是求解方法从网格差分变成了粒子插值。1.2 SPH算法的核心思想SPH的全称是Smoothed Particle Hydrodynamics中文一般叫光滑粒子流体动力学。它把连续介质离散成一组粒子每个粒子携带质量、密度、速度、应力、内能等物理量。粒子之间通过核函数进行加权插值某个粒子的物理量由其邻域内其他粒子的值按距离加权得到。所谓“光滑”就是用一个光滑长度h控制相互作用范围。你可以把粒子理解成一小团会变形的“面团”相邻面团之间有重叠重叠程度由光滑长度决定。计算中每个粒子不是和固定网格点交换数据而是动态搜索邻域粒子。正因为没有固定的单元拓扑材料随便怎么大变形、断裂、碎裂粒子都能重新排列不会出现单元畸变导致计算崩溃。在LS-DYNA中SPH粒子可以用节点和粒子单元定义核心关键字是SECTION_SPH后续用ELEMENT_SPH把节点转化为粒子。前处理阶段往往还是先画网格再把网格节点转成粒子或者直接用一个填充工具在几何体内生成粒子。求解时LS-DYNA会为每个粒子搜索邻居粒子求解连续性方程、动量方程和能量方程更新密度、速度和应力。1.3 什么时候用SPH什么时候用网格法SPH不是银弹。它的计算精度依赖粒子分布边界附近的粒子因为邻居缺失容易出现边界效应应力波在边界处的传播精度也不如成熟的有限元。光滑长度和粒子间距的选取对结果影响很大粒子数量一多计算量也非常吃紧。所以我自己在实际项目里不追求“全SPH”而是看问题定方案。低速冲压、结构碰撞、静态强度分析这些场景Lagrange网格仍然是最稳的选择。爆炸驱动飞片、超高速穿甲、水下爆炸附近的结构响应、射流切割这类大变形和强界面问题SPH有明显优势。还有一种工程上很常用的做法是混合建模弹体或冲击源用SPH远处或受力较小的结构保留FEM网格两者之间用接触耦合。这样既能保住SPH的大变形能力又能控制模型规模。对第一次做超高速碰撞的人来说我建议先从全SPH模型入手。模型简单关键字少粒子间相互作用自动计算容易跑通。跑通之后再尝试SPHFEM混合理解接触卡片的作用。2. 动手前先理清K文件的结构2.1 K文件是什么K文件是LS-DYNA的输入文件后缀通常是.k或.dyn本质是一堆关键字数据卡。每一个关键字对应一种定义比如材料卡、单元卡、接触卡、控制卡。求解器从上往下读卡按卡片的命令组装模型。一个K文件哪怕少了一个空格、多了一个逗号都可能造成模型读入失败或者更糟——参数错乱但能算结果却完全不对。所以强烈建议不要纯手写K文件去定义几千上万个节点。正确的做法是用LS-PrePost、HyperMesh或ANSA做前处理生成完整K文件骨架再手动修改关键卡片。比如调整材料参数、状态方程系数、加载速度这些手工编辑比点鼠标快得多。我自己习惯把同一套K文件复制多份分别做参数扫描比反复打开GUI操作高效很多。K文件的基本骨架包括*KEYWORD开头*TITLE标题*CONTROL_*控制参数*SECTION_*单元属性*MAT_*材料模型*EOS_*状态方程*PART零件定义*NODE节点*ELEMENT_单元然后是边界条件、初始速度、接触、输出控制最后END结束。理解这个顺序后你就能在几百行的K文件里快速定位问题。2.2 单位制最容易被忽略的致命错误LS-DYNA本身不规定单位所有数字必须自洽。很多人第一次算SPH穿靶结果靶板纹丝不动或者弹体直接消失在屏幕上八成是单位制出了问题。我最常用的是kg-mm-ms这套单位制碰撞冲击问题换算起来非常方便速度单位mm/ms刚好和m/s一一对应应力天然是GPa。物理量kg-mm-ms单位制长度mm质量kg时间ms速度mm/ms1 mm/ms 1 m/s密度kg/mm³应力/压力GPa力kN再看一组常用材料参数在这个单位制下的写法铝靶板密度2.70e-6 kg/mm³弹性模量68.9 GPa钢弹体密度7.85e-6 kg/mm³弹性模量210 GPa。Johnson-Cook模型里的A、B参数也必须换算成GPa状态方程里的冲击波速度C则用mm/ms比如铝合金C大约5.328 mm/ms对应实际物理速度5328 m/s。最容易错的就是把文献里g-cm-μs单位制下的参数原样抄进kg-mm-ms模型结果差出几个数量级。2.3 整体模型与装配思路模型设定为一个直径6mm的钢球弹体以2000m/s速度垂直撞击一块40mm×40mm×2mm的6061-T6铝合金靶板。靶板厚度方向是重点考察区域因为穿孔过程在这个方向上演化最快。弹体用SPH粒子靶板也用SPH粒子边界设置固定约束。整体模型可以看成“弹体团”对准“靶板团”以高速迎面接触。前处理阶段最简单的方式是在LS-PrePost里分别建立六面体网格然后使用SPH Generation工具把网格节点转成SPH粒子。转成粒子后原有网格单元可以不保留。弹体和靶板各自指定不同的Part ID这样后处理时可以通过Part区分粒子也能单独设置初始速度。也有人用TrueGrid或Python脚本直接生成粒子节点坐标再写K文件但入门阶段没必要。LS-PrePost生成节点和粒子单元后我们再打开K文件手工加入材料、状态方程、初始速度、固定边界和输出控制比从零手写稳得多。3. SPH模型搭建与K文件模板3.1 粒子间距与模型尺寸计算SPH粒子间距是决定计算精度和成本的核心参数。粒子间距太大冲击波空间分辨率不够穿孔直径偏大碎片形态粗糙粒子间距太小粒子数量爆炸式增长普通电脑跑不动。工程上有一个常见经验粒子间距不超过最小结构特征尺寸的十分之一厚度方向至少要有6到10层粒子。靶板厚度2mm如果取10层粒子粒子间距h0.2mm。靶板平面方向上40mm除以0.2mm就是200个粒子厚度方向10层靶板粒子总数约为200×200×1040万。弹体直径6mm体积约113mm³每个粒子体积0.008mm³弹体粒子约1.4万。加起来超过41万粒子对入门教学来说偏大四核笔记本会非常吃力。我建议教学模型把靶板缩小到30mm×30mm×2mm粒子间距取0.25mm厚度方向8层。这样靶板粒子数是120×120×8115200颗弹体粒子约8600颗总粒子数约12万普通台式机还能接受。如果你只是想验证流程可以把靶板再缩小到20mm×20mm×2mm粒子数降到5万左右几分钟能算完。计算粒子数量时有一点要注意SPH粒子的体积按h³估算但实际光滑长度和粒子间距有关。如果粒子分布不均匀局部粒子过密会导致时间步骤降。用六面体网格均匀生成粒子能最大程度避免这个问题。3.2 弹体和靶板的材料与状态方程超高速碰撞中材料在高压、高应变率、高温下同时工作单靠一个理想弹塑性材料远远不够。最常用的组合是Johnson-Cook强度模型加Gruneisen状态方程。Johnson-Cook描述屈服应力随等效塑性应变、无量纲应变率和温度的变化方程里包含A、B、n、C、m五个主要参数。Gruneisen状态方程描述材料在冲击压缩下的压力-体积关系把冲击波速度和粒子速度的线性关系引入计算适合冲击波传播问题。下面这组参数是教学默认值来源于公开文献但做了适当简化。实际工程项目中必须核对材料牌号、热处理状态和实验条件不能直接拿来用于最终设计。参数钢弹体近似6061-T6铝靶板近似密度 kg/mm³7.83e-62.70e-6剪切模量 GPa77.026.0弹性模量 GPa210.068.9泊松比0.300.33A GPa0.7920.265B GPa0.5100.426n0.260.34C0.0140.015m1.031.00参考温度 K293293熔化温度 K1793925冲击波速度C mm/ms4.5695.328S11.491.338Gamma02.172.00还需要在材料卡里考虑失效。Johnson-Cook材料卡可以定义D1到D5失效参数描述失效应变对压力、应力三轴度和应变率的依赖。如果不需要精确失效也可以在*MAT_ADD_EROSION里增加一个最大等效塑性应变阈值超过后粒子被删除。穿孔模拟一定要有失效机制否则靶板只会无限变形而不出洞。3.3 接触、边界和初始速度的设置经常有人问SPH粒子之间要不要设置接触这一点要分清模型类型。如果弹体和靶板全部由SPH粒子构成粒子与粒子之间的相互作用由SPH控制方程自动完成不需要额外定义接触关键字。粒子彼此靠近时会产生排斥力发生碰撞时通过状态方程和本构模型产生应力响应。如果采用SPHFEM混合模型比如弹体是SPH、靶板保留FEM网格那么必须在SPH粒子和FEM单元表面之间定义接触。常用*CONTACT_AUTOMATIC_NODES_TO_SURFACESPH粒子作为从面FEM网格作为主面。接触设置不当SPH粒子会直接穿透FEM网格看到“粒子哗啦啦全飞走”的诡异结果多半就是接触卡没写对。边界条件方面靶板边缘需要固定。SPH粒子的边界条件直接施加在边界粒子节点上用*BOUNDARY_SPC_SET定义节点集约束X、Y、Z三个方向的位移。要注意的是如果边界区域内粒子间距不均匀约束边界粒子时可能会产生初始应力建议边界粒子区域保持均匀排列。初始速度设置用*INITIAL_VELOCITY。弹体全部节点归入一个SET_NODE_LIST速度值按kg-mm-ms单位制写成2000.0方向为X正向。如果弹体初始位置和靶板之间有间隙速度会先让弹体飞行一小段距离再接触这符合实际穿靶过程。3.4 K文件模板下面给出一个最小化模板定义了材料、状态方程、SPH单元属性、Part、初始速度和输出控制。节点和*ELEMENT_SPH数据量太大模板里用注释代替实际建模时由LS-PrePost自动生成。*KEYWORD *TITLE SPH steel projectile impact on aluminum plate $ 单位制: kg-mm-ms; 应力单位: GPa; 速度单位: mm/ms *CONTROL_SPH $ ncbs boxi boxj boxk boxm boxn 1 *CONTROL_TERMINATION $ endtim endcyc dtmin endeng endmas 0.020 *CONTROL_TIMESTEP $ dtinit tssfac isdo tslimt dt2ms lctm erode 0.0 0.6 *SECTION_SPH $ secid icform iaxis icont 1 1 2 1 *MAT_JOHNSON_COOK $ 靶板 6061-T6 铝教学参数 $ mid rho G E PR A B N C M TM TR EPSO 1 2.70e-6 26.0 68.9 0.33 0.265 0.426 0.34 0.015 1.00 925 293 1.0 *EOS_GRUNEISEN $ eosid C S1 S2 S3 gama0 a E0 1 5.328 1.338 0.0 0.0 2.00 0.0 0.0 *MAT_JOHNSON_COOK $ 弹体 钢教学参数 2 7.83e-6 77.0 210.0 0.30 0.792 0.510 0.26 0.014 1.03 1793 293 1.0 *EOS_GRUNEISEN 2 4.569 1.490 0.0 0.0 2.17 0.0 0.0 *PART $ pid secid mid eosid 1 1 1 1 2 2 2 2 *SET_NODE_LIST $ sid 1 $ 弹体粒子节点号由LS-PrePost生成 $ 此处会列出弹体所有节点编号 *INITIAL_VELOCITY $ nsid nsidex icid 1 0 $ vx vy vz vxr vyr vzr 2000.0 0.0 0.0 0.0 0.0 0.0 *SET_NODE_LIST $ sid 2 $ 靶板边界节点号由LS-PrePost生成 *BOUNDARY_SPC_SET $ nsid cid dofx dofy dofz 2 1 1 1 *DATABASE_BINARY_D3PLOT $ dt lcdt beam npltc psetid 0.001 *END注意*BOUNDARY_SPC_SET里面需要确认DOF约束关键字在LS-DYNA对应版本中的格式老版本直接用DOFX1, DOFY1, DOFZ1。不同版本略有差异读入K文件后先查看d3hsp信息再确认没有报错。这个模板并不完整节点数据、粒子数据、单元编号都需要前处理生成后填入。4. 提交计算与穿孔结果后处理4.1 求解器参数设置与监控K文件准备好之后提交计算的方式很简单命令行下输入ls-dyna iimpact.k ncpu8 memory200m也可以用LS-Run或作业管理工具提交。SPH粒子模型的计算规模大ncpu尽量多分配。算例运行期间LS-DYNA会生成d3hsp文件这是文本格式的计算信息记录了模型概要、单位制识别、材料卡读取结果、初始时间步长和每一步的耗时。我每次提交后第一件事就是打开d3hsp检查有没有关键字被跳过材料参数有没有读入异常。终止时间设置为0.02ms也就是20微秒。对于直径6mm弹体撞击2mm厚靶板这个时间足够完成穿孔和碎片云初步扩展。时间步长受粒子间距和声速限制SPH模型的时间步通常非常小可能在纳秒量级。如果计算时间太长可以在*CONTROL_TIMESTEP里设置质量缩放因子但质量缩放会引入人为质量增加超高速碰撞这种高动态问题不建议用太激进的缩放。计算过程中可以用LS-PrePost实时监控d3plot文件查看粒子运动是否正常。如果一开始就出现粒子疯跑、速度暴涨说明参数有问题趁早停掉改K文件不要浪费时间跑完再排查。4.2 铝靶板穿孔的典型结果判读打开d3plot后第一眼看整体变形。正常穿靶过程会经历几个阶段弹体接触靶板表面靶板前表面产生压缩波和崩落碎片靶板背面开始鼓包内部出现层裂与微裂纹随后靶板厚度方向被击穿形成穿孔弹体和碎片继续向后飞散形成碎片云。SPH模型可以清晰看到这些阶段因为粒子可以分离、飞溅不会像网格单元那样僵硬地被删除。建议查看的物理量是Von Mises应力、等效塑性应变和材料状态。等效塑性应变可以反映材料失效区域失效粒子被删除后会形成孔洞。穿孔直径的量测可以直接在LS-PrePost里用测量工具取孔洞边界粒子之间的最大距离。如果孔径显著大于弹体直径尤其要怀疑粒子间距偏大或失效参数设置过松。残余速度也是穿靶分析的重要指标。在靶板后方设置一个截面或跟踪弹体中心粒子的速度得到剩余速度随时间变化曲线。弹道极限速度通常定义为刚好穿透靶板时对应的初始撞击速度。用LS-DYNA做穿甲分析输出弹体速度曲线后可以画一张“初始速度-剩余速度”曲线这是后处理里最有工程价值的图之一。4.3 常见问题与调试心得下面整理几种我实际踩过的问题按出现频率排序现象大概率原因解决方向所有粒子瞬间飞散单位制不一致或材料参数太小核对密度、弹性模量、速度单位靶板怎么都穿不透失效参数未设置或过保守检查JC失效参数或*MAT_ADD_EROSION粒子严重穿透靶板粒子间距太大或接触定义缺失纯SPH模型检查核函数参数混合模型检查接触时间步长一直掉粒子初始重叠或局部过密重新生成均匀粒子检查初始位置结果只在局部振荡光滑长度过小或邻居粒子数不足调整*CONTROL_SPH中的光滑长度比例K文件读入报错关键字格式错误优先用LS-PrePost生成再手改经验之谈不要一上来就追求0.1mm粒子间距。40mm靶板配0.1mm间距粒子数几百万笔记本散热风扇直接起飞。先用0.3mm间距跑通确认模型行为合理再加密到0.2mm比较。SPH本来就是计算密集型方法合理取舍精度和规模是基本功。另外一个容易被忽视的问题是光滑长度。*CONTROL_SPH里的光滑长度或接触搜索参数会影响邻居粒子数量和应力波传播。默认值在某些版本里偏保守可以适当调整但每次调整后都要重新看时间步和能量曲线。SPH对参数很敏感一次只改一个变量做好记录这是调试这类模型最有效的方法。5. 这个能力还能用到哪里5.1 从靶板穿孔到防护设计超高速碰撞仿真不仅仅是看一个弹丸打穿一块板。航天微陨石防护、装甲抗弹、爆炸容器碎片防护、汽轮机叶片外物撞击这些问题背后的核心物理是同一类冲击破坏行为。你学会了铝靶板穿孔建模等于掌握了一套通用的高速冲击问题建模框架。比如Whipple防护结构前板把高速弹体打碎形成碎片云后板再阻挡碎片云。前板穿孔和后板损伤的仿真用SPH比用网格容易很多关键是碎片云的扩散角、粒子速度和后板损伤分布。这些指标可以直接指导防护结构设计。参数扫描也很容易做改变弹体直径、撞击速度、入射角度、板厚和板间距批量跑K文件提取残余速度和穿孔直径形成设计曲线。还可以把材料换成钛合金、复合材料、装甲钢研究不同材料在高应变率下的失效差异。JC模型参数需要重新查文献失效参数尤其要用实验标定。SPH模型的粒子间距也要随材料特性调整脆性材料的失效形态和韧性材料差异很大模拟时不能一刀切。5.2 LS-DYNA同一个软件的不同玩法LS-DYNA除了拉格朗日网格和SPH粒子还有ICFD不可压缩流求解器、CESE可压缩流求解器、电磁求解器等。有意思的是同一个LS-DYNA里既能做高能冲击也能做完全不同的工程问题比如降落伞折叠展开和充气模拟。很多人以为降落伞仿真只能用专门的膜结构软件其实LS-DYNA通过控制体法或ICFD也能模拟伞衣展开过程伞绳、伞衣和流场的耦合还能用类似SPH粒子的处理手段。虽然这篇不写降落伞但理解SPH的无网格思想后再看那些流固耦合大变形仿真会轻松很多。回到超高速碰撞。我个人的体会是SPH入门最忌讳一直调参数却不跑通也忌讳一上来就图粒子越细越好。先做一个小模型哪怕只有几万粒子把穿孔过程完整跑出来再对照文献里的弹道极限速度。很多困惑其实是在结果里看出来的不是看书看出来的。多试几组速度你就知道不同撞击速度下穿孔形态差距有多大也会慢慢建立起“什么是合理结果”的感觉。如果你第一次做建议沿用本文的模型尺寸和参数跑完20微秒之后先看靶板是否出洞再看弹体是否破碎。有条件的话找一篇公开的SPH穿靶文献复制它的材料参数和实验条件把仿真结果和文献对一下。对上了再开始尝试调整粒子间距或靶板厚度。这个方法我用了很多年最稳。