0.   引言

大量的断裂事故表明,构件断裂都是由于其内存在的各种类型裂纹,使得结构的承载力在一定程度上被削弱,从而影响工程结构的质量和安全。因此,研究裂纹扩展规律,对工程设计、施工具有巨大的指导意义。基于以上目的,很多学者采用数值方法模拟裂纹扩展,并取得大量成果。有限元法是最常用的断裂分析的数值方法之一,但由于裂纹不连续位移场的特殊性及裂尖端应力、应变场的奇异性问题,使得对有限元网格的划分有特殊要求,如宏观裂纹表面要保证与单元边界一致,并依靠有限单元边界来勾画其位置和形状；为描述裂尖应力奇异区,甚至需要在裂尖区域细化单元或改用特殊形式的缝端奇异单元；随着裂缝扩展,有限元网格需要不断地调整和重新划分,产生极大的前处理工作量,不利于具有复杂几何结构的构件的裂纹分析。

为了解决上述困难,在传统有限元基础上,涌现出了许多新的适合于不连续分析、特别是断裂分析的数值方法,如数值流形方法、广义有限元法、无单元法、扩展有限元法（XFEM）等。XFEM允许裂纹沿着单元内部的任何位置起裂,而不仅仅在单元的边界起裂,因此,扩展有限元已经成功解决了许多裂纹问题。然而,当位移不连续时,扩展有限元仍然需要引入外部准则,而且在解决三维裂纹问题时遇到了一系列困难。边界元方法缺点是系数阵为非对称满阵,难以精确模拟裂尖应力场,对三维问题这个缺陷尤为突出。DEM和DDA不能求解块体内部的应力应变,裂纹只能沿块体边界扩展。由石根华博士提出的数值流形法（NMM）,为复杂裂纹问题的分析提供了一种全新的思路和方法。在NMM中,定义了两套独立的数学网格和物理网格,来统一求解连续和不连续问题。因其数学覆盖系统独立于物理网格定义,对各类连续和不连续问题的分析,NMM的前处理需求可以十分简单,特别适于求解各类裂纹问题。张国新等^[1]首次应用NMM模拟裂纹扩展,通过模拟坝基断裂过程展现了NMM处理大变形,大位移的能力。至今为止,NMM已经有了很多进展^[2-3],从基本及理论^[4]到工程应用^[5]。但鲜有研究展示数值流形法同时模拟室内拉伸、剪切试验和工程尺度的断裂模块。所以本文结合改进后的断裂准则,采用C语言编程开发了用于模拟裂纹扩展的数值流形法程序,通过模拟多个算例,验证了程序的准确性和适应性。

1.   NMM基本原理

NMM最大的特点在于它包含两个覆盖系统,数学覆盖和物理覆盖系统。物理覆盖是由数学覆盖和物理域（裂纹,边界等不连续区域）叠加构成的。流形单元是由物理覆盖相互叠加形成的区域。上述过程如图1所示,图1（a）中M₁,M₂和M₃为数学片,其分别与物理域叠加后形成如图1（b）所示的物理片P₁₁,P₁₂,P₂₁、P₂₂,P₃物理片；这5个物理片叠加分别形成流形单元E₁和E₂。

图  1  流形单元形成过程

Figure  1.  Formation process of manifold elements

下载: 全尺寸图片 幻灯片

NMM的物理覆盖被几何或物理的边界分开,因此即使在非常复杂的不连续条件下,NMM依然可以自然的对不连续体进行建模。NMM的加权函数特征取决于数学覆盖的形状和位移场的阶次。全局位移函数是多个物理覆盖的公共部分上局部独立覆盖函数的加权平均值。对于流形单元e,有3个物理片P₁,P₂和P₃,则流形单元e的位移可表示为

式中,权函数w_ei(x,y)是形函数,其要求与有限元中的形函数一致,满足刚体平移,轴对称和变形协调等。NMM采用单纯形积分计算刚度矩阵和荷载矩阵。流形单元可以是任意形状的区域,甚至是凹多边形,这不同于有限元中要求单元形状必须是凸多边形。所以在数值流形法中,前处理只需要选择覆盖的形状和尺寸,就可以生成计算所需的网格,不要求使用者考虑边界和网格布置。

