📚 课程目录第3章 LS-DYNA显式动力学基础

3. LS-DYNA显式动力学基础

跌落仿真之所以几乎全部采用显式动力学求解，根本原因在于跌落过程本身：包装件从约0.8 m高度自由落体，触地瞬间速度可达约4 m/s，冲击持续时间通常只有几毫秒到十几毫秒，期间接触力、加速度和应力在极短时间内剧烈变化。隐式方法需要反复组装和求解刚度矩阵并迭代收敛，在这种强非线性、短时、高瞬态问题中代价极高；显式方法则通过中心差分直接推进，无需迭代求解联立方程，天然适合这类问题。本节围绕显式时间积分、条件稳定、时间步长控制、沙漏模式、能量平衡与跌落适用性展开。

显式时间积分：中心差分推进

LS-DYNA的显式求解基于运动方程 M·a = F_ext − F_int，其中M为集中质量矩阵。由于采用集中质量（lumped mass），M是对角阵，求逆只需逐节点做除法，无需矩阵分解。时间推进采用中心差分格式：

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)
a^n       = M^(-1) · (F_ext^n − F_int^n)

可见每一步只需已知当前时刻的外力与内力，即可显式算出下一时刻的位移，不存在迭代和收敛判断。内力由单元应力积分得到，应力又由应变率和应变更新，因此整个求解是“一步接一步”的推进过程。这也意味着显式方法对每一步的计算量小、可高度并行，但步数极多——跌落仿真动辄几十万到上百万步。

条件稳定与临界时间步长

显式中心差分是条件稳定的：时间步长必须小于临界值，否则数值解会指数发散。临界步长由最小单元尺寸和材料波速决定：

Δt_crit = L_e / c
c = sqrt(E / ρ)          （杆/壳单元，弹性波速）
Δt = α · Δt_crit         （α 为时间步缩放因子，默认约0.9）

其中L_e为单元特征长度，c为材料声速，E为弹性模量，ρ为密度。以钢材为例，E≈210 GPa，ρ≈7850 kg/m³，c≈5170 m/s；若最小单元尺寸为2 mm，则Δt_crit≈0.39 μs，乘以0.9后约0.35 μs。若跌落仿真物理时间为10 ms，则需约2.9万步。若模型中存在0.2 mm的细小单元，临界步长会骤降到0.035 μs，步数增加十倍，计算时间随之暴涨。

质量缩放（*CONTROL_TIMESTEP中的DT2MS）通过人为增大密度来放大临界步长，可显著减少步数，但会引入附加惯性。跌落仿真中应控制质量增加比例（建议不超过1%~2%），并在能量平衡中检查其影响。

时间步长控制实践

在LS-DYNA中，时间步长由*CONTROL_TIMESTEP卡片控制。常用做法是：先用DT2MS = 0运行，通过d3hsp文件查看最小单元和实际步长；若步数过多，再设置DT2MS为负值（如−1e-6），表示仅对步长小于该值的单元进行质量缩放，避免全局增质。同时应开启*DATABASE_BINARY_D3PLOT控制输出频率，避免输出文件过大。跌落仿真中还需注意接触搜索和穿透，步长过大会导致接触力计算失真。

沙漏模式

显式求解常用单点积分（缩减积分）单元以降低计算量，但单点积分单元存在零能模式，即沙漏模式：单元在无应变能增加的情况下发生非物理的锯齿状变形。沙漏能（hourglass energy）若不受控，会污染结果甚至使计算崩溃。控制手段包括：使用*HOURGLASS卡片设置沙漏控制类型（如类型4、5）和系数（通常0.03~0.1）；对关键部件采用全积分或选择性缩减积分；细化网格。判断标准是沙漏能应小于总能量的5%，理想情况小于1%。

能量平衡：结果可信度的核心判据

显式动力学必须检查能量平衡。总能量应满足：

E_total = E_kinetic + E_internal + E_hourglass + E_slide + E_damping + E_external

跌落仿真中，初始只有重力势能，触地后动能转化为内能（变形能）和接触滑移能。若沙漏能或滑移能异常增大，说明模型有问题。典型判据：沙漏能/总能量 < 5%；滑移能不应为负且不应超过内能的10%~15%；总能量在无外力做功时应基本守恒。能量曲线可通过glstat文件或LS-PrePost的ASCII面板查看。

跌落仿真适用性判断

显式动力学并非万能。判断是否适用，可看三点：一是过程是否短时高瞬态（跌落、碰撞、冲击、爆炸适用；准静态压溃、长期蠕变不适用）；二是是否强非线性（大变形、接触、材料失效适用）；三是是否需要精确的长期动力响应（模态、频响不适用）。跌落包装仿真恰好满足前两点，因此显式方法是标准选择。

实际建模时，包装件通常包含EPS泡沫、瓦楞纸板、塑料卡扣等，泡沫材料需用低密度泡沫本构（如*MAT_LOW_DENSITY_FOAM），并注意其应变率效应。接触方面，包装与地面、包装与产品之间需定义*CONTACT_AUTOMATIC_SURFACE_TO_SURFACE，并设置合理的摩擦系数。初始速度通过*INITIAL_VELOCITY施加，重力通过*LOAD_BODY_Z施加。求解时间一般取冲击持续时间的1.5~2倍，确保能量充分转化。

建议在提交计算前做一次“空跑”检查：用极小终止时间（如0.1 ms）运行，查看d3hsp中的最小时间步、最小单元和初始穿透信息，确认模型无异常后再正式计算。

综上，显式时间积分以中心差分实现高效推进，条件稳定要求步长受最小单元和声速约束，沙漏模式和能量平衡是判断结果可信度的关键。掌握这些基础，才能在后续跌落仿真中合理设置参数、正确解读结果。