基于XFEM的混凝土三点弯曲梁开裂数值模拟研究

胡少伟，鲁文妍

(1.南京水利科学研究院，江苏 南京210024;2.水文水资源与水利工程科学国家重点实验室，江苏 南京210024)

混凝土是复杂的非线性脆性材料，其断裂问题受到广泛关注. 自1961年混凝土断裂力学诞生以来，国内外学者开展了广泛的研究，提出许多混凝土断裂模型，如虚拟裂缝模型［1］、等效裂缝模型［2］、尺寸效应模型［3］、双K 断裂模型［4］等. 这些模型在描述混凝土断裂特性方面各具特色，同时也有其各自的局限性.近年来，由Belytschko 和Black 提出的扩展有限元法［5］得到不断发展，它将非连续富集技术应用到有限元中，以单位分解法［6］来确定裂缝扩展路径，具有需网格细化、允许裂缝在连续单元内部扩展等优点.笔者介绍了扩展有限元法的基本原理，给出了扩展有限元法分析断裂问题的实现方法，并基于最大主拉应力牵引损伤开裂准则和线性软化损伤模型，应用扩展有限元法模拟了含初始裂纹的混凝土三点弯曲梁的开裂扩展过程，展示了扩展有限元法在混凝土开裂扩展分析中的独特优势.

1 扩展有限元法的基本原理

扩展有限元法(Extended Finite Element Method，XFEM)利用裂纹富集技术来考虑裂缝尖端奇异应力场和裂缝面两侧位移的不连续性，并引入水平集函数来跟踪裂缝扩展路径和确定裂缝尖端位置［7］.它最基本的思想是用一些局部附加函数来改进传统有限元，使得裂缝几何独立于计算网格.XFEM 富集模式如图1所示.

图1 XFEM 富集模式

在裂缝所在区域(图1)定义了Heaviside 函数和裂缝尖端富集函数，其位移表达式为

式中:Ni为i 点常规有限元形函数;ui为常规有限元节点位移;H(x)为Heaviside 阶跃函数;Nj为j 点与Heaviside 函数相关的形函数;bj为被裂缝贯穿的单元结点附加自由度(图1中方点所示);Fl(x)为裂缝尖端富集函数;Nk为k 点与裂缝尖端弹性富集函数相关的形函数;ck为裂缝尖端所在单元的结点附加自由度(图1中圆点所示);I 为网格中所有节点的集合;J 为裂缝贯穿单元节点的集合;K 为裂缝尖端所在单元节点的集合.

对于没有加强的区域，式(1)位移函数退化为常规有限元位移函数:

裂缝贯穿单元只有Heaviside 阶跃函数加强，其位移函数为

裂缝尖端所在单元只有裂缝尖端函数加强，其位移函数为

扩展有限元应用水平集函数来追踪裂缝扩展迭代中裂缝面和裂缝尖端位置. Heaviside 阶跃函数通过裂缝面水平集函数φ(x)来描述裂缝面单元内部位移的不连续性，具体公式为

式中φ(x)为裂缝面水平集函数(如图2所示)，代表Heaviside 函数加强点相对于裂缝面的垂直距离(包含正负).

裂缝尖端水平集函数ψ(x)(如图2所示)代表裂缝尖端加强点相对于裂缝前沿边界的垂直距离(包含正负).裂缝尖端加强函数定义了裂缝尖端单元的不连续性，具体公式为

式中r 和θ 为局部裂纹尖端场坐标系统中的极坐标.

图2 裂缝水平集函数

2 裂缝扩展模拟方法

一般认为混凝土材料是脆性材料，裂纹萌生准则采用最大主应力准则(即当最大主应力值达到材料屈服应力值时裂纹萌生)，表达式为

应用最大周向应力准则来确定裂缝每一步扩展的方向，根据这一准则，裂缝扩展角计算如下:

式中KⅠ，KⅡ分别为Ⅰ型和Ⅱ型断裂应力强度因子.

3 单边切口混凝土三点弯曲梁模型数值分析

3.1 模型建立

混凝土三点弯曲梁的几何尺寸、受力及边界条件如图3所示.试件跨中位置设置一预制裂缝，初始裂缝长度为80 mm，即初始缝高比α0(α0=a0/h)为0.4.混凝土强度等级为C35，材料物理参数见表1，混凝土采用最大主拉应力牵引损伤开裂准则和基于能量的线性软化损伤模型.

图3 混凝土三点弯曲梁几何形式、受力及边界条件(单位:mm)

表1 混凝土材料物理力学参数