2.   裂纹扩展准则

经典的断裂力学中的断裂准则多适用于裂纹较少的材料扩展问题中。对于复杂裂纹扩展,且发生裂纹相互贯穿时,断裂准则很难适用^[6]。为此,很多学者在裂纹扩展中采用强度准则作为裂纹是否扩展的依据,其中莫尔–库仑准则更是广为应用^[7]。其主要原因是：①在很多工程中,M–C完全可以满足精度要求,并且反应材料性能^[8]；②更重要的特点是,M–C只需要3个材料常数,这些参数可以通过简单室内试验测试得到^[9]；③表达式简洁,便于工程人员使用。

在岩土力学中,以压应力为正,考虑一个与小主应力σ₃方向夹角为β的斜面,如图2所示。

图  2  莫尔-库仑准则

Figure  2.  Mohr–Coulomb criterion

下载: 全尺寸图片 幻灯片

根据莫尔–库仑准则,该斜面上的剪切强度为

式中,$c$为材料的黏聚力,$\phi$为内摩擦角,$\sigma$为法向应力。斜面上的法向应力计算式为

切向应力计算式为

存在一个临界平面,在该平面上剪应力最先达到剪切强度。如图3所示,该临界面与小主应力${\sigma }_3$方向夹角为

图  3  拉伸和剪切破坏示意图

Figure  3.  Schematic diagram of tensile and shear failures

下载: 全尺寸图片 幻灯片

当如图2所示的莫尔–库仑包络线外推到σ₁=0时,与σ₃坐标轴存在一个交点,该交点对应材料的单轴拉伸强度。然而试验观察到的材料的单轴拉伸强度通常低于这个值。因此需要将莫尔–库仑包络线截断（如图2所示）以反映材料发生拉伸破坏的情况。

式中,$T_0$表示材料的抗拉强度。

按照断裂力学计算,裂尖点处的应力是无穷大的,这在实际问题中是不可能存在的。即无法具体求得裂纹尖端处的应力。无法直接使用莫尔–库仑准则确定裂纹走向,为此,An等^[10]采用一种平均化的应力来得到裂尖点的应力,裂尖应力可以表示如下：

式中,$A_k$为单元$k$的面积,m是平均的单元的数量。

在裂纹扩展过程中,可能会出现拉伸断裂或者剪切断裂。拉伸断裂的裂纹扩展方向只有两个,但是剪切断裂的裂纹方向有4个,如图3所示。

图3（a）所示的是拉伸破坏,理论上,拉伸破坏面与小主应力σ₃的夹角是90°,从裂尖点起裂,有两个方向。一般情况下取图中1的方向,因为方向1与裂纹的夹角比较小^[10]。图3（b）所示的是剪切破坏,剪切破坏裂纹面与σ₃的夹角分别是±(π/4+φ/2),但是从裂尖点起裂,裂纹面分别有4个方向。一般情况下,与原始裂纹夹角最小的扩展方向是最终的扩展方向^[10],即图3（b）中的1方向。本文对其改进,加入最大应力准则。即,判断所有上述图3中可能的裂纹扩展方向上的应力情况,按照拉应力最小的方向进行扩展。

3.   裂纹扩展实例

3.1   NSCB试验

NSCB试验是测试脆性材料拉伸强度的有效手段之一,尤其是测试岩石材料的拉伸强度。除了采用室内试验,一些学者也采用连续的数值方法^[10]和非连续数值方法^[11]来研究NSCB拉伸试验。AYATOLLAHI等^[12]对不同角度的$\beta$进行试验,得到不同角度下的拉伸试验结果,并且采用网格重画分的方法模拟了裂纹扩展路径。随后XIE等^[13]基于应力强度因子准则,采用XFEM对AYATOLLAHI等^[12]的部分试验进行重新模拟。网格重画分模拟裂纹扩展时,计算比较耗时,而且网格继承时会损失计算精度^[14]。采用XFEM计算裂纹扩展时,容易产生病态矩阵,导致计算不容易收敛^[15]。

