基于ABAQUS标准扩展有限元法的不良裂缝构型影响分析

方淑君, 王 涛，聂念从，张凌瑞

(中南大学 土木工程学院，湖南 长沙 410075)

混凝土被广泛应用于建筑行业，在施工中混凝土水化热、养护不当、配筋过多等诸多因素都有可能导致裂缝或者空洞的产生，而结构运营期间钢筋锈蚀、疲劳载荷等不利载荷作用下也会产生裂缝.在实际工程中混凝土结构产生裂缝的因素多、裂缝构型多样、裂缝处应力情况复杂[1].传统有限元方法研究断裂问题时有着很多固有缺陷[2],例如，裂缝处网格需要加密，计算量太大；裂缝的动态扩展需要不断地网格重构等.Belytschko等[3-5]提出了扩展有限元法(XFEM),其原理是在单位分解思想下，添加了加强函数对不连续边界进行描述.这种方法在常规有限元基础上，独立于网格计算裂缝这类强不连续问题.方修君等[6]将XFEM应用于混凝土梁复合型开裂过程的数值模拟，证实XFEM是一种模拟准脆性材料复合开裂的有效方法.胡少伟等[7]利用XFEM较好地模拟了不同参数混凝土三点弯曲梁裂缝起裂、稳定扩展与失稳破坏的全过程.杨刚等[8]基于XFEM下的混凝土弹塑性本构关系，模拟了钢筋混凝土柱的荷载-位移曲线和开裂过程.已有的研究论证了将混凝土作为准脆性材料采用XFEM进行断裂分析结果良好.

标准扩展有限元方法利用了单位分解思想[9-10]，而单位分解插值会造成线性相关性，富集节点位移函数分解插值计算会使得总体刚度矩阵高度病态[11-14].目前，ABAQUS采用的是标准拓展有限元法，在对实际工程进行裂缝模拟时，裂缝发展可能会发展至网格边角处，刚度矩阵病态造成模拟结果失真.本文的裂缝研究基于线弹性断裂力学，研究XFEM的Ⅰ型开裂对网格的敏感性.模型选取三点弯曲梁试验进行模拟分析,通过模拟结果与试验数据对比，分析网格尺寸与不良裂缝构型两个因素对ABAQUS标准拓展有限元法的影响.为了对比不同的裂缝发展路径，引入离散弗雷歇距离定量地进行比较分析.

1 XFEM基本理论与离散弗雷歇距离

扩展有限元法是在常规的有限单元法的位移函数中加入一些加强函数以反映间断性.当加强基函数为多项式时矩阵奇异，而且随网格加密或基函数阶次提高，病态趋于严重，稳定性更差.

1.1 节点富集函数

用于增加节点自由度的加强函数通常包括裂纹尖端附近渐进函数和间断函数.富集节点的位移向量函数u表示为

(1)

1.2 空间不连续场积分

对包含裂缝的单元采取高斯积分法.单元被裂缝分割后，x，y方向的最小长度被称为特征长度，在单元内均匀分布积分点，而且保证在特征长度上至少有两个积分点.如图1所示，图上的积分点个数是7×8，对于不同的单元，积分点的布置须分别计算[15].

1.3 离散弗雷歇距离

弗雷歇距离考虑了时序空间，它的定义如图2所示.假设有两人分别按照各自的路径规划行走，两人之间有一根长短可以伸缩的绳子.弗雷歇距离就是在两个人都能够按照各自轨迹走到终点的过程中计算满足条件的最长的绳子长度[16].

假设曲线R有r个数据点，曲线Q有q个数据点.σ(R) 和σ(Q)分别表示两轨迹点按特定规律排序的集合，则有σ(R)=(u1,…,ut)和σ(Q)=(ν1,…,νt).可以得出以下点对L： (u1,ν1),(u2,ν2),…,(ut,νt).其中,0

(2)

离散弗雷歇距离表示为

δdF(R,Q)=min‖L‖.

