一种矩形板中心斜裂纹疲劳扩展寿命预测的简化方法

摘要

一种矩形板中心斜裂纹疲劳扩展寿命预测的简化方法，它有五大步骤：步骤一：参数化建模：步骤二：划分网格：步骤三：有限元计算及结果处理：步骤四：通过循环调用实现裂纹自动扩展模拟；步骤五：预测裂纹扩展寿命。本发明基于Python语言，进行ABAQUS二次开发，采用参数化建模，编制了建模与分析的命令流程序包，可以实现对带中心斜裂纹的矩形板裂纹扩展的自动模拟，功作量小，精度高。

说明书

技术领域

本发明属于结构与强度及其相关的技术领域，具体涉及一种矩形板中心斜裂 纹疲劳扩展寿命预测的简化方法。

背景技术

传统疲劳裂纹扩展分析通常将裂纹类型假设为纯I型裂纹，在裂纹扩展过程 中，裂纹扩展方向不发生改变，可通过事先建立应力强度因子与裂纹尺寸之间的 经验关系，从而简化裂纹扩展分析计算。而在工程实际中，复合型裂纹是常见形 式，对此分析则要复杂得多。首先，裂纹扩展率模型更为复杂，涉及多个类型的 裂纹尖端参数。其次，复合型裂纹扩展方向事先未知，且在扩展过程不断变化， 无法事先建立应力强度因子与裂纹尺寸、形状等参数之间的经验关系。因此对复 合型疲劳裂纹扩展过程的分析，数值模拟是一种有效的方法。

常用的疲劳裂纹扩展预测模型是基于Paris公式的。以二维情况为例，通过 对裂纹体的有限元分析，计算裂纹尖端应力强度因子幅值，通过指定一个小的疲 劳周数增加量△N获得该裂纹扩展方向上的裂纹长度增量△a，以复合型断裂准 则—最大周向应力准则为基础计算裂纹的扩展方向，重新定义一个新的裂纹尖 端，最后相对于这一新裂纹尖端的有限元网格被重新划分，重新计算新应力强度 因子幅值，这样一步一步地跟踪裂纹的扩展。

对于常用的通用型有限元软件如ABAQUS，由于并不能实现裂纹的扩展模 拟，按照上述预测模型的思路，需要对其进行二次开发，利用其脚本语言Python 编程来达到模拟的目的。

复合型裂纹的精确分析计算往往是十分复杂的，根据（参考文献：邢文珍， 刘雪惠.I-II复合型裂纹疲劳扩展规律研究[J].固体力学学报.1983）中的投影法 虽然方法简单，但其精度偏于保守。

发明内容

针对矩形板中心斜裂纹在复合载荷下疲劳裂纹扩展问题，通过对裂纹扩展路 径及寿命预测的对比分析，提出了一种简化方法来模拟疲劳裂纹扩展并预测疲劳 扩展寿命，借助Python语言进行参数化建模，编制了建模与分析的命令流程序， 更新模型重新划分网格，逐步实现疲劳裂纹的扩展模拟。

本发明是一种矩形板中心斜裂纹疲劳扩展寿命预测的简化方法，通过Python 语言编程序实现。程序的功能主要由一系列子程序实现：创建模型子程序用来参 数化几何建模；定义裂纹子程序的功能是定义裂纹，在模型中自动选取裂尖及裂 缝；网格划分子程序用于为模型分网，将裂尖前缘划分为奇异单元来反映裂纹的 奇异性，周围划分成用于计算围线积分的环状单元；结果处理子程序来实现对计 算结果数据如应力强度因子的提取及处理；计算新裂尖子程序，用来计算新的裂 纹尖端坐标。

该方法具体步骤如下：

