拓冰网站定制 
LS-DYNA这个名字很多朋友是从安全气囊展开、降落伞充气这类强非线性问题开始接触的。同一个求解器里还有一个专门对付超高速碰撞的算法家族SPH就是其中最常用的一种。今天这篇不聊理论推导直接用一块40mm见方、3mm厚的铝靶板和一枚4mm铝球在LS-DYNA里把3km/s垂直撞击穿孔的SPH仿真完整跑一遍从算法原理讲到K文件每个关键卡片最后给出一套可以改参数直接用的模板。适合已经会跑一般显式动力学模型、但对SPH还停留在“听说过”阶段的读者。先说清楚一个基本判断如果你只是做常规低速冲击比如50m/s的落锤或汽车碰撞传统Lagrange网格完全够用。但一旦弹丸速度到了2km/s以上材料行为、失效模式、靶板响应都和低速工况完全不同。这时候SPH不是炫技而是绕开网格畸变的成熟工程方案。1. 超高速工况下网格为什么“扛不住”SPH的适用边界1.1 传统Lagrange网格在3km/s下的三个死穴低速冲击里结构变形是主要的能量吸收方式单元虽然会畸变但通过自适应网格或侵蚀失效还能继续算。到了3km/s这个量级弹靶接触区的压力可以到几十GPa甚至上百GPa材料强度相对于流体静水压力已经成了次要项——换句话说金属在这个状态下表现得像流体。Lagrange网格这时候会遇到三件麻烦事单元畸变弹丸和靶板接触区的单元被极度压扁Jacobian变成负值计算直接崩掉。即使不崩精度也完全不可信。侵蚀删除导致质量丢失很多低速碰撞模型用*MAT_ADD_EROSION删单元但超高速撞击中被删除的单元里包含大量动量和能量删多了弹道极限、剩余速度全都不对。接触穿透高速下节点速度太快接触搜索步长跟不上出现节点穿透后后续结果基本报废。ALE方法能解决一部分畸变问题但超高速撞击伴随着自由表面、碎片云、裂纹扩展ALE网格在处理这些不连续现象时同样吃力。1.2 SPH把连续体拆成“一堆有质量的粒子”SPH的核心思想不复杂把连续的弹丸和靶板离散成一系列粒子每个粒子携带质量、速度、密度、压力、应力这些物理量粒子之间通过一个叫“光滑核函数”的权函数相互作用。你可以把它理解成把一整块豆腐切成很多小丁。每个小丁自带速度和密度撞击的时候小丁可以被撞飞、碎裂、重新堆积但不存在“网格扭曲”这回事。因为SPH的计算不依赖固定的单元拓扑粒子之间是“邻居关系”而不是“连接关系”粒子靠近了就算邻居飞远了就断开联系。这个特性让它天然适合碎片云、孔洞、开坑这类强不连续问题。也正因为如此SPH能捕捉到很多有限元方法很难模拟的现象着靶瞬间的冲击波传播、靶板背面崩落spall、碎片云的扩散形态。1.3 SPH不是万能药三个需要提前知道的短板SPH在超高速撞击里很好用但它有代价。拉伸不稳定材料处于拉伸应力状态时SPH粒子可能自发聚团产生非物理的断裂。这在低压区尤其明显需要通过足够的粒子分辨率、合理的核函数选择来压制。边界缺陷靠近自由表面的粒子其核函数积分区域被截断密度和压力计算会偏低。这不是bug而是SPH数学性质决定的。工程上通常靠加密表面附近粒子或使用修正光滑函数来缓解。计算成本SPH需要每一时间步搜索邻居粒子粒子间距减半粒子数量翻8倍成本涨得非常快。4万粒子还属于小规模百万粒子级别的模型就要认真考虑HPC资源了。所以SPH最适合的场景就是大变形、碎裂、自由表面演化、材料失效模式复杂。超高速碰撞正好全占。方法网格畸变碎片/裂纹计算成本适用场景Lagrange严重靠单元删除近似低低速冲击、结构碰撞ALE可接受需要界面重构中高流体结构耦合、爆炸SPH不涉及自然模拟高超高速碰撞、碎裂、穿甲2. 弹靶模型设计与单位制先把4cm见方的铝板量成K文件里的数字2.1 几何尺寸与粒子间距这个演示模型我刻意做成一个能跑得动又不会失真太多的规模靶板40mm × 40mm × 3mm铝合金弹丸直径4mm的铝球撞击速度3km/s垂直入射粒子间距0.5mm在LS-DYNA里我习惯用cm-g-μs单位制。40mm就是4.0cm3mm就是0.3cm弹丸半径0.2cm粒子间距0.05cm。靶板体积4.0×4.0×0.34.8cm³除以单个粒子占有的体积0.05³0.000125cm³靶板大约需要38400个粒子弹丸球体积约0.0335cm³折算下来粒子数大约270个。整个模型不到39000个粒子普通笔记本用多核求解器跑几十微秒的物理时间半小时内能出结果。粒子间距这个参数直接决定精度和成本。间距从0.5mm加密到0.25mm粒子数变8倍但孔径和剩余速度的收敛性会明显改善。先跑0.5mm的摸底看趋势正常再做加密。2.2 单位制换算最容易出错的0.5328cm-g-μs单位制下的关键换算关系很多第一次接触SPH的人都会在这里栽跟头物理量单位长度cm质量g时间μs速度cm/μs压力/应力Mbar1Mbar100GPa密度g/cm³几个实测参数换算过来就是3km/s 弹速 0.3cm/μs铝的声速约5328m/s 0.5328cm/μs铝合金弹性模量70GPa 0.7Mbar6061-T6铝合金的初始屈服强度324MPa 0.00324Mbar如果你把324填成0.324意味着屈服强度放大了100倍靶板会硬得像一块不可穿透的装甲如果你把0.5328填成5328冲击波声速比实际大了四个数量级临界时间步会急剧变小计算卡到怀疑人生。2.3 材料模型为什么Johnson-Cook要配EOS超高速碰撞中材料强度仍然重要但压力—密度关系同样关键。所以材料卡片要分两层第一层是强度模型我用Johnson-Cook。它把流动应力拆成初始屈服、应变强化、应变率强化和温度软化四部分。对6061-T6铝合金常用参数如下参数数值说明A0.00324 Mbar初始屈服应力B0.00114 Mbar应变强化系数n0.42应变强化指数C0.002应变率敏感系数m1.34温度软化指数TM925K熔化温度TR294K室温第二层是状态方程EOS我用Mie-Grüneisen也就是经典的冲击Hugoniot关系。它的作用是把密度和压力联系起来让材料在强烈压缩时表现出正确的流体动力学行为。铝合金参数大约取C0.5328cm/μsS11.338GAMMA02.0。纯粹的低速碰撞模型可以只给材料强度、不给EOS但超高速碰撞不给EOS压力峰值会差得离谱整个仿真没有意义。3. 核心K文件逐段拆解SPH控制、Johnson-Cook和Grüneisen状态方程3.1 K文件模板总览下面是这个模型的骨架模板。节点和SPH粒子列表需要你用LS-PrePost从六面体网格转换生成或者用脚本按规则网格生成后粘贴到*NODE和*ELEMENT_SPH段。控制卡片本身可以直接抄。*KEYWORD *TITLE SPH hypervelocity impact demo $ Units: cm, g, us, Mbar *NODE $ 节点列表由脚本或LS-PrePost生成后粘贴到此处 $ NID X Y Z TC RC ... *ELEMENT_SPH $ 粒子列表生成后粘贴到此处 $ EID PID N1 ... *PART $ PID SECID MID EOSID HGID GRAV ADPORT OPTID $ projectile part 1 1 1 1 $ target plate part 2 1 1 1 *SECTION_SPH $ SECID ELFORM CST ENT 1 1 0 0 $ HMIN HMAX IVIS START MAXV CONT 0.025 0.10 1 *MAT_JOHNSON_COOK $ MID RO G E PR DTF VP RATE 1 2.700 0.269 0.700 0.330 $ A B N C M TM TR EPSO 0.00324 0.00114 0.4200 0.0020 1.3400 925.0 294.0 1.0e-6 $ PC SPALL IT -0.003 1 *EOS_GRUNEISEN $ EOSID C S1 S2 S3 GAMAO A E0 1 0.5328 1.338 0.0 0.0 2.0 0.77 0.0 $ V0 1.0 *INITIAL_VELOCITY_GENERATION $ ID STYP OMEGA VX VY VZ IVATN ICID IVAT 1 2 0.0 0.0 0.0 -0.3 *CONTROL_TERMINATION $ ENDTIM ENDCYC DTMIN ENDNEG ENDMAS 20 *CONTROL_TIMESTEP $ DTINIT TSSFAC ISDO TSLIMT DT2MS LCTM ERODE MS1ST 0 0.900 *CONTROL_SPH $ NCBS BOXID DT ICRM FORM START MAXV CONT 0 0 0 0 0 0 0 0 *DATABASE_BINARY_D3PLOT $ DT LCDT BEAM NPLTC PSETID 0.05 *DATABASE_GLSTAT $ DT 0.10 *END节点生成我多说一句。最简单的方式是在LS-PrePost里先建一个规则六面体网格然后选中网格执行SPH粒子生成功能程序会自动把每个节点转成一个SPH粒子并写出对应的*ELEMENT_SPH卡片。手动写K文件时每个SPH元素只需要填EID、PID和参与计算的第一个节点N1后面的N2到N8留空。3.2 SPH相关卡片逐条说明*PART卡片里弹丸和靶板用同一个材料ID、同一个EOSID但PID必须分开。这样后处理里才能分别观察弹丸残体和靶板碎片的运动也方便后续给弹丸和靶板赋予不同材料。*SECTION_SPH是SPH的核心配置之一。SECID对应part里的section IDELFORM是SPH的粒子公式编号。初学者不要乱改这个编号不同公式对核函数和边界处理方式有影响先用默认。后面的HMIN和HMAX是光滑长度的最小和最大缩放系数一般取粒子基本间距的一半到两倍能覆盖碰撞过程中密度变化引起的邻居粒子重构。*CONTROL_SPH里的字段我习惯先全部保持默认。等模型跑通、确认粒子分布正常之后再去研究NCBS、FORM这些参数。一上来就调这些往往会把正常模型调到粒子乱飞。另外要注意SPH与SPH之间不需要接触定义。粒子之间的相互作用是依靠邻居搜索和核函数完成的弹丸粒子进入靶板粒子的支持域内自然就会产生力。如果你在模型里看到有人加*CONTACT_NODES_TO_SURFACE那通常是SPH粒子与FEM壳或实体单元之间需要耦合的情况不是SPH-SPH的常规做法。3.3 边界条件与初始速度靶板外边界需要约束否则整块靶板会在冲击反力作用下飞出去。在LS-PrePost里选择靶板外边界一圈粒子建成一个*SET_NODE_LIST然后约束这组粒子1到6自由度。这个操作依赖具体节点编号没法在通用模板里写死但一定不要省。初始速度用*INITIAL_VELOCITY_GENERATIONSTYP2表示按PART ID选择目标PART ID1是弹丸。速度方向为Z负方向大小0.3cm/μs。注意这里填的是速度不是加速度OMEGA保持0。4. 从提交计算到排查报错粒子飞散、拉伸不稳定和“看起来没穿”4.1 命令行提交与资源占用K文件准备好之后可以用SMP版求解器提交ls-dyna isph_impact.k ncpu8 memory500m39000个粒子对这种求解器来说是很小的规模500MB内存完全够。如果后续把粒子间距加密到0.25mm粒子数会超过30万内存建议加到2GB以上ncpu可以放到16。SPH计算的核心瓶颈在邻居搜索和光滑长度更新。网格模型一个单元最多接触相邻单元但SPH粒子要和半径内的所有邻居粒子交互。粒子数上去之后并行效率下降得比FEM明显所以不要指望SPH和FEM一样能线性扩展。4.2 最常见的四种异常现象我把自己跑SPH排查过的坑整理成了一张表可以当排查手册用现象最可能的原因排查顺序粒子还没接触靶板就乱飞单位制或材料参数错误压力/密度差数量级先查密度、声速、JC参数弹丸穿透靶板靶板无反应粒子质量未正确赋值或PID/MID对应错误检查*ELEMENT_SPH的PID查看单个粒子质量临界时间步骤降计算极慢粒子间距过小或某处粒子被压缩到极高密度观察最密区域检查EOS参数粒子聚集形成非物理裂纹拉伸不稳定或光滑长度设置不合理检查HMIN/HMAX适当增加背压这里重点说单位制。很多人从网上复制了一段别的模型的材料参数直接贴进来密度可能还是kg/m³压力还是Pa算出来结果自然离谱。SPH对初始质量特别敏感。一个铝粒子的质量是密度乘以初始体积即2.7×0.05³0.0003375g。如果你在预处理里看到单个粒子质量差了几个数量级先回头核对单位制不要急着调算法参数。4.3 怎么判断“穿孔”是否真的发生SPH仿真的穿孔不是看单元删除而是看靶板粒子是否被推开并形成贯通通道。在LS-PrePost里打开d3plot设定Z方向视图逐帧播放弹丸撞击靶板正面产生高压冲击波粒子向四周飞散靶板背面出现隆起和崩落弹丸残体穿过背表面碎片云向后喷射靶板上留下一个近似圆形的孔洞如果靶板背面粒子只有微小位移、没有形成通道那说明弹丸速度低于弹道极限没有穿孔。3km/s对3mm铝板来说肯定穿孔但如果速度降到某个阈值以下就会变成“嵌在靶板里”。5. 后处理里怎么量化穿孔d3plot切片、能量曲线和孔径统计5.1 d3plot输出设置与切片模板里*DATABASE_BINARY_D3PLOT的DT0.05意思每0.05μs输出一帧。3km/s弹丸穿过3mm靶板只需要大约1μs如果DT设成1μs你只会看到撞击前和撞击后的两个状态根本看不见穿孔过程。对这类模型输出间隔控制在临界时间步的1到2倍比较合适也就是0.05μs到0.1μs。在LS-PrePost里打开d3plot后可以用平面裁剪功能切掉靶板的一半这样能清楚看到粒子通道内部的密度分布和碎片形态。SPH粒子默认显示成点可以把粒子显示尺寸调大一点但记住这只是显示效果不代表真实粒子大小。5.2 孔径的统计口径孔径是超高速碰撞最常用的输出指标但不同论文的统计口径差别很大。我建议至少统计两个值瞬时孔径撞击后某个时刻比如10μs靶板中心区域的开放孔径反映动态过程。残余孔径等碎片云基本散去、孔洞形态稳定后的孔径通常比瞬时孔径小因为靶板回弹会让孔洞部分收缩。量孔径时可以提取靶板原始中面附近一层粒子的坐标剔除飞散的碎片粒子对剩余粒子做圆拟合取拟合圆直径。如果粒子飞散严重就直接取孔的凸包直径。实测下来同一模型用不同统计口径孔径差异可能有10%到15%所以对比文献时一定要看好对方用的是哪种定义。5.3 用能量曲线验证模型合理性后处理不要只看动画能量曲线才是验证模型是否正常的硬指标。初始状态下系统总能量等于弹丸动能。撞击开始后弹丸动能减小一部分转化为靶板内能一部分转化为碎片动能。理想情况下总能量基本守恒。如果总能量曲线明显上涨通常意味着时间步控制不当、粒子质量或接触计算有问题。如果总能量大幅下降可能原因是粒子飞出了计算域或者能量被数值耗散。总能量漂移在5%以内算正常超过这个数就需要回去检查模型了。6. 复现时的参数敏感度与个人经验6.1 粒子间距的收敛性测试把粒子间距从0.1cm扫到0.05cm再到0.025cm你会看到结果逐渐收敛孔径和弹丸剩余速度变化幅度逐级减小。如果0.05cm和0.025cm的结果差距还很大说明0.05cm分辨率不够必须加密。这个收敛性测试在正式发论文、出报告之前一定要做。很多审稿人看SPH论文第一件事就是问你的结果对粒子分辨率是否敏感对工程应用来说一般找到计算结果变化小于5%的间距就够用了没必要一味追求更小粒子。6.2 材料参数对结果的支配程度我做了一组经验对比后发现穿孔孔径主要由EOS参数和密度决定Johnson-Cook强度参数对碎片形态、崩落区大小影响更明显对最终孔径影响相对弱。这符合超高速碰撞的物理规律——高压阶段材料行为由流体动力学控制强度模型只在低压膨胀阶段起作用。所以你调试模型时如果孔径差得远优先查EOS的C和S1如果碎片形态不对再查Johnson-Cook的A、B和n。不要一上来同时改七八个参数那只会让你永远找不到原因。6.3 我自己保留的几个习惯跑这类模型我始终保留三个例行检查一是跑正式工况前先拿一个单粒子或薄靶小模型验证材料参数单位制没搞错二是每次改完参数先跑2μs短时程打开d3plot看粒子行为是否合理再跑完整20μs三是后处理必看总能量曲线确认没有数值异常。SPH的优势在大变形和碎裂问题里确实无出其右但它对参数、单位制和数值控制的敏感度也比FEM高不少。按照上面的流程走一遍不敢说马上变成SPH高手但至少能少踩一半的坑。剩下的一半等你看到自己的粒子动画里出现漂亮碎片云的时候自然就有动力继续钻了。 LS-DYNA SPH超高速碰撞仿真：3km/s铝球撞击铝板完整案例
LS-DYNA这个名字很多朋友是从安全气囊展开、降落伞充气这类强非线性问题开始接触的。同一个求解器里还有一个专门对付超高速碰撞的算法家族SPH就是其中最常用的一种。今天这篇不聊理论推导直接用一块40mm见方、3mm厚的铝靶板和一枚4mm铝球在LS-DYNA里把3km/s垂直撞击穿孔的SPH仿真完整跑一遍从算法原理讲到K文件每个关键卡片最后给出一套可以改参数直接用的模板。适合已经会跑一般显式动力学模型、但对SPH还停留在“听说过”阶段的读者。先说清楚一个基本判断如果你只是做常规低速冲击比如50m/s的落锤或汽车碰撞传统Lagrange网格完全够用。但一旦弹丸速度到了2km/s以上材料行为、失效模式、靶板响应都和低速工况完全不同。这时候SPH不是炫技而是绕开网格畸变的成熟工程方案。1. 超高速工况下网格为什么“扛不住”SPH的适用边界1.1 传统Lagrange网格在3km/s下的三个死穴低速冲击里结构变形是主要的能量吸收方式单元虽然会畸变但通过自适应网格或侵蚀失效还能继续算。到了3km/s这个量级弹靶接触区的压力可以到几十GPa甚至上百GPa材料强度相对于流体静水压力已经成了次要项——换句话说金属在这个状态下表现得像流体。Lagrange网格这时候会遇到三件麻烦事单元畸变弹丸和靶板接触区的单元被极度压扁Jacobian变成负值计算直接崩掉。即使不崩精度也完全不可信。侵蚀删除导致质量丢失很多低速碰撞模型用*MAT_ADD_EROSION删单元但超高速撞击中被删除的单元里包含大量动量和能量删多了弹道极限、剩余速度全都不对。接触穿透高速下节点速度太快接触搜索步长跟不上出现节点穿透后后续结果基本报废。ALE方法能解决一部分畸变问题但超高速撞击伴随着自由表面、碎片云、裂纹扩展ALE网格在处理这些不连续现象时同样吃力。1.2 SPH把连续体拆成“一堆有质量的粒子”SPH的核心思想不复杂把连续的弹丸和靶板离散成一系列粒子每个粒子携带质量、速度、密度、压力、应力这些物理量粒子之间通过一个叫“光滑核函数”的权函数相互作用。你可以把它理解成把一整块豆腐切成很多小丁。每个小丁自带速度和密度撞击的时候小丁可以被撞飞、碎裂、重新堆积但不存在“网格扭曲”这回事。因为SPH的计算不依赖固定的单元拓扑粒子之间是“邻居关系”而不是“连接关系”粒子靠近了就算邻居飞远了就断开联系。这个特性让它天然适合碎片云、孔洞、开坑这类强不连续问题。也正因为如此SPH能捕捉到很多有限元方法很难模拟的现象着靶瞬间的冲击波传播、靶板背面崩落spall、碎片云的扩散形态。1.3 SPH不是万能药三个需要提前知道的短板SPH在超高速撞击里很好用但它有代价。拉伸不稳定材料处于拉伸应力状态时SPH粒子可能自发聚团产生非物理的断裂。这在低压区尤其明显需要通过足够的粒子分辨率、合理的核函数选择来压制。边界缺陷靠近自由表面的粒子其核函数积分区域被截断密度和压力计算会偏低。这不是bug而是SPH数学性质决定的。工程上通常靠加密表面附近粒子或使用修正光滑函数来缓解。计算成本SPH需要每一时间步搜索邻居粒子粒子间距减半粒子数量翻8倍成本涨得非常快。4万粒子还属于小规模百万粒子级别的模型就要认真考虑HPC资源了。所以SPH最适合的场景就是大变形、碎裂、自由表面演化、材料失效模式复杂。超高速碰撞正好全占。方法网格畸变碎片/裂纹计算成本适用场景Lagrange严重靠单元删除近似低低速冲击、结构碰撞ALE可接受需要界面重构中高流体结构耦合、爆炸SPH不涉及自然模拟高超高速碰撞、碎裂、穿甲2. 弹靶模型设计与单位制先把4cm见方的铝板量成K文件里的数字2.1 几何尺寸与粒子间距这个演示模型我刻意做成一个能跑得动又不会失真太多的规模靶板40mm × 40mm × 3mm铝合金弹丸直径4mm的铝球撞击速度3km/s垂直入射粒子间距0.5mm在LS-DYNA里我习惯用cm-g-μs单位制。40mm就是4.0cm3mm就是0.3cm弹丸半径0.2cm粒子间距0.05cm。靶板体积4.0×4.0×0.34.8cm³除以单个粒子占有的体积0.05³0.000125cm³靶板大约需要38400个粒子弹丸球体积约0.0335cm³折算下来粒子数大约270个。整个模型不到39000个粒子普通笔记本用多核求解器跑几十微秒的物理时间半小时内能出结果。粒子间距这个参数直接决定精度和成本。间距从0.5mm加密到0.25mm粒子数变8倍但孔径和剩余速度的收敛性会明显改善。先跑0.5mm的摸底看趋势正常再做加密。2.2 单位制换算最容易出错的0.5328cm-g-μs单位制下的关键换算关系很多第一次接触SPH的人都会在这里栽跟头物理量单位长度cm质量g时间μs速度cm/μs压力/应力Mbar1Mbar100GPa密度g/cm³几个实测参数换算过来就是3km/s 弹速 0.3cm/μs铝的声速约5328m/s 0.5328cm/μs铝合金弹性模量70GPa 0.7Mbar6061-T6铝合金的初始屈服强度324MPa 0.00324Mbar如果你把324填成0.324意味着屈服强度放大了100倍靶板会硬得像一块不可穿透的装甲如果你把0.5328填成5328冲击波声速比实际大了四个数量级临界时间步会急剧变小计算卡到怀疑人生。2.3 材料模型为什么Johnson-Cook要配EOS超高速碰撞中材料强度仍然重要但压力—密度关系同样关键。所以材料卡片要分两层第一层是强度模型我用Johnson-Cook。它把流动应力拆成初始屈服、应变强化、应变率强化和温度软化四部分。对6061-T6铝合金常用参数如下参数数值说明A0.00324 Mbar初始屈服应力B0.00114 Mbar应变强化系数n0.42应变强化指数C0.002应变率敏感系数m1.34温度软化指数TM925K熔化温度TR294K室温第二层是状态方程EOS我用Mie-Grüneisen也就是经典的冲击Hugoniot关系。它的作用是把密度和压力联系起来让材料在强烈压缩时表现出正确的流体动力学行为。铝合金参数大约取C0.5328cm/μsS11.338GAMMA02.0。纯粹的低速碰撞模型可以只给材料强度、不给EOS但超高速碰撞不给EOS压力峰值会差得离谱整个仿真没有意义。3. 核心K文件逐段拆解SPH控制、Johnson-Cook和Grüneisen状态方程3.1 K文件模板总览下面是这个模型的骨架模板。节点和SPH粒子列表需要你用LS-PrePost从六面体网格转换生成或者用脚本按规则网格生成后粘贴到*NODE和*ELEMENT_SPH段。控制卡片本身可以直接抄。*KEYWORD *TITLE SPH hypervelocity impact demo $ Units: cm, g, us, Mbar *NODE $ 节点列表由脚本或LS-PrePost生成后粘贴到此处 $ NID X Y Z TC RC ... *ELEMENT_SPH $ 粒子列表生成后粘贴到此处 $ EID PID N1 ... *PART $ PID SECID MID EOSID HGID GRAV ADPORT OPTID $ projectile part 1 1 1 1 $ target plate part 2 1 1 1 *SECTION_SPH $ SECID ELFORM CST ENT 1 1 0 0 $ HMIN HMAX IVIS START MAXV CONT 0.025 0.10 1 *MAT_JOHNSON_COOK $ MID RO G E PR DTF VP RATE 1 2.700 0.269 0.700 0.330 $ A B N C M TM TR EPSO 0.00324 0.00114 0.4200 0.0020 1.3400 925.0 294.0 1.0e-6 $ PC SPALL IT -0.003 1 *EOS_GRUNEISEN $ EOSID C S1 S2 S3 GAMAO A E0 1 0.5328 1.338 0.0 0.0 2.0 0.77 0.0 $ V0 1.0 *INITIAL_VELOCITY_GENERATION $ ID STYP OMEGA VX VY VZ IVATN ICID IVAT 1 2 0.0 0.0 0.0 -0.3 *CONTROL_TERMINATION $ ENDTIM ENDCYC DTMIN ENDNEG ENDMAS 20 *CONTROL_TIMESTEP $ DTINIT TSSFAC ISDO TSLIMT DT2MS LCTM ERODE MS1ST 0 0.900 *CONTROL_SPH $ NCBS BOXID DT ICRM FORM START MAXV CONT 0 0 0 0 0 0 0 0 *DATABASE_BINARY_D3PLOT $ DT LCDT BEAM NPLTC PSETID 0.05 *DATABASE_GLSTAT $ DT 0.10 *END节点生成我多说一句。最简单的方式是在LS-PrePost里先建一个规则六面体网格然后选中网格执行SPH粒子生成功能程序会自动把每个节点转成一个SPH粒子并写出对应的*ELEMENT_SPH卡片。手动写K文件时每个SPH元素只需要填EID、PID和参与计算的第一个节点N1后面的N2到N8留空。3.2 SPH相关卡片逐条说明*PART卡片里弹丸和靶板用同一个材料ID、同一个EOSID但PID必须分开。这样后处理里才能分别观察弹丸残体和靶板碎片的运动也方便后续给弹丸和靶板赋予不同材料。*SECTION_SPH是SPH的核心配置之一。SECID对应part里的section IDELFORM是SPH的粒子公式编号。初学者不要乱改这个编号不同公式对核函数和边界处理方式有影响先用默认。后面的HMIN和HMAX是光滑长度的最小和最大缩放系数一般取粒子基本间距的一半到两倍能覆盖碰撞过程中密度变化引起的邻居粒子重构。*CONTROL_SPH里的字段我习惯先全部保持默认。等模型跑通、确认粒子分布正常之后再去研究NCBS、FORM这些参数。一上来就调这些往往会把正常模型调到粒子乱飞。另外要注意SPH与SPH之间不需要接触定义。粒子之间的相互作用是依靠邻居搜索和核函数完成的弹丸粒子进入靶板粒子的支持域内自然就会产生力。如果你在模型里看到有人加*CONTACT_NODES_TO_SURFACE那通常是SPH粒子与FEM壳或实体单元之间需要耦合的情况不是SPH-SPH的常规做法。3.3 边界条件与初始速度靶板外边界需要约束否则整块靶板会在冲击反力作用下飞出去。在LS-PrePost里选择靶板外边界一圈粒子建成一个*SET_NODE_LIST然后约束这组粒子1到6自由度。这个操作依赖具体节点编号没法在通用模板里写死但一定不要省。初始速度用*INITIAL_VELOCITY_GENERATIONSTYP2表示按PART ID选择目标PART ID1是弹丸。速度方向为Z负方向大小0.3cm/μs。注意这里填的是速度不是加速度OMEGA保持0。4. 从提交计算到排查报错粒子飞散、拉伸不稳定和“看起来没穿”4.1 命令行提交与资源占用K文件准备好之后可以用SMP版求解器提交ls-dyna isph_impact.k ncpu8 memory500m39000个粒子对这种求解器来说是很小的规模500MB内存完全够。如果后续把粒子间距加密到0.25mm粒子数会超过30万内存建议加到2GB以上ncpu可以放到16。SPH计算的核心瓶颈在邻居搜索和光滑长度更新。网格模型一个单元最多接触相邻单元但SPH粒子要和半径内的所有邻居粒子交互。粒子数上去之后并行效率下降得比FEM明显所以不要指望SPH和FEM一样能线性扩展。4.2 最常见的四种异常现象我把自己跑SPH排查过的坑整理成了一张表可以当排查手册用现象最可能的原因排查顺序粒子还没接触靶板就乱飞单位制或材料参数错误压力/密度差数量级先查密度、声速、JC参数弹丸穿透靶板靶板无反应粒子质量未正确赋值或PID/MID对应错误检查*ELEMENT_SPH的PID查看单个粒子质量临界时间步骤降计算极慢粒子间距过小或某处粒子被压缩到极高密度观察最密区域检查EOS参数粒子聚集形成非物理裂纹拉伸不稳定或光滑长度设置不合理检查HMIN/HMAX适当增加背压这里重点说单位制。很多人从网上复制了一段别的模型的材料参数直接贴进来密度可能还是kg/m³压力还是Pa算出来结果自然离谱。SPH对初始质量特别敏感。一个铝粒子的质量是密度乘以初始体积即2.7×0.05³0.0003375g。如果你在预处理里看到单个粒子质量差了几个数量级先回头核对单位制不要急着调算法参数。4.3 怎么判断“穿孔”是否真的发生SPH仿真的穿孔不是看单元删除而是看靶板粒子是否被推开并形成贯通通道。在LS-PrePost里打开d3plot设定Z方向视图逐帧播放弹丸撞击靶板正面产生高压冲击波粒子向四周飞散靶板背面出现隆起和崩落弹丸残体穿过背表面碎片云向后喷射靶板上留下一个近似圆形的孔洞如果靶板背面粒子只有微小位移、没有形成通道那说明弹丸速度低于弹道极限没有穿孔。3km/s对3mm铝板来说肯定穿孔但如果速度降到某个阈值以下就会变成“嵌在靶板里”。5. 后处理里怎么量化穿孔d3plot切片、能量曲线和孔径统计5.1 d3plot输出设置与切片模板里*DATABASE_BINARY_D3PLOT的DT0.05意思每0.05μs输出一帧。3km/s弹丸穿过3mm靶板只需要大约1μs如果DT设成1μs你只会看到撞击前和撞击后的两个状态根本看不见穿孔过程。对这类模型输出间隔控制在临界时间步的1到2倍比较合适也就是0.05μs到0.1μs。在LS-PrePost里打开d3plot后可以用平面裁剪功能切掉靶板的一半这样能清楚看到粒子通道内部的密度分布和碎片形态。SPH粒子默认显示成点可以把粒子显示尺寸调大一点但记住这只是显示效果不代表真实粒子大小。5.2 孔径的统计口径孔径是超高速碰撞最常用的输出指标但不同论文的统计口径差别很大。我建议至少统计两个值瞬时孔径撞击后某个时刻比如10μs靶板中心区域的开放孔径反映动态过程。残余孔径等碎片云基本散去、孔洞形态稳定后的孔径通常比瞬时孔径小因为靶板回弹会让孔洞部分收缩。量孔径时可以提取靶板原始中面附近一层粒子的坐标剔除飞散的碎片粒子对剩余粒子做圆拟合取拟合圆直径。如果粒子飞散严重就直接取孔的凸包直径。实测下来同一模型用不同统计口径孔径差异可能有10%到15%所以对比文献时一定要看好对方用的是哪种定义。5.3 用能量曲线验证模型合理性后处理不要只看动画能量曲线才是验证模型是否正常的硬指标。初始状态下系统总能量等于弹丸动能。撞击开始后弹丸动能减小一部分转化为靶板内能一部分转化为碎片动能。理想情况下总能量基本守恒。如果总能量曲线明显上涨通常意味着时间步控制不当、粒子质量或接触计算有问题。如果总能量大幅下降可能原因是粒子飞出了计算域或者能量被数值耗散。总能量漂移在5%以内算正常超过这个数就需要回去检查模型了。6. 复现时的参数敏感度与个人经验6.1 粒子间距的收敛性测试把粒子间距从0.1cm扫到0.05cm再到0.025cm你会看到结果逐渐收敛孔径和弹丸剩余速度变化幅度逐级减小。如果0.05cm和0.025cm的结果差距还很大说明0.05cm分辨率不够必须加密。这个收敛性测试在正式发论文、出报告之前一定要做。很多审稿人看SPH论文第一件事就是问你的结果对粒子分辨率是否敏感对工程应用来说一般找到计算结果变化小于5%的间距就够用了没必要一味追求更小粒子。6.2 材料参数对结果的支配程度我做了一组经验对比后发现穿孔孔径主要由EOS参数和密度决定Johnson-Cook强度参数对碎片形态、崩落区大小影响更明显对最终孔径影响相对弱。这符合超高速碰撞的物理规律——高压阶段材料行为由流体动力学控制强度模型只在低压膨胀阶段起作用。所以你调试模型时如果孔径差得远优先查EOS的C和S1如果碎片形态不对再查Johnson-Cook的A、B和n。不要一上来同时改七八个参数那只会让你永远找不到原因。6.3 我自己保留的几个习惯跑这类模型我始终保留三个例行检查一是跑正式工况前先拿一个单粒子或薄靶小模型验证材料参数单位制没搞错二是每次改完参数先跑2μs短时程打开d3plot看粒子行为是否合理再跑完整20μs三是后处理必看总能量曲线确认没有数值异常。SPH的优势在大变形和碎裂问题里确实无出其右但它对参数、单位制和数值控制的敏感度也比FEM高不少。按照上面的流程走一遍不敢说马上变成SPH高手但至少能少踩一半的坑。剩下的一半等你看到自己的粒子动画里出现漂亮碎片云的时候自然就有动力继续钻了。