满足断裂过程区裂纹张开位移条件应力函数的半解析解法*

侯永康1, 段树金1,2, 安蕊梅1

(1. 石家庄铁道大学 土木工程学院, 石家庄 050043; 2. 石家庄铁道大学 道路与铁道工程安全保障省部共建教育部重点实验室, 石家庄 050043)

引 言

1959年Barenblatt[1]提出了“内聚力”模型,认为在裂纹尖端存在一个微小的内聚区,在这个区域内,内聚力和张开位移间存在着σ=f(w)的依存关系作为内聚力模型的一个特例,Dugdale[2]提出了条形塑性区模型,认为塑性区集中在裂纹尖端前缘沿裂纹方向长为ρ、高为t的狭带上,并且假设塑性区内应力均为σ,解除位移约束,通过对应力积分给出了问题的解答,得到了塑性区的长度或裂纹尖端张开位移

随后在混凝土类材料的研究过程中,发现由于混凝土类材料结构的非均质性,其在断裂过程中的变形机制和力学特性与金属材料有着本质的区别,应力与应变的本构关系表现出了软化特性[3-4]1976年,Hillerborg等[5]在B-D模型的基础上提出了虚拟裂纹模型,对混凝土类材料的内聚区特性进行了描述,并引入了“断裂过程区(FPZ)”这个新概念他们认为材料虽已产生受拉损伤但骨料间仍具有“桥联”作用,即在断裂过程区内分布着阻碍裂缝扩展的内聚力,并且内聚力与虚拟裂纹张开位移(COD)之间存在着某种反比关系,即拉应变软化关系,曲线下的面积被定义为断裂能GF,但是没有给出解析解答在此之后,国内外学者专家对混凝土类材料的拉应变软化曲线关系和断裂能的确定进行了深入研究[6-8]

目前,国内外确定混凝土类材料拉应变软化曲线的方法主要有:直接拉伸法、J积分法和逆分析法,分别描述如下

直接拉伸法[9]:通过单轴拉伸试验获得拉应变软化曲线在直拉试验中,丁晓唐等[10]为减小偏心影响,在试件和螺栓传力连接杆之间加设铰接头,进行了混凝土直接拉伸试验,得到了拉应变软化曲线Zhao等[11]采用三峡大坝混凝土的施工配合比制作试块,进行直拉断裂试验,确定了混凝土的软化关系表达式

J积分法[12]:在裂尖区域,约束应力σ是张开位移w的函数,J积分表示为

J=σ(w)dw

(1)

将式(1)对张开位移w求偏微分可得

(2)

将所求出的σ-w曲线去掉弹性变形,只考虑断裂过程区变形,可得出拉应变软化曲线

逆分析法:首先假定拉应变软化关系为某些参数控制的函数表达式,采用数值分析的方法计算特定条件下材料的力学行为,当计算结果与测试结果误差达到最小时,即确定的拉应变软化关系为假定的函数形式赵志方等[13]采用Levenberg-Marquardt优化算法,根据断裂试验获得的荷载-裂缝口张开位移曲线做逆分析,获得了混凝土的软化曲线其在之后的研究中[14],基于改进进化算法的逆推法确定混凝土软化曲线,验证了用逆推法确定软化曲线的可行性冯孝杰[15]采用大坝及其湿筛混凝土楔入劈拉试件的试验数据进行基于进化算法的混凝土软化曲线的逆分析求解,并与直拉法的结果进行了对比研究Su等[16]采用IDCM法测定了荷载-张开位移曲线,由逆分析法推求出拉应变软化曲线

但是在实际应用中,直接拉伸法对试验机的刚度要求很高,且很难保证试验过程中试件与试验机准确的几何对中与力学对中,所以精准测得拉应变软化曲线的难度较高;而J积分法虽然可由简单试验直接求得,但是要精确测得形成裂纹所消耗的能量是非常困难的;逆分析法基于数值计算考虑,通常将软化曲线简化成为线性关系,尽管这在工程应用中对构件开裂、荷载-位移曲线的求解起到了简化的作用,但在某种程度上并不能正确反映断裂过程区内应力、位移的分布及其相互之间的关系,只能从宏观角度对过程区的影响进行实用意义上的评价