(3)

通过固定位移荷载的值，提取裂缝曲线的离散数据点，并计算整个时程内裂缝曲线之间的δdF.通过δdF描述不同模拟方法导致的裂缝构型的差异大小.

1.4 不良裂缝的定义

2 试验设置与数值模拟

2.1 三点弯曲混凝土梁试验

本文试验选用三点弯曲混凝土梁，详细数据与试验结果采用王青原等[18]与唐宇翔等[19]发表的试验数据.试验选用材料为硅酸盐水泥CEM I 52.5，骨料的最大粒径为10 mm.混凝土配合比(水∶水泥:细骨料∶粗骨料)=0.58∶1∶2.73∶2.73.梁的尺寸(长×高×厚)为710 mm×150 mm×80 mm，梁跨度S为600 mm，预制的裂缝长度为45 mm，宽度为2 mm.试验梁的混凝土弹性模量为23.8 GPa.加载方式为位移控制加载，其速率为0.01 mm/min.在切缝口放置位移计测量裂缝口张开位移(CMOD).试验分为三组，分别为A-1,A-2,A-3.三维裂缝通过0,1,2号截面转化为二维平面对裂缝形状进行描述，试验布置与截面见图3.

2.2 数值模拟

采用ABAQUS的XFEM对三点弯曲混凝土梁数值模拟，模型尺寸与试验梁尺寸一致.混凝土材料参数见表1，模型设置简支约束.加载采用跨中位移加载，竖直向下位移为0.8 mm.单元类型采用8节点等参减缩积分单元(C3D8R).开裂准则采用最大主应力准则，采用基于能量的线性软化本构关系表示断裂过程的软化特性.初始裂缝接触选用硬接触且无摩擦.在部分收敛困难的模型中采用增加最大迭代次数、引入黏结阻尼等方法改善模型收敛性，ABAQUS数值建模见图4.

表1 混凝土材料参数

2.3 模型设置

通过改变不良系数u，分析不良裂缝位置的影响而设置模型组1-1～1-5.为分析网格尺寸的影响设置模型组1-1～5-1.模型网格信息见表2，网格尺寸数据见图5.

表2 模型网格参数信息

3 数据分析

3.1 网格尺寸影响分析

在裂缝构型较好，即不良系数u=1时，讨论网格尺寸在裂缝计算中的影响.

3.1.1 裂缝发展路径对比

试验裂缝与模型裂缝发展的路径见图6.

试验中，在拉应力作用下，跨中裂缝为Ⅰ型裂缝，总体趋势是竖直向上，但遇到硬度较大的骨料会发生弯折.在不良系数u=1时，各模型0，1，2号截面裂缝均是竖直向上发展的，原因是模型假设为匀质实体没有骨料的作用，裂缝不会弯折.以模型1-1的裂缝发展时程曲线为基准，求解离散弗雷歇距离δdF，见表3.分析数据可见，大尺寸网格会使裂缝长度增值较大，数据波动较大求得δdF也较大.但是在裂缝发展趋势相同情况下，需要控制网格尺寸以对δdF值进行控制.

表3 裂缝路径δdF数据

通过离散弗雷歇距离可以定量地描述裂缝发展路径的差异，摆脱了定性描述实际裂缝与模拟裂缝相似程度的模糊性.

3.1.2P-δ曲线与P-CMOD曲线

将不同网格尺寸且u=1的模型数据与试验数据进行比较，反力曲线与裂缝加载位移曲线见图7～图8.

在反力P-δ曲线图与P-CMOD曲线图中，试验数据曲线与模型模拟曲线走势均相同.但是在接近混凝土梁完全破坏时，模拟曲线下降更快、裂缝发展更快.原因是，在模拟时将混凝土梁假设为各向同性的匀质固体，而实际混凝土中的骨料咬合作用会在一定程度上抑制裂缝发展，即曲线下降得较为平缓；实际裂缝路径存在弯折，在裂缝面上存在法向应力和切应力，而在模拟时裂缝路径竖直向上，裂缝面不存在法向应力与切向应力.

