📚 课程目录第2章显式动力学基础

2. 显式动力学基础

在汽车碰撞仿真中，LS-DYNA 求解的是结构在极短时间（通常 50～150 ms）内经历大变形、大转动和复杂接触的高度非线性动力学问题。显式动力学之所以能胜任这类问题，核心在于它不组装整体刚度矩阵、不做迭代求解，而是用中心差分法在时间域上逐步推进。理解这套时间积分机制，是判断模型能否算得稳、算得准的前提。

中心差分法：显式推进的数学内核

中心差分法把加速度用位移的二阶中心差分近似。设当前时刻为 t_n，时间步长为 Δt，则节点加速度可写为：

a_n = (u_{n+1} - 2u_n + u_{n-1}) / Δt²

代入运动方程 M·a + C·v + F_int = F_ext（M 为质量矩阵，C 为阻尼矩阵，F_int 为内力，F_ext 为外力），并利用显式质量矩阵通常为对角阵这一特点，可解出下一时刻位移：

u_{n+1} = f(u_n, u_{n-1}, F_int, F_ext, M, Δt)

由于 M 是对角阵，求逆只需逐分量相除，无需矩阵求逆和迭代。这就是显式方法单步计算量小、天然适合并行和 GPU 加速的根本原因。整个求解过程是“已知前两步、推出下一步”的递推：先由位移差得到速度与加速度，再由单元应力更新内力，最后推进位移。每一步都不需要收敛判断，因此不存在隐式方法中常见的迭代不收敛问题，但代价是必须严格控制时间步长。

条件稳定性与临界时间步长

中心差分法是条件稳定的：只有当时间步长小于某个临界值时，数值解才不会发散。这个临界值由结构中最小单元的特征尺寸和材料波速决定。对一维杆单元，临界时间步长近似为：

Δt_crit = L_e / c

其中 L_e 是单元特征长度，c 是材料中的声速。对线弹性材料，c = √(E/ρ)，E 为弹性模量，ρ 为密度。LS-DYNA 实际使用的稳定步长还要乘以一个小于 1 的安全系数（默认约 0.9），即 Δt = 0.9·Δt_crit。这意味着：网格越小、材料越硬、密度越低，允许的时间步长就越小。例如钢材 E≈210 GPa、ρ≈7.85×10⁻⁹ t/mm³，声速约 5170 mm/ms；若最小单元尺寸为 5 mm，则 Δt_crit≈0.97 μs，实际步长约 0.87 μs。一个 100 ms 的碰撞过程需要约 11.5 万步，这正是显式碰撞仿真计算量巨大的原因。

时间步长控制：从被动接受到主动管理

在 LS-DYNA 中，时间步长并非完全由求解器自动决定而不可干预。常用的控制手段包括：

质量缩放（Mass Scaling）：通过人为增大密度来提升波速分母中的 ρ，从而增大 Δt_crit。LS-DYNA 提供 *CONTROL_TIMESTEP 中的 DT2MS 参数，可对小于目标步长的单元自动缩放质量。但质量缩放会引入附加惯性，必须监控质量增加百分比（通常控制在 1%～5% 以内），否则会改变碰撞响应。
网格尺寸控制：在碰撞仿真中，最小单元往往出现在圆角、倒角或薄壁特征处。一个 0.5 mm 的单元会把全局步长拖到 0.1 μs 量级。因此建模阶段就要有意识地控制最小特征尺寸，必要时用局部细化代替全局细化。
子循环（Subcycling）：LS-DYNA 默认启用子循环，让不同部件按各自稳定步长推进，再在同步点交换信息。这能缓解局部小单元对全局步长的拖累，但同步开销和接触搜索频率仍需关注。
初始时间步长设置：通过 *CONTROL_TIMESTEP 的 DTINIT 显式指定初始步长，或调整 TSSFAC 安全系数（默认 0.9，可适当降低至 0.8 以增强稳定性）。

实际工程中，建议在求解前用 LS-DYNA 的 *CONTROL_TIMESTEP 配合 *DATABASE_BINARY_D3PLOT 输出步长信息，并在 message 文件中检查 “smallest timestep” 和 “mass scaling” 报告。若发现某几个单元反复主导最小步长，应回到网格阶段处理，而不是一味依赖质量缩放。

沙漏能与能量平衡：结果可信度的守门人

显式动力学常用单点积分单元以降低计算量，但单点积分会带来沙漏模式（Hourglassing）——单元在无能量消耗的情况下产生非物理的锯齿状变形。LS-DYNA 通过沙漏控制（如 *HOURGLASS 关键字，类型 1～6）施加人工刚度或粘性阻尼来抑制。但沙漏控制本身会引入附加能量，即沙漏能。

判断沙漏能是否可接受，通常看两个指标：

沙漏能绝对值：一般要求沙漏能不超过内能的 5%。若超过 10%，结果通常不可信。
能量比曲线：在 glstat 文件中观察 “hourglass energy / internal energy” 随时间的变化。若在碰撞后期该比值持续上升，说明沙漏控制不足或网格质量差。

能量平衡是另一个关键检查点。显式求解中，总能量应近似守恒（考虑外力做功和阻尼耗散）。若总能量出现异常跳变，往往意味着：接触穿透导致滑移能突变、质量缩放过大引入虚假动能、或单元负体积导致内能计算崩溃。在 LS-DYNA 中，通过 *DATABASE_GLSTAT 和 *DATABASE_MATSUM 输出能量曲线，是每次碰撞仿真后必须执行的检查步骤。

工程判断示例：某轿车正面 100% 重叠刚性墙碰撞模型，初始动能约 120 kJ。求解完成后检查 glstat：内能峰值约 95 kJ，沙漏能峰值约 3.2 kJ，沙漏能/内能≈3.4%，滑移能约 8 kJ，总能量波动小于 2%。此时可认为能量平衡正常、沙漏控制有效。若沙漏能/内能升至 12%，则应优先检查前纵梁和保险杠区域的网格长宽比与沙漏控制类型，而非直接接受结果。

建议在每次提交计算前，用 *CONTROL_ENERGY 开启 HGEN 和 RWEN 选项，确保沙漏能和滑移能被正确统计。计算完成后，先看 message 文件中的时间步长摘要，再看 glstat 的能量曲线，最后才看变形和加速度结果——这个顺序能帮你快速排除数值问题，避免在错误结果上做物理分析。

总结来说，中心差分法决定了显式求解的递推结构，条件稳定性决定了时间步长的上限，时间步长控制决定了计算成本与精度的平衡，而沙漏能与能量平衡则是判断结果是否可信的最终依据。这四者在汽车碰撞仿真中环环相扣：网格质量影响步长，步长影响计算规模，沙漏控制影响能量分布，能量分布又反过来验证模型是否合理。掌握这条逻辑链，才能在 LS-DYNA 中既算得动、又算得对。