基于累积剪切应变损伤变量的微观疲劳裂纹扩展模拟方法

技术领域

该发明涉及疲劳数值模拟技术领域，尤其涉及一种基于累积剪切应变损伤变量的微观疲劳裂纹扩展模拟方法。

背景技术

疲劳是影响结构寿命与可靠性的重要因素。据统计，80％以上的工程结构破坏与疲劳有关。十九世纪初期，发生于铁路机车、铁轨以及桥梁的疲劳断裂最早引起了人们对疲劳问题的关注。近年来，由于航空技术的迅速兴起，对材料疲劳失效问题进行深入的研究愈加迫切。

材料疲劳失效行为的演化过程往往由应力集中、损伤、微观疲劳裂纹萌生、微观疲劳裂纹扩展、宏观疲劳裂纹扩展和失效组成，随着研究的深入，人们发现仅凭宏观尺度下的研究还不足以完全理解疲劳行为，微观尺度下的疲劳裂纹扩展行为逐渐成为研究者们关注的焦点。

对微观尺度下的裂纹扩展行为展开研究时，比较常用的方法有实验法和有限元法，但疲劳实验具有成本高、耗费长等缺点，且采用实验的方法时，往往需要从结构中进行取样，这对结构本身就是一种损害。随着计算能力的不断提高，有限元法越来越显示出其独特的优越性，尤其是晶体塑性有限元技术，为复杂条件下的疲劳研究提供了重要支撑，特别对于材料的微观力学行为和显微组织影响下的变形机理贡献尤为显著。

裂纹扩展行为的研究过程中，裂纹尖端网格的处理作为关注的焦点所在，常用的方法有节点释放技术、嵌入非连续模型、网格重划分技术、虚拟裂纹闭合技术、和内聚力模型方法。这些方法往往对于网格的质量要求过高，且在不连续力学问题中易发生不收敛情况。随着计算机技术的发展，无网格法和扩展有限元法应运而生，由于扩展有限元法不依赖于网格的划分，在求解不连续力学问题中较为高效，可实现裂纹的全自动扩展模拟，故在对疲劳裂纹扩展行为的研究中具有独特的优势。

在现有的、基于扩展有限元法的、微观尺度疲劳裂纹扩展的研究中，已实现了对单晶材料疲劳裂纹扩展的数值模拟，但基于单晶材料预测的疲劳裂纹扩展路径呈一条直线，不能反映晶粒取向影响下的裂纹偏转现象，也不能反映滑移和位错对裂纹扩展速率的影响。

基于此，为了准确描述多晶材料在微观尺度下的疲劳裂纹扩展行为，本发明提供了一种结合晶体塑性有限元技术和扩展有限元技术的、基于累积剪切应变损伤变量的、多晶材料微观疲劳裂纹扩展数值模拟方法。

发明内容

该发明拟解决的技术问题是针对以上不足，提供一种基于累积剪切应变损伤变量的微观疲劳裂纹扩展模拟方法，本方法可通过数值模拟实现多晶材料的微观疲劳裂纹扩展，并获得与微观疲劳裂纹扩展相关的一系列关键指标，如循环周次、疲劳裂纹扩展速率等。

为解决以上问题，该发明采用的技术方案如下：

一种基于累积剪切应变损伤变量的微观疲劳裂纹扩展模拟方法，包括以下步骤：

S1、通过Matlab语言生成多晶材料的Voronoi图和对应的随机晶粒取向；

S2、在ABAQUS中通过布尔运算对多晶材料建立代表性体积单元，通过Python 语言读取S1中Voronoi图的拓扑信息和每个晶粒的取向信息，随同材料参数批量赋予到每个晶粒中；

S3、编写用户自定义损伤子程序UDMGINI，将其嵌入至用户自定义材料子程序UMAT中；

S4、将UMAT子程序嵌入ABAQUS软件，进行有限元计算，生成每一增量步对应的应力增量和塑性应变增量，并基于此对当前时刻的应力应变状态值进行更新；

S5、UMAT子程序将应力应变状态、增量信息和材料参数作为输入变量，生成并旋转得到全局弹性矩阵，进一步迭代求解得到晶体塑性结果参量，并将其储存至雅可比矩阵DDSDDE和可共享的状态变量组STATEV中，完成材料响应的生成；

S6、UDMGINI子程序对S5中包含结果参量的STATEV组进行调用，并对裂纹尖端前方单元的累积剪切应变及最活跃滑移系进行计算与判定，累积剪切应变为：

式中：γ为累积剪切应变，是塑性剪切应变率

S7、若最活跃滑移系a上的累积剪切应变γ

S8、UDMGINI子程序将S7中的扩展信息反馈至ABAQUS求解器中；

S9、ABAQUS求解器根据反馈信息对各个单元的刚度矩阵进行组装；

S10、对已组装的刚度矩阵进行有限元求解；

S11、对求解结果的收敛性进行判定：若结果收敛，则根据单元平均解判断并更新裂纹集合，成功实现裂纹扩展，并将结果储存以进行下一个增量步的迭代运算；若不收敛，则自动采用较小增量步重复以上步骤直至收敛或达到最小增量步长为止。

优选的，所述S3中，UMAT子程序框架基于晶体塑性理论，采用基于滑移的幂函数型晶体塑性本构方程，硬化方式采用各向同性硬化；UDMGINI子程序通过 FORTRAN语言编写，通过可共享的状态变量组STATEV与UMAT子程序完成衔接工作。

优选的，所述S5中，通过Newton-Raphson迭代法实现隐式求解并解决塑性滑移增量的迭代问题。