本文在NMM方法中,加入修正的M–C强度准则来模拟SBC试验的裂纹扩展。根据AYATOLLAHI等^[12]的试验,计算模型图4所示,其中$\beta$分别为0°、10°、20°、30°、40°、43°、47°和50°,并且S/R=0.72。数值模型是基于改进后的M–C强度准则的2D-NMM模型,材料各向同性。

图  4  半圆盘弯曲拉伸试验模型

Figure  4.  Bending tensile tests of half disc

下载: 全尺寸图片 幻灯片

图5为半圆盘弯曲拉伸试验和的NMM数值模拟结果对比。当β=0°时,为纯Ⅰ型拉伸问题,裂纹从裂尖开始沿着拉力最大的中轴线上扩展,直至破坏。随着β的角度的增加,模型变成Ⅱ型占优势的裂纹扩展问题。在微观尺度上,在向试样中心线的曲线方向上,产生微裂纹,随着荷载增加,裂纹向沿着曲线向试样中心线处扩展。由图5（b）中展示了数值模拟得到的试样裂纹扩展结果,对比发现,各个不同初始角度的裂纹扩展路径与试验结果图5（a）吻合的非常好。对于应力强度因子以及不同长度的初始裂纹对结果的影响,XIE等^[13]已经采用XFEM方法详细讨论,这里不再赘述。

图  5  半圆盘弯曲拉伸试验结果和NMM模拟对比

Figure  5.  Comparison of bending tensile test results of half disc and NMM simulation

下载: 全尺寸图片 幻灯片

3.2   四点双边剪切梁

四点双边剪切梁（double-edge notched（DEN）concrete beams loaded in four-point-shear）是混凝土结构测试中的重要试验之一,早在90年代BOCCA等^[16]对四点双边剪切梁进行试验并采用数值方法模拟裂纹开裂过程；随后GEERS^[17]采用损伤模型对四点双边剪切梁试验进行模拟；ZHU和TANG^[18]基于有限元微观尺度对的四点双边剪切梁进行模拟。本文根据BOCCA等^[16]的试验模型,NMM模拟尺寸如图6所示。

图  6  四点双边剪切梁模型

Figure  6.  Four-point bilateral shear beam model

下载: 全尺寸图片 幻灯片

假设材料为2D均匀材料,其单轴压缩强度和弹性模量分别为37.9 MPa和27.9 GPa。如图6所示,试样的长度为400 mm,宽为150 mm,预制裂纹的宽度为5 mm,长度为25 mm。在数值模型中,裂纹宽度为0 mm。为了使得在数值模拟过程中整个模型是平衡的,顶部两个荷载的比值应该是15∶1。数值结果如图7所示。

图  7  NMM模拟结果与试验结果对比

Figure  7.  Comparison between NMM simulation and test results

下载: 全尺寸图片 幻灯片

图7四点双边弯曲剪切梁流形网格图以及数值模拟结果与试验结果对比,试验结果来源于BOCCA等^[16]的试验。由于模型的对称性,破坏结果称中心对称形状。在荷载的作用下,预制的裂纹之间存在很大的拉应力,拉应力的分布在裂尖附近比较大,所以裂纹从裂尖开始扩展。由图7可知,模拟得到的裂纹扩展路径与试验高度一致,并且不需要采用复杂的网格重画分加密网格技术^[19],只需要少量的网格就可以得到比较好的结果。Zhu等^[18]采用基于损伤的、非均匀的有限元模拟的四点弯曲剪切梁结果,详细分析了该模型的受力情况,本文不再赘述。通过本例结果可以看出NMM在模拟剪切裂纹扩展时也同样适用。

3.3   重力坝开裂

在这部分内容中,将研究全尺寸的重力坝模型,见图8,这个模型已经被研究过很多次^[20-21]。在本模型中将在NMM采用修改的MC的强度准则来模拟重力坝在静态水荷载作用下的坝体开裂情况。假设水坝的材料本构符合线弹脆性本构,并且在距离水坝顶端36.5 m处有一个初始裂纹。水坝初始受重力作用,并且假设水体刚好到达坝顶时坝体此时的受力状态。由于受坝体已经漫顶流作用,水坝左边还受均匀的漫顶压力作用。水坝假设为均匀材料,弹性模量为250 GPa,泊松比为0.2,拉伸强度为10 MPa。水坝底边为固定边界,其他边界为自由边,结果如图8所示。