Duan等提出的加权积分法[17-20],成功地消除了裂纹尖端应力应变的奇异性,在裂纹尖端构成了内聚力作用和闭合形裂纹相并存的过程区,且二者之间存在确定的反比关系,但是,这些解答依赖于假定的权函数因此,本文以Duan-Nakagawa模型为基础,采用加权积分法,提出了一种求解满足断裂过程区裂纹张开位移条件的应力函数的半解析解法;并以带板对称边裂纹I型问题为例,对此方法进行了应用;根据试验测试数据,反推出权函数,得到满足试验测试数据的应力函数和位移函数,以及相应的拉应变软化曲线和断裂能

1 求解基本思路

在进行试件的断裂试验时,无法直接测得裂纹尖端区域的应力分布,但可以测得断裂过程区的形状、大小和相应的应变分布根据虚裂纹模型的定义,即可得到断裂过程区的长度和张开位移的分布所以此方法以可测得的张开位移作为边界条件,得到相应的应力函数和位移函数下面以对称边裂纹带板为例,说明问题求解的基本思路

1.1 问题的一般解

设一任意含对称边裂纹的无限大板(如图1所示),其中,有限连续区间为2(ak+bk),bk为断裂过程区长度,ak为裂纹尖端到模型中心的距离,σy0k为作用的外荷载,处于Ⅰ型状态(以下简称为裂纹元)采用加权积分法,权函数取ρ(t)=1/bk时,其应力函数和位移函数的一般解为[17]

(3)

(4)

式中,为剪切模量,G=E/[2(1+ν)];对于平面应变问题κ=3-4ν,对于平面应力问题κ=(3-ν)/(1+ν);ν为Poisson(泊松)比其沿x轴的应力和位移如图2所示此解在断裂过程区的内聚应力和张开位移依赖于假定的权函数

图1 对称边裂纹无限大板 图2 沿x轴正应力、张开位移分布
Fig. 1 The infinite plate with double edge cracks Fig. 2 The distribution of normal stress and opening displacement along x-axis

图3 对称边裂纹带板
Fig. 3 The double edge notched plate

1.2 求解满足断裂过程区张开位移条件的特解

现取一含对称边裂纹的带板,受σy0作用,有限连续区间为2(a+b),如图3所示已知在断裂过程区的k点处满足给定的张开位移Uk现用选点法求满足这一条件的特解,这一过程可视为一般解(裂纹元解)的叠加

取带板对称边裂纹模型的一半,如图4所示,将其断裂过程区分为n份,定义n条分别以1,2,3,…,n为断裂过程区尖端的裂纹元,σy01σy02,…,σy0n为各条裂纹元所施加的外荷载各裂纹元在外荷载作用下的内聚力分布如图5所示,分别以σy1σy2,…,σyk,…,σyn表示,将各条裂纹元进行叠加可得

(5)

图4 裂纹元叠加模型
Fig. 4 The superposition model of element cracks

图5 各裂纹元正应力分布
Fig. 5 The normal stress distribution of each element crack

图6为各裂纹元的张开位移分布示意图,图中,uk为第k条裂纹元的张开位移曲线将各裂纹元的张开位移在断裂过程区内叠加,需满足

(6)

式中,uk,j为第k条裂纹元在j点处的张开位移

图6 各裂纹元张开位移分布
Fig. 6 The opening displacement distribution of each element crack

根据方程组(6), 结合式(3)、(4), 得到满足给定裂纹张开位移边界条件的权函数, 再进行加权积分计算, 即可得到相应的应力函数上述思路的关键是求出满足给定张开位移的权函数

2 反推权函数的基本步骤

根据给定裂纹断裂过程区张开位移,求解权函数的主要步骤描述如下:

1) 设加权积分法构造应力函数所采用的权函数为ρ(t)=1/bkbk为裂纹元断裂过程区的长度

2) 给定裂纹断裂过程区长度为b,将其分为n段,沿断裂过程区设n+1个点

3) 首先,选择断裂过程区尖端在n点的裂纹元n,其断裂过程区长度为bn=b,则其所采用的权函数为ρ(t)=1/bn由图6可得,在n-1点处给定裂纹的张开位移即为裂纹元nn-1点处的张开位移,即Un-1=un,n-1而裂纹位移已经给定,将Un-1代入式(4),可得裂纹元n所需的外加荷载σy0n通过外加荷载σy0n,推出裂纹元n在断裂过程区其他各点的张开位移,即un,n-2un,n-3,…,un,1un,0