试验与模拟结果部分数据见表4.对选取的参数最大反力PMax、临界裂缝张开位移CMODc和临界加载位移δc，模型1-1拟合最好且相对误差很小，而模型2-1与3-1拟合较差、相对误差较大，3.1.3节对其误差来源进行分析.

表4 试验与模拟部分数据表

3.1.3 初始裂缝端点与网格边界关系

模型2-1网格长度a=10 mm裂缝端点落在单元内部，而模型4-1和5-1调整网格长度使裂缝端点落在单元边界上.相关曲线见图9～图10，部分模型数据见表5.

模型4-1与模型5-1的最大反力PMax比模型2-1拟合更好，而前两个模型与模型2-1的区别仅为网格划分的差异.相较于模型2-1，模型4-1与5-1的网格尺寸变化较小，对结果的主要影响因素为裂缝端点位置.由于ABAQUS的XFEM方法有一定的简化，在进行裂缝动态分析时，裂缝产生即贯穿整个单元.当裂缝端点落在单元内时，该单元被认为完全破坏，在原有基础上削减了梁承载力.

表5 部分模型模拟数据

网格尺寸需考虑初始裂缝端点位置，不能使端点落入单元内.此外，网格尺寸大小也会影响裂缝发展的时程情况.

3.2 不良裂缝构型的影响分析

3.2.1 裂缝发展路径

在不良系数减小时裂缝发展路径会产生显著变化，且δdF与δ′dF值均增大，在u=0时，达到最大值.在构件、边界与荷载均关于xoy平面对称的条件下，δ′dF值过大且裂缝面并不对称.这种情况验证了随不良系数减小，单元刚度矩阵病态更严重，微小的节点位移差异导致两个截面裂缝发展差异巨大.同样也说明基于ABAQUS平台的XFEM网格划分会影响裂缝发展的结果，不良系数u值越接近0，模拟的裂缝发展路径失真情况越严重.

表6 裂缝路径δdF数据表

3.2.2P-δ曲线与P-CMOD曲线

将相同网格尺寸而不良系数u不同的模型数据与试验数据进行比较，反力曲线与裂缝张开位移曲线见图13～图14.部分模型模拟数据见表7.

在P-δ曲线图与P-CMOD曲线图中，试验数据曲线与模型模拟结果曲线在达到反力最大值段均拟合较好，但在梁接近破坏时表现出明显差异.当不良系数减小至0时，出现应力锁死现象，产生了虚假的承载力.在初始裂缝起裂阶段与初始裂缝位置无关，则最大反力值各模型值均误差较小.因此基于ABAQUS的XFEM在不良裂缝构型时仍可以得出较为准确的承载力数值.

表7 部分模型模拟数据

3.2.3 黏结阻尼的影响

表8 裂缝路径δdF数据表

4 结 论

1) 通过离散弗雷歇距离可以定量地描述裂缝发展路径的差异，摆脱了定性描述实际裂缝与模拟裂缝相似程度的模糊性.

2) 应用基于ABAQUS的XFEM分析裂缝扩展时，初始裂缝端点必须设置在单元边界处，否则会削弱构件的承载力.

3) 网格划分尺寸和裂缝构型会影响裂缝的发展.网格尺寸过大会影响裂缝时程发展.不良系数u值越接近0，模拟的裂缝发展路径失真情况越严重, 不良系数u值要控制在0.25～1内.

4) 裂缝构型不影响计算裂缝起裂荷载，在不良系数u较小时仍可以得到较准确的最大反力值.但当不良系数u接近0时，模型在构件即将完全破坏的阶段会存在应力锁死现象，产生虚假的承载能力.

5) 通过引入黏结阻尼，在提高收敛性的同时有效降低刚度矩阵病态引起的数值震荡，可以一定程度上提高裂缝发展路径的准确性.