3. LS-DYNA显式动力学基础与求解控制

薄壁结构溃缩仿真本质上是一个高度非线性的瞬态动力学问题：压溃过程中伴随大变形、塑性铰形成、接触状态频繁切换以及局部屈曲，隐式方法在每一步都要组装并求解除刚度矩阵，遇到接触突变和单元畸变时极易不收敛。LS-DYNA 采用显式中心差分时间积分，用"小步快跑"的方式绕开迭代求解，这正是它适合吸能盒、前纵梁、薄壁圆管这类溃缩问题的根本原因。本节围绕显式积分原理、稳定时间步长、质量缩放、沙漏控制和能量平衡检查五个控制点展开，并给出可直接落到关键字文件中的设置方式。

显式中心差分积分原理

显式方法不求解位移方程组，而是由当前时刻的加速度直接推进速度与位移。以节点自由度为例，半离散后的运动方程为 M·a = F_ext − F_int，其中 M 为集中质量矩阵（对角阵），因此加速度可逐节点显式求出：a_n = M⁻¹(F_ext − F_int)。随后按中心差分更新：

v_(n+1/2) = v_(n-1/2) + a_n · Δt_n
d_(n+1)   = d_n + v_(n+1/2) · Δt_(n+1/2)

速度定义在整步之间的半时刻，位移定义在整时刻，这种交错格式使求解完全不需要矩阵求逆和迭代，单步计算量正比于单元数。代价是条件稳定：时间步长必须小于临界值，否则数值解会指数发散。因此显式求解的核心不是"能不能算"，而是"步长是否受控"。

稳定时间步长与临界步长估算

临界步长由结构中最小单元的特征尺寸和材料波速决定。对壳单元常用估算式：

Δt_crit = L_e / c ,  c = sqrt(E / ρ)

其中 L_e 为单元特征长度，c 为材料声速。以典型车身用钢为例，E = 210 GPa，ρ = 7.85e-9 t/mm³，则 c ≈ 5170 mm/ms。若最小壳单元尺寸为 3 mm，Δt_crit ≈ 3/5170 ≈ 5.8e-4 ms。LS-DYNA 实际使用的时间步还会乘以安全系数（默认 0.9 左右），即约 5.2e-4 ms。若总物理时间为 30 ms，则需约 5.8 万步。可见单元越小、步数越多，这是显式仿真成本的主要来源。

控制手段有三条：一是避免在非关键区域出现极小单元，网格划分时把最小尺寸控制在 3～5 mm；二是用 *CONTROL_TIMESTEP 设置 DTINIT 或 TSSFAC 调整安全系数；三是当局部小单元无法避免时，启用质量缩放。需要强调的是，稳定步长由全局最小单元决定，一个 0.5 mm 的畸变单元就能把整模型步长拖低一个数量级。

质量缩放：用可控的质量代价换步长

质量缩放通过人为增大密度来降低波速、放大临界步长。LS-DYNA 提供两种方式：*CONTROL_TIMESTEP 中的 DT2MS 为固定步长缩放，DT2MS 取负值时表示只对步长小于 |DT2MS| 的单元缩放，正值则对所有单元统一缩放。工程上更常用负值方式，避免对本来满足步长要求的单元无谓增重。

*CONTROL_TIMESTEP
$  DTINIT    TSSFAC    ISDO     TSLIMT     DT2MS      LCTM
   0.0       0.9       0        0.0      -5.0e-4      0

上例表示把全局步长目标设为 5.0e-4 ms，仅对不满足该步长的单元进行缩放。缩放后必须检查附加质量占比：在 d3hsp 文件中查看 "added mass" 与总质量之比，一般要求小于 1%～5%。对吸能结构而言，附加质量会直接改变动能与吸能量，若占比过大，压溃力和吸能量都会被污染。因此质量缩放应作为"局部补救"，而不是全局提速手段。

沙漏控制