4) 同理,断裂过程区尖端在n-1点的裂纹元n-1,其断裂过程区长度为bn-1,则其所采用的权函数为ρ(t)=1/bn-1由图6可知,在n-2点处给定裂纹的张开位移为Un-2=un,n-2+un-1,n-2,则由给定的Un-2和上一步计算所得的un,n-2可推算出un-1,n-2,即裂纹元n-1在n-2点处的张开位移un-1,n-2代入式(4),可得裂纹元n-1所需的外加荷载σy0(n-1),通过外加荷载σy0(n-1)推出裂纹元n-1在断裂过程区其他各点的张开位移,即un-1,n-3un-1,n-4un-1,1un-1,0

5) 依次类推

6) 断裂过程区尖端在1点的裂纹元1,其断裂过程区长度为b1,则其所采用的权函数为ρ(t)=1/b1由图6可得,在0点处给定裂纹的张开位移U0=u1,0+u2,0+u3,0+…+un-1,0+un,0,则由给定的U0和上述步骤依次推出的u2,0u3,0,…,un-2,0un,0,可推算出u1,0所以将u1,0代入式(4),可得裂纹元1所需的外加荷载σy01

7) 最后得到分别以1,2,3,…,n-1,n为断裂过程区尖端裂纹元的外加荷载将各个裂纹元的外荷载叠加,可得满足给定位移边界条件时,各裂纹元断裂过程区尖端所需荷载的值对荷载进行拟合,得到的曲线方程即为满足断裂过程区各点张开位移的权函数,问题得解

3 算 例

以图3所示带板对称边裂纹Ⅰ型问题为例(取模型中的右半部分),对此方法进行应用

设模型中的基本参数为a=40 mm,b=10 mm,初始裂纹长度为10 mm,弹性模量E=2.1×104MPa,抗拉强度ft=5 MPa,Poisson比ν=0.1,外荷载为σy0

根据表1所给定的断裂过程区张开位移Uk,按照第2节反推权函数的求解步骤,取裂纹元断裂过程区的长度分别为:10,9.5,9,8.5,8,7,6,4,2 mm,计算可得各裂纹元相应的荷载,对所求荷载进行叠加并拟合,从而得出满足张开位移条件的权函数如式(7)所示,拟合曲线如图7所示

ρ(t)=0.021 4t-0.842 9 (40≤t≤50)

(7)

图7 权函数拟合
Fig. 7 The fitting weighted function

表1断裂过程区张开位移分布

Table 1 The distribution of COD in FPZ

x/mm4040.54141.5424344464850(dCOD/2)/mm00.001 10.001 60.003 80.005 30.011 20.019 40.048 80.087 70.139

采用所求得的权函数ρ(t),进行加权积分计算,得到积分后的应力函数解析表达式,由此得到相应的正应力分量和位移分量解析表达式,分别如式(8)和式(9)所示:

,

(8)

(9)

式(9)可积,但由于公式较长,积分结果予以省略

图8为根据式(8)所得到的正应力分量分布图(由对称性,取模型的1/4)

图8 正应力分布 图9 位移分布 Fig. 8 The normal stress distribution Fig. 9 The displacement distribution

图9为根据式(9)得到的裂纹张开位移分布图(模型的1/4)

裂纹延长线上断裂过程区内的张开位移和内聚力分布如图10所示从而得出相应的拉应变软化曲线(如图11所示),而拉应变软化曲线与xy轴所围成面积的大小即为断裂能GFGF=297.238 N/m

图11中也分别给出当拉应变软化曲线类型选直线型[5]和折线型[21]时的情况由图可得,本文方法获得的裂纹最大发展宽度为0.137 4 mm,而直线型和折线型获得的裂纹最大发展宽度分别为0.129 mm、0.141 mm经比较,与直线型和折线形的裂纹最大发展宽度相差较小,说明该方法具有较高的精度

图10 断裂过程区内聚力和张开位移
Fig. 10 The cohesive stress and COD in FPZ
图11 拉应变软化曲线
Fig. 11 The tensile strain softening curves

4 结 论

本文提出了一种基于Duan-Nakagawa模型,根据给定的断裂过程区张开位移,推求满足张开位移条件的应力函数的半解析解法

1) 本方法可以导出满足断裂过程区裂纹张开位移条件的应力函数和位移函数的全场解,其中的关键是权函数的确定

