Abaqus中显示动力学分析

ABAQUS 中的显式动力学分析 摘要 ABAQUS 的分析模块主要由两个求解器完成 ABAQUS Stand 和 ABAQUS Explicit 前者主 要用于静态分析和简单的准静态以及接触分析 而后者用于高速碰撞 复杂接触以及材料的 失效和退化仿真分析 ABAQUS Explicit 进行动力学分析最重要的是求解算法 稳定时间极 限和能量平衡问题 它们直接关系到计算的时间和计算所占用的内存以及计算结果的准确性 本文通过理论分析和实际举例验证和分析了 ABAQUS Explicit 理论在实际中的应用 关键词 ABAQUS Explicit 求解算法 稳定时间极限 能量平衡 张敬东 成焕波 合肥工业大学 机械与汽车工程学院 合肥 230009 0 前言 The Explicit Dynamic 分析方法起初专门用于对高速动力学的分析 主要原因是隐式求解方 法特别占用内存而且极其耗时 比如在爆炸分析时载荷在瞬间增大而且持续时间非常短 精 确的捕捉结构响应的应力波是极其重要的 而应力波与系统的最高阶固有频率和最小时间增 量步有直接的关系 The Explicit Dynamics 对于结构的瞬态响应和随后的震荡接触尤其适 用 而对于结构的刚度在分析的过程中发生显著变化的复杂后屈曲处理以及高度非线性的准 静态分析和材料的失效和退化 该方法也能以及其独特的算法和极小的时间增量步完成 并 且不存在隐式求解中经常遇到的收敛问题 1 求解算法 ABAQUS Explicit 在整个运算过程中采用中心差分的算法求解整体平衡方程 运用动态方程 在每个微小的增量步内计算出下一个增量步的动态状况 假定节点的质量矩阵 M 加速度 u 外力作用为 P 单元的内力 I 因此节点的平衡方程是 M u PI 根据方程可以计算出节点在当前增量步 t 时的节点加速度为 u t M 1 P I t 由于显示求解器所用的质量矩阵全部为对角矩阵 而无需求解方程组 节点的加速度完全由 质量矩阵和施加在节点上的外力决定 因此求解十分节省资源 在极短的时间内假定加速度为常值就可以运用中心差分对加速度积分计算出节点的速度方 程 u t 2 t u t 2 t 2 tu ttttt 再次运用中心差分对速度积分求解出位移方程 u t t u t t t t u t 2 t 由 方程可知模型节点的位移 速度和加速度只与初始增量步的大小有关 因此 ABAQUS Explicit 采用极小的初始时间增量和很多的增量步完成中心差分的运算 由于对于不同节点的网格类型和网格的形状 后处理时查看的结果可能是节点应力也可能是 单元应力 假定单元应变增量为 d 单元的应变率为 单元应力为 由材料力学中的应力应变关系可以推导出单元应力应变的关系 dtftt 注 上述方程中t 为初始时间增量步 在 ABAQUS Stand 中进行非线性分析时 每一个时间增量步都要求解大量的非对角矩阵平衡方 程组 及其耗时而且在涉及接触时最容易出现不收敛的现象 因而采用 ABAQUS Explicit 进 行动力学分析具有省时 节约内存和精确性较高的有点等 2 稳定时间极限 stable t ABAQUS Explicit 的算法表明时间增量步t 的大小直接关系到计算的准确性 如果该时间大 于最大的允许时间会导致数值不稳定和求解的无限大现象 因此在分析之前必须尽可能精确 地估算出稳定时间极限 要求到数量级 稳定时间极限由系统最高阶频率 max 决定 在无阻尼状态下 max 2 stable t 有阻尼状态下 1 2 2 max stable t 最高阶频率的临界阻尼系数 但如果要确定系统的最高阶频率必须对系统进行线性摄动分析提取多阶模态值 过程较为复 杂 但是系统的最高阶频率与每个单元的最高阶频率紧密相关 并且单元的最高阶频率略高 于整体的最高阶频率 为此 ABAQUS Explicit 采用更为简便的方法估算稳定极限时间 即计 算每个单元的最高阶频率对应的稳定时间极限 stable t d e C L L e 单元在各个尺寸方向上的最小值 Cd 模型材料的波速 这里只适用于弹塑性材料和线弹性材料 超弹性完全不可压缩材料的 波速趋近于无穷 由式可以得出稳定时间极限与单元的最小尺寸成正比 与材料的波速成反比 因而网格越 细 所需的稳定时间极限越小 但是不同材料的波速不同 Cd E E 材料的弹性模量 材料的密度 由此可见材料越刚硬波速越大 稳定时间极限越小 而材料的密度越大波速越小 稳定时间 极限越大 增大稳定时间极限有利于分析的顺利进行和数值稳定同时也使计算的准确性下降 因而在实 际运算中可以先估计出稳定时间极限 用较大的时间增量步先进行分析计算 在运算出结果 之后进行细化网格精确运算 在控制稳定时间极限时经常采用的方法是用细化网格 质量缩 放控制材料的波速 3 能量平衡 能量平衡作为 ABAQUS 后处理历程数据 是评估分析精确性的重要组成部分 通过对比各种能 量的比值可以判断结构相应的准确性 总的能量平衡方程 EI EV EFD EKE EW EPW ECW EMW ETota l Costant EI 内能 EV 粘性耗散能 EFD 摩擦耗散能 EKE 动能 EW 外力做功 EPW 接触罚函数做功 ECW 约束罚函数做功 EMW 质量缩放功 ETotal 总能量 总能量为常值 100 在分析过程中允许的最大误差为 1 对于验证能量平衡主要采用伪应变能占总内能的比例 EI EE EP ECD EA EE 可恢复弹性应变能 EP 非弹性 如塑性 耗散能 ECD 蠕变或者粘弹性耗散能 EA 伪应变能 主要是由沙漏因素和壳单元与梁单元的横向剪应变储存的能量 按照 ABAQUS Explicit 能量评估的标准 伪应变能在总的内能力所占的比例不得高于 5 在分析完成之后先单一输出各个能量 再利用 ABAQUS 后处理中的 Operate on XY Date 对 能量进行叠加运算即可 4 举例分析 1 模型 如图所示为一方形钢板 中间均匀分布 6 道加强筋 在极短的瞬间受到爆炸产生的冲击载荷 采用 ABAQUS Explicit 解决爆炸载荷问题 方板模型和几何尺寸 方板几何尺寸 1 1m 加强筋之间的距离均为 0 25m 方板的厚度 0 012m 加强筋厚度 0 006m 高度 0 005m 2 材料 钢板为弹塑性材料 基本参数如下 7800Kg m 3 弹塑性数据 E 2 1 10 11 0 3 真是应力 Pa 真实塑性应变 3 10 8 0 3 5 10 8 0 025 3 75 10 8 0 1 3 94 10 8 0 2 4 10 8 0 35 表 1 塑性应力应变数据 表 2 爆炸载荷与时间的关系 时间 载荷大小 Pa 0 0 1 10 3 4 10 6 1 10 2 4 10 6 2 10 2 0 5 10 2 0 由上述数据可知 当实际应力大于 3 10 8Pa 时 产生塑性应变 PEEQ 值大于零 因此塑性应 变耗散能大于零 3 受载情况 钢板的四边六个自由度全部限制 上表面受到瞬间爆炸载荷 爆炸载荷随时间变化的大小如 表 2 受载模型 网格模型 4 网格划分 对于承受冲击载荷的情况只能选择线性积分单元 由于该模型的形状相对规则且为壳单元 故使用线性四边壳单元 S4R 先粗划网格 方板网格布置 20 20 加强筋网格布置 20 2 并根据网格的尺寸大小计算稳定时间极限 L e 0 05 C d E 5 2 10 3 stable t d e C L 9 6 10 6s 5 分析结果如下 应力云图 能量平衡图 分析过程的十个临界单元增量步时间 Most critical elements Element number Rank Time increment Increment ratio Instance name 5 1 4 093231E 06 1 000000E 00 PART 1 1 9 2 4 093231E 06 1 000000E 00 PART 1 1 10 3 4 093231E 06 1 000000E 00 PART 1 1 19 4 4 093231E 06 1 000000E 00 PART 1 1 20 5 4 093231E 06 1 000000E 00 PART 1 1 25 6 4 093231E 06 1 000000E 00 PART 1 1 39 7 4 093231E 06 1 000000E 00 PART 1 1 40 8 4 093231E 06 1 000000E 00 PART 1 1 45 9 4 093231E 06 1 000000E 00 PART 1 1 49 10 4 093231E 06 1 000000E 00 PART 1 1 由上述的时间增量步可知 增量步时间的数量级为 10 6s 这与理论上要求的数量级保持一致 而且根据能量平衡理论的要求 ALLAE ALLIE 5 由能量平衡图可看出 ALLAE ALLIE 在 3 4 之间 符合能量平衡要求 在这里为了验证网格的大小对增量步时间和能量平衡的影响 进一步细化网格 网格加倍方 板 40X40 加强筋 40X4 L e 0 025 C d E 5 2 10 3 stable t d e C L 4 8 10 6s 分析结果如下图 细网格模型 细化网格后应力云图 细化网格的能量平衡图 伪应变能占总内能的比值 本文主要讨论了 ABAQUS 显式动力学的基本理论和求解的问题以及具体的分析步骤和注意 要点 2 在进行显式动力学分析时 材料的波速 网格单元的尺寸大小都会影响计算结果的准确性 和临界时间 因而必要时可以采用质量缩放和细化网格等方式加速求解 3 能量平衡理论主要用于检验沙漏和横向剪应变产生的伪应变能所占内能的比例 进一步验 证分析的准确性 4 通过具体的爆炸载荷实例说明以上理论的正确性和诸