显式方法普遍采用单点积分壳单元以节省计算量，但单点积分会产生零能变形模式，即沙漏。沙漏模式在视觉上表现为网格出现交替的锯齿状畸变，能量上表现为沙漏能异常增长。控制方式有两类：一是粘性沙漏控制，通过 *HOURGLASS 设置 IHQ 与系数 QM，常用 IHQ=4（类型 4）配合 QM=0.03～0.05；二是刚度沙漏控制，适用于需要更硬响应的场合，但可能引入虚假刚度。

*HOURGLASS
$   IHQ       QM
     4      0.05

判断沙漏是否受控，关键看沙漏能占内能的比例。经验判据是沙漏能应小于内能的 5%，严格场合要求小于 2%。若超标，优先检查网格质量与加载方式，其次再调整沙漏控制系数，切忌一味加大 QM，否则会人为提高结构刚度、压低压溃位移。

能量平衡检查

能量平衡是判断显式仿真是否可信的总闸门。总能量应满足：

E_total = E_kinetic + E_internal + E_hourglass + E_sliding + E_damping − E_external

在 glstat 文件中逐时刻输出各项能量。健康溃缩仿真的典型特征是：初始动能或外力做功逐步转化为内能（塑性变形与吸能），动能先升后降，沙漏能与滑移能保持小量，总能量在无外力输入时基本守恒。若总能量出现明显漂移，通常意味着接触穿透、单元负体积或时间步失控。工程上常设定三条红线：沙漏能/内能 < 5%，滑移能/内能 < 5%，总能量偏差 < 5%。

溃缩仿真中的参数设置示例

把上述控制点落到一个薄壁方管轴向压溃算例上。方管截面 60 mm × 60 mm，壁厚 1.5 mm，长度 200 mm，材料为 DP780，密度 7.85e-9 t/mm³，弹性模量 210 GPa，屈服强度约 450 MPa。网格最小尺寸控制在 3 mm，采用 Belytschko-Tsay 壳单元。关键字设置如下：

*CONTROL_TERMINATION
$   ENDTIM
    30.0

*CONTROL_TIMESTEP
$  DTINIT    TSSFAC    ISDO     TSLIMT     DT2MS      LCTM
   0.0       0.9       0        0.0      -5.0e-4      0

*CONTROL_ENERGY
$  HGEN      RWEN     SLNTEN    RYLEN
    2         2         2         2

*HOURGLASS
$   IHQ       QM
     4      0.05

*DATABASE_GLSTAT
$   DT
    0.1

*DATABASE_MATSUM
$   DT
    0.1

其中 *CONTROL_ENERGY 的 HGEN=2 表示计算沙漏能，RWEN=2 表示计算刚性墙能量，SLNTEN=2 表示计算滑移能，这些开关是后续能量核查的前提，必须打开。加载方式可用刚性墙以恒定速度 5 m/s 压溃，或给方管一端施加 10 m/s 初速度撞击刚性墙，两种方式都能在 30 ms 内完成充分压溃。

结果判读与常见问题

仿真结束后，先看 glstat 的能量曲线：内能应单调上升并趋于平台，动能先升后降，沙漏能与滑移能保持低位。再看 matsum 确认附加质量占比。最后看压溃力—位移曲线，典型薄壁方管应呈现周期性折叠峰谷，峰值力不宜超过平均压溃力的 1.5～2 倍，否则说明初始触发或网格过渡不合理。

常见问题与对应处理：若沙漏能持续增长，检查是否使用了单点积分单元且网格过于规则，可改用 IHQ=4 并适当提高 QM，或局部改用全积分单元；若滑移能异常偏高，检查接触厚度与穿透，调整 *CONTACT 的 SOFT 与 PENALTY 参数；若总能量漂移，优先排查负体积单元与时间步是否被局部小单元拖垮，必要时对畸变区域做网格重划或删除失效单元。质量缩放带来的附加质量若超过 5%，应回到网格阶段减小最小单元尺寸，而不是继续加大 DT2MS。