2) 由此得到了一种确定材料拉应变软化曲线和断裂能的新方法, 其理论基础坚实, 简单易行

3) 由于断裂过程区尖端张开位移变化趋于剧烈,所以选点时在断裂过程区尖端附近应布置得相对密集些

4) 采用Dugdale方式只能给出断裂过程区的内聚力为恒定或线性这样一些简单分布时的解,而本方法可以得到满足断裂过程区张开位移条件的任意特解

今后的研究工作可以进行各种含切口试件的断裂破坏试验,根据断裂过程区的变形特征,用半解析解法得到相应的位移函数和应力函数,并验证拉应变软化曲线和断裂能作为材料断裂性能参数的有效性

参考文献(References):

[1] BARENBLATT G I. The formation of equilibrium cracks during brittle fracture, general ideas and hypotheses, axially-symmetric cracks[J].Journal of Applied Mathematics&Mechanics, 1959,23(3): 434-444.

[2] DUGDALE D S. Yielding of steel sheets containing slits[J].Journal of the Mechanics and Physics of Solids, 1960,8(2): 100-104.

[3] RUSCH H, HILSDORF H. Deformation characteristics of concrete under axial tension[R]. Vorunterschungen 44, Munich, Bericht, 1963.

[4] EVANS R H, MARATHE M S. Microcracking and stress-strain curves for concrete in tension[J].Material and Structures, 1968,1(1): 61-64.

[5] HILLERBORG A, MODEER M, PETERSON P. Analysis of crack formation and crack growth in concrete by means of fracture mechanics and finite elements[J].Cement and Concrete Research, 1976,6(6): 773-782.

[6] ELICES M, GUINEA G V, GMEZ J, et al. The cohesive zone model: advantages, limitations and challenges[J].Engineering Fracture Mechanics, 2002,69(2): 137-163.

[7] 卿龙邦, 李庆斌, 管俊峰, 等. 基于虚拟裂缝模型的混凝土断裂过程区研究[J]. 工程力学, 2012,29(9): 112-116.(QING Longbang, LI Qingbin, GUAN Junfeng, et al. Study of concrete fracture process zone based on fictitious crack model[J].Engineering Mechanics, 2012,29(9): 112-116.(in Chinese))

[8] BAZANT Z P. Concrete fracture models: testing and practice[J].Engineering Fracture Mechanics, 2002,69(2): 165-205.

[9] LI Q B, ANSARI F. High-strength concrete in uniaxial tension[J].ACI Materials Journal, 2000,97(1): 49-57.

[10] 丁晓唐, 丁鑫, 刘海霞, 等. 混凝土直拉试验和三点弯曲断裂试验确定的软化曲线的比较[J]. 水电能源科学, 2014,32(1): 116-118.(DING Xiaotang, DING Xin, LIU Haixia, et al. Comparison study of softening curves of concrete by direct tension test and three-point bending fracture test[J].Water Resources and Power, 2014,32(1): 116-118.(in Chinese))

[11] ZHAO Z, ZHANG J, ZHOU H, et al. Two methods for determining softening relationship of dam concrete and wet-screened concrete[J].Advances in Structural Engineering, 2016,15(15): 1125-1138.

[12] RICE J R. A path independent integral and the approximate analysis of strain concentration by notches and cracks[J].Journal of Applied Mechanics, 1968,35(2): 379-386.

[13] 赵志方, 李铭, 赵志刚. 逆推混凝土软化曲线及其断裂能的研究[J]. 混凝土, 2010(7): 4-7.(ZHAO Zhifang, LI Ming, ZHAO Zhigang. Research on reversing softening curve and fracture energy of concrete[J].Concrete, 2010(7): 4-7.(in Chinese))

[14] 赵志方, 王刚, 周厚贵, 等. 混凝土拉伸软化曲线确定方法的对比研究[J]. 浙江工业大学学报, 2015,43(4): 455-459.(ZHAO Zhifang, WANG Gang, ZHOU Hougui, et al. A comparative study of the methods for determining the tensile softening curve of concrete[J].Journal of Zhejiang University of Technology, 2015,43(4): 455-459.(in Chinese))

[15] 冯孝杰. 大体积混凝土软化曲线的新确定方法[D]. 硕士学位论文. 杭州: 浙江工业大学, 2012.(FENG Xiaojie. A new method of determining softening curve of mass concrete[D]. Master Thesis. Hangzhou: Zhejiang University of Technology, 2012.(in Chinese))