优选的，所述S6中，面心立方结构模型的随机晶粒中包含4个滑移面和12 个滑移系，4个滑移面在循环加载工况下均被激活且在不同区域占据主导地位， UDMGINI子程序通过调取UMAT子程序中分别代表12个滑移系的状态变量： STATEV(9*NSLPTL+1)至STATEV(10*NSLPTL)，通过比较12个滑移系的独立累积剪切应变并选择包含最大独立累积剪切应变的滑移系作为最活跃滑移系。

该发明采用以上技术方案，与现有技术相比，具有以下优点：

1.该模型兼顾了单晶材料和多晶材料，能够描述具有面心立方结构、体心立方结构和六方密排结构材料的微观疲劳裂纹扩展行为。

2.模型预测的疲劳裂纹扩展路径不再只是一条直线，而是能准确描述晶粒取向影响下的裂纹偏转现象，并能反映滑移和位错对裂纹扩展速率的影响。

3.模拟寿命与实际寿命之间可通过缩放系数进行关联，可以在更少的模拟周期下完成对实际寿命的准确预测。

下面结合附图和实施例对该发明进行详细说明。

附图说明

图1为该发明实施例中多晶材料微观疲劳裂纹扩展数值模拟方法流程图；

图2为该发明实施例中代表GH4169合金微观组织的代表性体积单元；

图3为该发明实施例中SENT试样几何形状和边界条件示意图；

图4为该发明实施例中SENT试样疲劳载荷波形示意图；

图5为该发明实施例中SENT试样裂纹扩展图。

具体实施方式

下面将结合说明书附图和具体的实施例对本发明所述的基于累积剪切应变损伤变量的微观疲劳裂纹扩展模拟方法做进一步的详细说明。

实施例：

本实施例以对含裂纹SENT试样进行微观疲劳裂纹扩展数值模拟为例，对本发明进行详细说明。

一种基于累积剪切应变损伤变量的微观疲劳裂纹扩展模拟方法，如图1所示，包括以下步骤：

步骤一、通过Matlab语言生成多晶材料的Voronoi图和对应的随机晶粒取向。

步骤二、在ABAQUS中通过布尔运算对多晶材料建立代表性体积单元，通过 Python语言读取S1中Voronoi图的拓扑信息和每个晶粒的取向信息，随同材料参数批量赋予到每个晶粒中。

其中代表性体积单元的示意图如图2所示，整体尺寸为200μm×200μm，共包含150个晶粒和22534个元素。SENT试样的几何形状和边界条件如图3所示，通过对试样上下表面施加相同的位移载荷控制加载，疲劳载荷加载波形如图4所示，加载波形为三角波，其中每2个time为一个循环。预制初始裂纹长度a＝20μm，为了消除刚体位移的影响，对左侧边界沿x方向的位移进行约束，且允许其在y方向上自由移动。

对SENT试样进行网格化：网格类型采用四节点平面应变单元CPE4，最小网格尺寸为0.0012μm。

SENT试样材料选用GH4169合金，本模型中用到的GH4169合金的化学成分及晶体塑性本构参数如下表1和表2所示。

表1.本模型中用到的GH4169合金的化学成分表

表2.本模型中用到的GH4169合金的晶体塑性本构参数

步骤三、编写用户自定义损伤子程序UDMGINI。

UDMGINI子程序通过Fortran语言编写，ABAQUS用户手册中提供了编写该 UDMGINI子程序的固定格式，如下所示：

将编写好的UDMGINI子程序嵌入至用户自定义材料子程序UMAT中。UMAT子程序的固定格式如下所示：

UMAT子程序框架基于晶体塑性理论，采用基于滑移的幂函数型晶体塑性本构方程，硬化方式采用各向同性硬化，可通过临界剪切应力的演化实现：

h

式中，g

步骤四、将UMAT子程序嵌入ABAQUS软件，进行有限元计算，生成每一增量步对应的应力增量和塑性应变增量，并基于此对当前时刻的应力应变状态值进行更新。

步骤五、UMAT子程序将应力应变状态、增量信息和材料参数作为输入变量，生成并旋转得到全局弹性矩阵，进一步迭代求解得到晶体塑性结果参量，并将其储存至雅可比矩阵DDSDDE和可共享的状态变量组STATEV中，完成材料响应的生成。

步骤六、UDMGINI子程序对S5中包含结果参量的STATEV组进行调用，并对裂纹尖端前方单元的累积剪切应变及最活跃滑移系进行计算与判定，累积剪切应变为：

式中：γ为累积剪切应变，是塑性剪切应变率

本模型选用材料GH4169合金为面心立方结构，模型随机晶粒中包含4个滑移面和12个滑移系，4个滑移面在循环加载工况下均被激活且在不同区域占据主导地位，UDMGINI子程序通过调取UMAT子程序中分别代表12个滑移系的状态变量：STATEV(9*NSLPTL+1)至STATEV(10*NSLPTL)，通过比较12个滑移系的独立累积剪切应变并选择包含最大独立累积剪切应变的滑移系作为最活跃滑移系。

步骤七、若最活跃滑移系a上的累积剪切应变γ

步骤八、UDMGINI子程序将S7中的扩展信息反馈至ABAQUS求解器中。

步骤九、ABAQUS求解器根据反馈信息对各个单元的刚度矩阵进行组装。

步骤十、对已组装的刚度矩阵进行有限元求解。

步骤十一、对求解结果的收敛性进行判定：若结果收敛，则根据单元平均解判断并更新裂纹集合，成功实现裂纹扩展，并将结果储存以进行下一个增量步的迭代运算；若不收敛，则自动采用较小增量步重复以上步骤直至收敛或达到最小增量步长为止。

最后应说明的是：以上所述仅为该发明的优选实施例而已，并不用于限制该发明，尽管参照前述实施例对该发明进行了详细的说明，对于本领域的技术人员来说，其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分技术特征进行等同替换。凡在该发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在该发明的保护范围之内。