图  8  重力坝模型

Figure  8.  Gravity dam model

下载: 全尺寸图片 幻灯片

图9展示了NMM模拟重力坝不同时间步下网格图和裂纹扩展路径。整个坝体的破坏过程以拉伸破坏为主,与DIAS等^[20]中的结论一致,对于本算例,最重要的是验证不同的数值算法和裂纹扩展准则,而裂纹扩展路径是否一致是判断的重要标准。图10是NMM和FEM的模拟结果对比。

图  9  重力坝开裂过程

Figure  9.  Cracking process of gravity dam

下载: 全尺寸图片 幻灯片

图  10  NMM与FEM结果对比

Figure  10.  Comparison between NMM and FEM results

下载: 全尺寸图片 幻灯片

图10中,网格较粗的是本文采用NMM计算的网格及扩展路径,网格较细的结果是采用FEM^[20]的计算结果。对比发现,二者的结果高度一致。说明了本文的计算结果是正确的。但是本文的计算方法采用的网格较粗,且在裂纹扩展路径上形成结构非连续面,而不是FEM^[20]中得到的损伤路径。通过模拟重力坝开裂,验证了本文基于NMM提出的裂纹扩展算法模拟这类拉伸问题是准确的、可靠的。

3.4   边坡滑移

Camones等^[22]采用离散单元法模拟了边坡滑移过程,由于离散元模拟过程中需要采用大量粒子计算,其效率较低^[23]。由前面的算例可以知道,NMM在网格较粗的时候,对一些简单的拉伸问题和剪切问题都有较好的适应性。对于多裂纹扩展问题,例如边坡问题,可能同时包含多个拉伸和剪切裂纹,NMM算法也有较好的适应性^[10,24]。但文献[10, 24]中采用M–C强度准则,裂纹扩展角的方向是与原始裂纹面夹角最小的方向。本文采用前面修改后的M–C准则模拟边坡多裂纹扩展,计算模型见图11。

图  11  边坡尺寸

Figure  11.  Slope dimensions

下载: 全尺寸图片 幻灯片

图12展示了边坡在重力作用下,边坡内部发生破坏,裂纹扩展切割形成多个小块体。NMM的模拟结果与DEM的模拟结果相类似,说明了在初始阶段上,边坡上的裂纹扩展计算是没有问题。图13表示各个时步下整个边坡系统的流形单元和边界状况。在初始阶段,边坡受重力作用,裂纹主要受张拉作用扩展。随后形成的各个块体沿着斜坡向下加速滑移,最后冲击底边。本文的结果与文献[10,22,24]的结果一致,说明本文修正的M–C准则也可以应用于边坡复杂裂纹扩展问题。本算例为NMM计算边坡滑移及边坡滑移灾害评估打下基础。

图  12  边坡破坏结果

Figure  12.  Results of slope failure

下载: 全尺寸图片 幻灯片

图  13  边坡破坏滑移过程

Figure  13.  Failure and sliding process of slope

下载: 全尺寸图片 幻灯片

4.   结论

本文回顾了NMM的基本概念,通过改进莫尔–库仑准则和最大周向应力准则,作为新的裂纹扩展准则。采用提出的裂纹扩展准则编制C语言程序,通过模拟半圆盘弯曲拉伸试验和四点双边剪切梁试验,确认了提出的算法可以模拟拉伸和剪切破坏时的裂纹扩展,并且结果是可靠、可信的。随后,模拟了重力坝开裂问题,在流形元单元尺寸远大于有限元网格时,NMM’依然可以得到与有限元模拟结果一致的结果,说明NMM具有很好的适用性和计算精度。最后通过模拟边坡滑移问题,计算结果与DEM模拟结果具有高度一致性,说明所提出的裂纹扩展准则在模拟边坡裂纹扩展问题也有很好的适应性。本文结果不仅验证了提出的强度准则在模拟裂纹扩展问题是有效的,也为模拟工程中的裂纹扩展问题提供新工具。