[16] SU R K L, CHEN H H N, KWAN A K H. Incremental displacement collocation method for the evaluation of tension softening curve of mortar[J].Engineering Fracture Mechanics, 2012,88: 49-62.

[17] DUAN S J, NAKAGAWA K. Stress functions with finite stress concentration at the crack tips for a central cracked panel[J].Engineering Fracture Mechanics, 1988,29(5): 517-526.

[18] ZHU M, CHANG W V. An unsymmetrical fracture process zone model and its application to the problem of a radical crack with an inclusion in longitudinal shear deformation[C]//Proceedings of FRAMCOS-3/Fracture Mechanics of Concrete Structures. Freiburg, Germany, 1997: 1097-1106.

[19] 段树金, 前田春和, 藤井康寿, 等. 沿直线有多条裂纹的薄板弯曲问题[J]. 工程力学, 1999,16(3): 21-29.(DUAN Shujin, MAEDA H, FUJII K, et al. The problem of bending of a thin plate with a number of cracks along the line[J].Engineering Mechanics, 1999,16(3): 21-29.(in Chinese))

[20] DUAN S J, NAKAGAWA K. A mathematical approach of fracture macromechanics for strain-softening material[J].Engineering Fracture Mechanics, 1989,34(5): 1175-1182.

[21] ROELFSTRA P E, WITTMANN F H. Numerical method to link strain softening with failure of concrete[C]//Fracture Toughness and Fracture Energy in Concrete. Amsterdam, Netherlands, 1986: 163-175.

A Semi-Analytical Method for Stress Functions Meeting Crack Opening Displacements in Fracture Process Zones

HOU Yongkang1, DUAN Shujin1,2, AN Ruimei1

(1.School of Civil Engineering,Shijiazhuang Tiedao University,Shijiazhuang050043,P.R.China; 2.Key Laboratory of Roads and Railway Engineering Safety Control of Ministry of Education,Shijiazhuang Tiedao University,Shijiazhuang050043,P.R.China)

Abstract:Based on the Duan-Nakagawa model, with the weighted integral method, a semi-analytical method for stress functions meeting crack opening displacements in fracture process zones was proposed. The weighted function was determined by means of the boundary selected point method and the superposition of analytical functions with the same crack length but different fracture process zone lengths, to meet the given crack opening displacement in the fracture process zone, and then the final stress function and displacement function can be obtained with the weighted integral method. As an example, a special analytical solution for a double edge notched plate under Mode-I loading was derived, and the tensile strain softening curve and the fracture energy were obtained.

Key words:fracture mechanics; fracture process zone; semi-analytical method; stress function; crack opening displacement; tensile strain softening curve; double edge notched plate

引用本文/Cite this paper: 侯永康, 段树金, 安蕊梅. 满足断裂过程区裂纹张开位移条件应力函数的半解析解法[J]. 应用数学和力学, 2018,39(8): 979-988.HOU Yongkang, DUAN Shujin, AN Ruimei. A semi-analytical method for stress functions meeting crack opening displacements in fracture process zones[J].Applied Mathematics and Mechanics, 2018,39(8): 979-988.

文章编号1000-0887(2018)08-0979-10

ⓒ 应用数学和力学编委会,ISSN 1000-0887

*收稿日期2017-11-23;

修订日期:2018-01-15

基金项目河北省自然科学基金 (A2015210029);河北省教育厅青年基金项目(QN2014062);河北省研究生创新资助项目(CXZZBS2017132)

作者简介侯永康(1990—),男,博士生(通讯作者. E-mail: 15233255808@163.com).

摘要基于Duan-Nakagawa模型,采用加权积分法,提出了一种满足断裂过程区裂纹张开位移条件应力函数的半解析解法该方法结合边界选点法,通过叠加含有相同裂纹长度但断裂过程区长度不同的解析函数,得到满足给定裂纹张开位移的权函数,再进行加权积分得到相应的应力函数和位移函数以带板对称边裂纹I型问题为例,应用上述方法成功导出了特定的应力函数和位移函数,以及相应的拉应变软化曲线和断裂能

断裂力学; 断裂过程区; 半解析解法; 应力函数; 裂纹张开位移; 拉应变软化曲线; 对称边裂纹

中图分类号TU528; O346

文献标志码:A

DOI:10.21656/1000-0887.380296