仿真笔记——一文详解ABAQUS中冲击动力学问题的求解方法

冲击载荷随时间迅速变化。当物体的局部位置受到冲击时，所产生的扰动会逐渐传到未扰动的区域去，这种现象称为应力波的传播。当载荷作用时间短、变化快，且受力物体在加载方向的尺寸又足够大时，这种应力波的传播就显得特别重要。

研究动力学问题最终将简化为求解动力学平衡方程式：节点质量矩阵M乘以节点加速度 等于节点的合力（所施加的外力P与单元内力I之间的差值）：

   

1. 冲击动力学求解方法

 

ABAQUS/CAE（前后处理）

 

ABAQUS/Standard（隐式求解器）

 

ABAQUS/Explicit（显式求解器）

   

2. 动力学显式有限元方法

ABAQUS/Explicit应用中心差分方法对运动方程进行显式的时间积分，应用一个增量步的动力学条件计算下一个增量步的动力学条件。在增量步开始时，程序求解动力学平衡方程式：节点质量矩阵M乘以节点加速度等于节点的合力（所施加的外力P与单元内力I之间的差值），即式(2-1)：

                                    在当前增量步开始时（t时刻），计算加速度：

                (2-2)        

由于显式算法总是采用一个对角的或者集中的质量矩阵，所以求解加速度并不复杂，不必同时求解联立方程。采用节点的加速度完全取决于节点质量和作用在节点上的合力，使得节点计算的成本非常低。

对加速度在时间上进行积分采用中心差分法，在计算速度的变化时假定加速度为常数。应用这个速度的变化值加上前一个增量步中点的速度来确定当前增量步中点的速度：

               (2-3)

速度对时间的积分加上在增量步开始时的位移可以确定增量步结束时的位移：

                       (2-4)

这样，在增量步开始时提供了满足动力学平衡条件的加速度。得到了加速度，在时间上“显式地”前推速度和位移。所谓“显式”是指在增量步结束时的状态仅依赖于该增量步开始时的位移、速度和加速度。这种方法精确地积分常值的加速度。为了使该方法得到精确的结果，时间增量必须相非常小，这样在增量步中加速度几乎为常数。由于时间增量步必须很小，所以一个典型的分析需要成千上万个增量步。因为不必同时求解联立方程组，所以每一个增量步的计算成本很低。大部分计算成本消耗在单元计算上，以此确定作用在节点上的单元内力。单元的计算包括确定单元应变和应用材料本构关系（单元刚度）确定单元应力，从而进一步计算内力。

下面给出显式动力学方法的总结：

(1) 节点计算

1) 动力学平衡方程

                            (2-5)   

2) 对时间显式积分

           (2-6)  

                  (2-7)

(2) 单元计算

1) 根据应变速率，计算单元应变增量

2) 根据本构关系计算应力

             (2-8)

3) 集成节点内力I

(3) 设置时间t为，返回到步骤(1)。

                       (2-9)

在完全Newton迭代求解方法的基础上，ABAQUS/Standard使用自动增量步。在时刻 增量步结束时，Newton方法寻求满足动力学平衡方程，并计算出同一时刻的位移。由于隐式算法是无条件稳定的，所以时间增量 比应用显式方法的时间增量相对大些。对于非线性问题，每一个典型的增量步需要经过几次迭代才能获得满足给定容许误差的解答。每次Newton迭代都会得到对于位移增量 的修正值 。每次迭代需要求解的一组瞬时方程为

                      (2-10)

对于较大的模型，这是一个非常大的计算过程。有效刚度矩阵 是关于本次迭代的切向刚度矩阵和质量矩阵的线性组合。直到一些变量满足了给定的容许误差才结束迭代，如位移修正值等。

稳定性限制了ABAQUS/Explicit求解器所能采用的最大时间步长，这是ABAQUS/Explicit 进行计算的一个重要因素。应用显式方法，基于在增量步开始时刻 的模型状态，通过时间增量 前推到当前时刻的模型状态。使得状态能够前推并仍能保持对问题的精确描述的时间非常短。如果时间增量大于这个最大的时间步长，则时间增量已经超出了稳定性限制（stability limit）。超过稳定性限制的后果是数值不稳定，可能导致计算不收敛。由于一般不可能精确地确定稳定性限制，因而采用保守的估计值。为了提高计算效率，ABAQUS/Explicit选择时间增量，使其尽可能地接近而且又不超过稳定性限制。

在系统中，以最高频率（ ）的形式定义稳定性限制。无阻尼的稳定极限由下式定义：

                                   (2-11)有阻尼的稳定极限由下面的表达式定义：

 (4-14)                  (2-12)式中， 是最高频率模态的临界阻尼部分。由此可以看出，阻尼通常是减小稳定性限制的。

基于逐个单元的估算，稳定极限可以通过单元长度 和材料的波速 重新定义：

                             (2-13)

                           (2-14)

其中， 是杨氏模型； 是材料密度。由上式可以看出，材料的刚度越大，波速越大，导致越小的稳定极限；密度越大，波速越小，导致较大的稳定极限。

稳定极限是当膨胀波通过由单元特征长度定义的距离时需要的时间。如果知道最小的单元尺寸和材料的波速，就能够估算稳定极限。例如，如果最小单元尺寸是5mm，膨胀波速是5000m/s，那么稳定的时间增量就在 的数量级上。

   

3. 接触