步骤一：参数化建模：建立带中心斜裂纹矩形板的有限元模型。如图2所示， 首先调用创建模型子程序，建立宽为2W，高为2H，中心在坐标原点的矩形模型； 再调用定义裂纹子程序，定义长度为2a，与x轴夹角b，在矩形模型中心的裂纹， （首先计算出裂纹边缘两点坐标x₁=a·cosb，y₁=a·sinb，x₂=-a·cosb，y₂=-a·sinb， 将两点连成裂纹）。

步骤二：划分网格：由于裂纹的存在，裂尖附近的网格需要奇异化处理，以 反应裂尖应力场，以及实现应力强度因子的计算。利用网格划分子程序用于为模 型分网，主要分为裂尖区域以及远离裂尖区域，将裂尖前缘划分为奇异单元来反 映裂纹的奇异性，周围划分成用于计算围线积分的环状单元；远离裂尖区域由于 不影响裂尖应力场的计算，则划分为普通的三角形单元，以提高程序计算效率。

步骤三：有限元计算及结果处理：提交作业，完成有限元计算分析。利用 ABAQUS计算应力强度因子及所需要的数据，通过结果处理子程序来实现对计算结 果数据如应力强度因子的提取。

步骤四：通过循环调用实现裂纹自动扩展模拟。判断每次裂纹扩展时所得应 力强度因子K是否达到临界值，如果达到，转到步骤五；如果没有达到，利用计 算新裂尖子程序计算新的裂纹尖端坐标，继续步骤一。

在这里，由于采用本文提出的简化方法，将裂纹的后续扩展简化为一条直线， 故不再需要计算裂纹扩展方向，只需根据给定的裂纹扩展增量（步长）计算新的 裂尖坐标即可。为了提高计算分析效率且保证计算精度，计算新裂尖子程序中， 包含了变步长的处理方法，即在初始的裂纹扩展中给定较小的步长，而后面的由 于扩展寿命所占总扩展寿命的比重较小，故可采用较大的步长。

步骤五：预测裂纹扩展寿命。上述过程实现了复合载荷下裂纹扩展的数值模 拟，但我们的目的为了预测裂纹扩展寿命的预测。疲劳裂纹扩展服从Paris公式

da/dN＝C(Δk)^m

这里da/dN是疲劳裂纹扩展率，A和M是材料参数，ΔK是应力强度因子K的幅值。 写成差分形式，可以得到

$\mathrm{\Delta N}=\frac{\mathrm{\Delta \alpha }}{C{\left(\mathrm{\Delta K}\right)}^m}$

如果已知裂纹扩展增量Δa，则可以利用方程得到扩展此步长经历的裂纹的扩展寿 命ΔN。

于是，裂纹由初始长度a₀扩展到临界值a_c，所经过的循环周次N可由积分得 到

$N=\mathrm{\Sigma \Delta N}=\Sigma \frac{{\mathrm{\Delta \alpha }}_i}{C·\left\{{\left({\mathrm{\Delta K}}_i\right)}^n\right\}}$

Δa_i是某一扩展步长，ΔK_i是此步长所对应的应力强度因子幅值。通过步骤三 中，对有限元计算结果的处理，我们可得到每次扩展步长Δa_i的ΔK_i，从而计算得 到裂纹扩展寿命N。流程图如图4所示。

本发明方法的优点如下：

1）基于Python语言，进行ABAQUS二次开发，采用参数化建模，编制了 建模与分析的命令流程序包，可以实现对带中心斜裂纹的矩形板裂纹扩展的自动 模拟。

2）提出的简化方法模拟结果发现与未简化前相差很小（在1％左右），但大 大减小了工作量。如果不涉及对裂纹扩展具体轨迹的研究，对于计算裂纹扩展寿 命及相关精度是足够的。

附图说明

图1是带中心斜裂纹矩形板示意图

图2是裂纹扩展简化模拟示意图（W为模型半宽，H为模型半高，α为裂纹与 x方向夹角）

图3是不同裂尖单元示意图（a为裂纹长度的单元数）