利用大型有限元软件ABAQUS，建立含有初始裂缝的混凝土三点弯曲梁数值模型，如图4所示.模型左端采用固定支座，右端采用滑动支座，跨中位置施加一竖向位移荷载，并约束除竖向以外的其他自由度.混凝土梁采用C3D8R 实体单元，对梁中部裂纹扩展区域网格进行适当加密. 加载点下垫块作用为防止加载过程中梁上加载面及支座处出现应力集中，但若将垫块与梁绑定则人为引入了过多的约束，而将垫块与混凝土梁之间设置接触则存在收敛问题，因此建模时在支座和加载点采用耦合点位移加载，耦合点位于垫块中心轴线上.

图4 混凝土三点弯曲梁数值模型

3.2 结果分析

利用所建立的模型对混凝土三点弯曲梁进行数值模拟，将模拟结果与课题组所做同尺寸试件试验结果进行对比.图5为混凝土三点弯曲梁荷载－裂缝张口位移(F－CMOD)模拟结果曲线与试验结果曲线的对比，两曲线大体吻合.

图5 混凝土三点弯曲梁F－CMOD 曲线

由图5可知:在初始阶段，荷载与张口位移之间呈线性关系，加载至某一步使混凝土开裂后，荷载与裂缝口张开位移之间呈现出非线性关系;荷载增大到峰值荷载后开始卸载. 试验结果荷载峰值(失稳荷载)略高于模拟值，这主要是由于实际支座与试验台面之间的摩擦阻碍混凝土梁裂缝开裂，因此需要额外的力克服这部分摩擦力.此外，达到失稳荷载后，试验曲线迅速下降，而模拟曲线下降相对较平缓，这是由于实际试验过程中混凝土梁本身的重力加速了裂缝的开裂.

表2为三点弯曲梁模拟及试验计算结果的对比汇总.同样由于试验中重力及摩擦力的存在，使得模拟值与试验值在起裂韧度、失稳韧度等方面存在差异，但总体相差不大.

表2 三点弯曲梁计算结果

图6为混凝土三点弯曲梁最大主应力云图. 由图6可知，加载初期梁顶端中部加载部位为压应力区，裂缝扩展过程中在裂缝尖端区域拉应力较大，存在局部应力集中现象.

图6 混凝土三点弯曲梁最大主应力云图(单位:Pa)

初始裂缝经过10 步扩展在梁的整个高度方向 贯穿，每次扩展均贯穿一整个单元.以初始裂缝尖端相邻单元13683 为例进行分析，图7为单元13683应变与刚度退化标量(SDEG)和分析步关系曲线.在初始阶段SDEG =0，表明混凝土单元仍然完好;当t=0.022 4 时，SDEG 开始增大，此时随着载荷增加混凝土开始发生损伤，微裂缝产生，混凝土单元刚度逐渐退化;到t=0.031 9 时，SDEG 达到1，此时混凝土单元发生开裂，对应应变曲线的峰值，即一旦开裂，混凝土释放能量，应变开始逐渐减小.SDEG 由0到1 经历的时间非常短，即裂纹萌生后不久就发生开裂，这符合混凝土脆性材料的特性. 随着荷载增大，裂缝以单元为单位不断向前扩展.

图7 单元13683 应变与SDEG 和分析步的关系曲线

4 结 语

基于最大主拉应力牵引损伤开裂准则，应用扩展有限元技术对含初始裂缝的C35 混凝土梁进行了三点弯曲数值模拟，依据基于能量的线性软化损伤模型，模拟了混凝土梁在含有初始裂缝情况下从加载到失稳破坏的全过程，获得了裂纹萌生、扩展的过程及破坏形态，分析了裂缝以单元为单位向前扩展的损伤演化过程. 将模拟结果与相关试验进行了对比分析，二者具有较好的吻合性，分析认为存在的差异主要是由于混凝土梁自重、支座摩擦力，此外还可能是材料参数等原因造成的.

［1］李忠献，刘永光.基于虚拟裂缝模型求解混凝土等效断裂韧度的实用解析方法［J］. 工程力学，2006，23(11):91－97.

［2］王学志，张晓刚，宋玉普，等.多裂纹混凝土等效断裂模型研究［J］.水力发电，2006，32(2):20－22，44.

［3］Bazant Z P，Kim J K，Pfeiffer P A.Determination of fracture properties from size effect test［J］.Journal of Structural Engineering，1986，112(2):289－307.

［4］徐世烺，王建敏.水压作用下大坝混凝土裂缝扩展与双K 断裂参数［J］. 土木工程学报，2009，42(2):119－125.

［5］Belytschko T，Black T. Elastic crack growth in finite elements with minimal remeshing［J］.International Journal for Numerical Methods in Engineering，1999，45 (5):601－620.

［6］Melenk J M，Babuska I. The partition of unity finite elementmethod:basic theory and applications［J］. Computer Methods in Applied Mechanics and Engineering，1996，139(14):289－314.

［7］Stolarska M，Chopp D L，Moës N，et al. Modelling crack growth by level sets in the extended finite element method［J］.International Journal for Numerical Methods in Engineering，2001，51(8):943－960.