隐式与显式动力学的区别-弹性动力学有限元基本解法

结构系统的通用运动学方程为：

 （1）

求解该动力学振动响应主要有三类方法：（1）时域法（2）频域法（3）响应谱法

时域法又可分为：（1）直接积分法，（2）模态叠加法。直接积分法又可分为中心差分法（显式），Wilson （隐式）法以及Newmark（隐式）法等。

本文介绍中心差分法（显式）与Newmark（隐式）法。

1 中心差分法（显式）

假定0，，，…，时刻的节点位移，速度与加速度均为已知，现求解时刻的 结构响应。中心差分法对加速度，速度的导数采用中心差分代替，即为：

将（2）式代入（1）式后整理得到

 （3）

式（3）中

分别称为有效质量矩阵，有效载荷矢量。，，，为结构载荷，质量，阻尼，刚度矩阵。

求解线性方程组（3），即可获得时刻的节点位移向量，将代回几何方程与物理方程，可得时刻的单元应力和应变。

中心差分法在求解瞬时的位移时，只需时刻以前的状态变量和，然后计算出有效质量矩阵，有效载荷矢量，即可求出，故称此解法为显式算法。

中心差分法，在开始计算时，需要仔细处理。=0时，要计算，需要知道的值。因此应该有一个起始技术，因而该算法不是自动起步的。由于，，是已知的，由=0时的（2）式可知：

中心差分法中时间步长的选择涉及两个方面的约束：数值算法的稳定性和计算时间。中心差分法的实质是用差分代替微分，并且对位移和加速度的导数采用线性外插，这了的取值不可过大，否则结果可能失真过大。

可以证明：中心差分法是条件稳定的。即当时间步长必须小于由该问题求解方程性质所决定的一个时间步长的临界值。LS-DYNA中，采用“变时间步长法”，即每一时刻的步长由当前结构的稳定性条件来控制。具体算法为：计算每一个单元的极限时间步长，=1,2…，取为下一个时刻的时间步长。各种单元的计算方法如下。

（1）1D杆，梁单元

其中为时间步长因子，系统默认为0.9。为杆，梁单元的长度。材料声速。

（2）2D板，壳单元

其中为时间步长因子，系统默认为0.9。为 壳单元的最小单元边长度。为材料的声速。

（3）3D单元

其中

和为无量纲常数，默认=1.5，=0.06.

 为单元等效长度， 为单元体积，为单元最大侧面积。

为材料声速。

时间步因子可由用户设置，减小相当于减少时间步长。设置时间步长因子的关键字为*CONTROL_TIMESTEP，控制参数为TSSFAC。

另外，质量缩放可以人为控制时间步长。即调整单元密度来改变时间步长。以壳单元为例说明质量缩放改变时间步长。

 （假如=1）

由得到

LS-DYNA中有2种质量缩放方案，修改*CONTROL_TIMESTEP中的参数。

方案1：DT2MS为正的时间步

通过调整单元密度，使所有单元时间步相同，只用于惯性效应不重要的情况。

方案2：DT2MS为负的时间步

质量缩放只用于小于指定时间步长的单元。惯性效应应该通过变形体的动能与内能的比例进行衡量（一般应该小于10%等）。

2 Newmark法（隐式）

Newmark假定在时间间隔内，加速度线性变化，即采用如下的加速度，速度公式：

式中，为按积分的精度和稳定性要求可以调整的参数。

根据（4）式可给出和用，，表示的表达式，代入（1）式中整理得到

 （5）

其中

 称之为有效刚度矩阵和有效载荷矢量。由上式可以看出求解当前，需要用到当前时刻的，因此该算法为隐式算法。

当载荷历史全部已知时，为已知量，求解需要迭代实现。

可以证明，当参数，时，Newmark法是无条件稳定的，即的大小不影响数值稳定性。此时时间步长的选择主要根据解得精度确定。一般，Newmark法可以比中心差分法的时间步长大得多。

3．结论

比较两种算法，显式中心差分法非常适合研究波的传播问题，如碰撞、高速冲击、爆炸等。分析式（3）发现，显式中心差分法的M与C矩阵是对角阵，如给定某些有限元节点以初始扰动，在经过一个时间步长后，和它相关的节点进入运动，即U中这些节点对应的分量成为非零量，此特点正好和波的传播特点相一致。另一方面，研究波传播的过程需要微小的时间步长，这也正是中心差分法的特点。

而Newmark法更加适合于计算低频占主导的动力问题，从计算精度考虑，允许采用较大的时间步长以节省计算时间，同时较大的时间步长还可以过滤掉高阶不精确特征值对系统响应的影响。隐式方法要转置刚度矩阵，增量迭代，通过一系列线性逼近（Newton-Raphson）来求解。正因为隐式算法要对刚度矩阵求逆，所以计算时要求整体刚度矩阵不能奇异，对于一些接触高度非线性问题，有时无法保证收敛。