图4是本发明流程框图

图5扩展步长对裂纹扩展路径影响示意图

图6扩展步长对寿命预测的影响示意图

图7不同裂纹倾斜角的裂纹扩展路径示意图

图8 45°倾斜角裂纹模拟对比（1为未简化时的模拟结果，2为简化结果，3 为投影法结果）

图9 30°倾斜角裂纹模拟对比（1为未简化时的模拟结果，2为简化结果）

图10 60°倾斜角裂纹模拟对比（1为未简化时的模拟结果，2为简化结果）

具体实施方式

见图4，本发明一种矩形板中心斜裂纹疲劳扩展寿命预测的简化方法，通过 Python语言编写各个步骤的子程序实现，该方法具体步骤如下：

步骤一：参数化建模：建立带中心斜裂纹矩形板的有限元模型。如图2所示， 首先调用创建模型子程序，建立宽为2W，高为2H，中心在坐标原点的矩形模型； 再调用定义裂纹子程序，定义长度为2a，与x轴夹角b，在矩形模型中心的裂纹， （首先计算出裂纹边缘两点坐标x₁=a·cosb，y₁=a·sinb，x₂=-a·cosb，y₂=-a·sinb， 将两点连成裂纹）。

步骤二：划分网格：由于裂纹的存在，裂尖附近的网格需要奇异化处理，以 反应裂尖应力场，以及实现应力强度因子的计算。利用网格划分子程序用于为模 型分网，主要分为裂尖区域以及远离裂尖区域，将裂尖前缘划分为奇异单元来反 映裂纹的奇异性，周围划分成用于计算围线积分的环状单元；远离裂尖区域由于 不影响裂尖应力场的计算，则划分为普通的三角形单元，以提高程序计算效率。

步骤三：有限元计算及结果处理：提交作业，完成有限元计算分析。利用 ABAQUS计算应力强度因子及所需要的数据，通过结果处理子程序来实现对计算结 果数据如应力强度因子的提取和处理，即得到每次扩展步长Δa_i的ΔK_i。

步骤四：通过循环调用实现裂纹自动扩展模拟。判断每次裂纹扩展时所得应 力强度因子K是否达到临界值，如果达到，转到步骤五；如果没有达到，利用计 算新裂尖子程序计算新的裂纹尖端坐标，继续步骤一。

在这里，由于采用本文提出的简化方法，将裂纹的后续扩展简化为一条直线， 故不再需要计算裂纹扩展方向，只需根据给定的裂纹扩展增量（步长）计算新的 裂尖坐标即可。为了提高计算分析效率且保证计算精度，计算新裂尖子程序中， 包含了变步长的处理方法，即在初始的裂纹扩展中给定较小的步长，而后面的由 于扩展寿命所占总扩展寿命的比重较小，故可采用较大的步长。

步骤五：预测裂纹扩展寿命。上述过程实现了复合载荷下裂纹扩展的数值模 拟，但我们的目的为了预测裂纹扩展寿命的预测。疲劳裂纹扩展服从Paris公式

da/dN＝C(Δk)^m

这里da/dN是疲劳裂纹扩展率，A和M是材料参数，ΔK是应力强度因子K的幅值。 写成差分形式，可以得到

$\mathrm{\Delta N}=\frac{\mathrm{\Delta \alpha }}{C{\left(\mathrm{\Delta K}\right)}^m}$

如果已知裂纹扩展增量Δa，则可以利用方程得到扩展此步长经历的裂纹的扩展寿 命ΔN。

于是，裂纹由初始长度a_o扩展到临界值a_c，所经过的循环周次N可由积分得 到

$N=\mathrm{\Sigma \Delta N}=\Sigma \frac{{\mathrm{\Delta \alpha }}_i}{C·\left\{{\left({\mathrm{\Delta K}}_i\right)}^n\right\}}$

Δa_i是某一扩展步长，ΔK_i是此步长所对应的应力强度因子幅值。通过步骤三 中，对有限元计算结果的处理，我们可得到每次扩展步长Δa_i的ΔK_i，从而计算得 到裂纹扩展寿命N。

为了对上述模型进行验证，我们对一带中心斜裂纹的矩形板（如图1）进行 了疲劳裂纹扩展模拟。几何参数分别是初始裂纹长度为2a，矩形板宽2W，裂纹 倾斜角度α。材料参数如下：

E＝7250kg/(mm)²，υ＝0.321，

K_Ic＝116kg/(mm)²，σ＝15.3kg/(mm)²

A＝1.039×10^-10，m＝2.7438，R＝0.048

a=7mm    W=35mm

根据所述疲劳裂纹预测模型，其关键是应力强度因子的求解。为了考察有限 元计算应力强度因子的准确性和可靠性，有必要对有限元计算结果进行评估，并 确定划分网格尺寸对其影响。

表1是裂纹未扩展时所求的应力强度因子与文献：中国航空研究院.应力强 度因子手册.科学出版社(1981)中的手册值的比较，验证了ABAQUS求解应力 强度因子的准确性。由于应力强度因子手册中倾斜角有限，这里选取45°和75°时 进行比较。

表1应力强度因子有限元解与手册值

可见，ABAQUS可以精确地求解出相应的应力强度因子。

这里选取了四种网格尺寸（以裂尖网格单元特征长度与裂纹半长a之比）进 行裂纹尖端网格划分：1/2、1/5、1/10、1/20。计算模型选择中心裂纹（即裂纹 倾斜角为0°）进行计算。

计算模型选择中心裂纹（即裂纹倾斜角为0°）进行计算。计算结果如表2.

表2不同裂尖单元尺寸对应力强度因子的影响

采用有限元方法进行应力强度因子计算时，可以选取裂尖单元特征尺寸与裂 纹半长之比小于1/10进行计算，这时计算精度较为满意。

为了考察给定裂纹扩展步长对裂纹扩展数值模拟结果的影响，对于初始裂纹 倾斜角为45°情况，选择了几种步长进行计算(步长/半裂纹长度)：0.01、0.03、 0.05、0.10、0.15。裂纹扩展路径模拟对比结果如图5所示（图中只显示裂纹 扩展路径，初始裂纹未显示），扩展步长对寿命预测的影响如图6所示。从中可 以看出，取扩展步长为原裂纹半长的0.05以内时，裂纹扩展路径基本相同，对 寿命的影响也很小，因此，在模拟裂纹扩展时，裂纹扩展步长可以选取裂纹半长 的0.05之内。

从图7的裂纹扩展轨迹可以看出，斜裂纹总是倾向于垂直加载方向扩展，初 始倾斜角越小，裂纹扩展方向更倾向于水平方向，而且裂纹扩展轨迹除了开始的 部分变化较大，大部分基本是沿水平方向的。这与实验现象是吻和的。

由于裂纹总是倾向于垂直加载方向扩展，并且从上图的裂纹实际轨迹可以看 出裂纹扩展轨迹基本平行于水平方向，数值模拟结果显示，裂纹在初始扩展时其 偏转角较大但这一阶段在整个扩展过程中所占比例很小（只在扩展开始的几个增 量步内如此），而后续扩展角基本就是平行于水平方向，我们不妨简化裂纹就是 沿着水平方向扩展的，如图2所示。图3是不同裂尖单元示意图。

以45°倾斜角为例，我们在此简化基础上对裂纹扩展进行了模拟，并与另一 种简化方法“投影法”进行了对比，如图8所示（其中，1为未简化时的模拟结果， 2为简化结果，3为投影法结果。图9、图10为30°和60°时简化前后的对比）。

模拟结果表明，这种简化与直接模拟结果十分相近（误差均在1％左右），但 大大减小了工作量，可以在工程实际中应用。