复合裂纹的应力强度因子有限元分析

徐  慧^1,2，伍晓赞¹，程仕平²，李燕峰^1,2，龙朝辉³，邓超生¹

(1. 中南大学 物理科学与技术学院，湖南　长沙，410083；

2. 中南大学 材料科学与工程学院，湖南　长沙，410083；

3. 湖南师范大学 物理与信息科学学院，湖南　长沙，410081)

摘要：讨论裂纹尖端的应力应变与应力强度因子的关系, 建立计算复合型裂纹应力强度因子的有限元方法，应用有限元分析软件ANSYS计算I型裂纹和II型裂纹的应力强度因子以及裂纹尺寸和载荷对应力强度因子的影响。研究结果表明：ANSYS解与解析解很接近，误差很小，验证了复合裂应力强度因子计算方法的有效性；ANSYS解在裂纹较大和较小时误差相对较大，这主要与划分网格的精度有关，裂纹较大时网格不够密，裂纹较小时网格产生变形的影响，可以通过增加网格精度来减小计算误差。

关键词：复合型裂纹；有限元法；ANSYS软件；应力强度因子

中图分类号：O346.1；O242.21         文献标识码：A         文章编号：1672-7207(2007)01-0079-05

Finite element analysis of stress intensity factor in composite mode crack

XU Hui^1,2, WU Xiao-zan¹, CHENG Shi-ping², LI Yan-feng^1,2 , LONG Zhao-hui³, DENG Chao-sheng¹

(1. School of Physical Science and Technology, Central South University, Changsha 410083, China;

2. School of Materials Science and Engineering, Central South University, Changsha 410083, China;

3. School of Physics and Information Science, Hunan Normal University, Changsha 410081, China)

Abstract: The releationship between stress and stress intensity factor（SIF）was studied. Using the software ANSYS, the method of calculating SIF of composite mode cracks was set up. At the same time, the SIF of mode I and II crack was calculated and the influence of crack size and loading on the SIF was studied. The results show that the results of ANSYS is almost equal to analytical results, which proves that the method is feasible. When the crack is very big or very small, the result of ANSYS’ error is comparatively obvious, which is due to the problem of the precision of gridding. When the crack is too big, the dense of gridding is inadequate; when the crack is too small, the cell of gridding will produce distortion. The error can be decreased by increasing the precision of gridding.

Key words: composite mode cracks; finite element method; ANSYS software stress intensity factor

在断裂力学的工程应用中, 应力强度因子是一个判断含裂纹构件的断裂或计算裂纹扩展速率的重要参量。当前已有许多种计算应力强度因子的计算方法, 如解析法^[1]、边界配位法^[2-5]、有限单元法^[6-8]、边界元素法、体力法、权函数法^[9-12]、线弹簧模型和复变函数法^[13]，研究的对象有椭圆裂纹^[14]、片状裂纹、孔边裂纹等。由于有限单元法计算结果相对精确，所以，对于结构或裂纹形状复杂和受复杂载荷作用的结构件，一般采用有限单元法。目前，人们应用有限元方法计算I型裂纹的应力强度因子的研究较多^[15-17]，而对计算II型、III型裂纹的应力强度因子的研究较少。在此，本文作者在I型裂纹的基础上，扩展有限元软件ANSYS计算裂纹断裂参量的方法。

1　基本理论和方法

1.1　应力强度因子的求解方法

在进行裂纹分析时,应力强度因子是一个重要的指标。构件中的裂纹可以抽象成二维或三维模型。二维和三维裂纹模型的裂纹尖端的坐标分别如图1和图2所示。

图1　二维裂纹前缘任意点的局部坐标(法平面)和位移场

Fig.1　Displacement fields and local polar coordinates

at a arbitrary point in the 2-D crack front

图2　三维裂纹前缘任意点的局部坐标(法平面)和位移场

Fig.2　Displacement fields and local polar coordinates

at a arbitrary point in the 3-D crack front

对于三维裂纹模型，裂纹尖端附近在图2所示的局部坐标系(局部直角坐标系由坐标分量X，Y和Z构成, XY 面为裂纹前缘任一点的法平面, YZ 平面为其切平面。局部柱坐标系的坐标分量由r，θ和Z组成，r和θ所构成的平面为裂纹前缘任一点的法平面)下的位移场可表示如下：

1.2　复合裂纹理论

复合裂纹因其自身结构的复杂性，从而导致其受力的复杂性，因此，在计算中，不仅要考虑裂纹的结构，还要考虑其复杂的受力。

本文在计算复合裂纹的应力强度因子时采用的方法是将复杂的裂纹分解成Ⅰ，Ⅱ和Ⅲ型分量，然后进行分别计算。以试样的厚、宽、长3个方向取直角坐标XOYZ(见图3)，空间斜裂纹在试样正面(XOY)上与OY的夹角为α，在试样侧面(XOZ)上与OX的夹角为β，试样受铅垂力P的作用，P平行于Z轴，空间斜裂纹的裂纹面为平行四边形，外力P在裂纹面上可分解为 P_Ⅰ，P_Ⅱ和P_Ⅲ 3个分力。其中：P_Ⅰ为垂直于裂纹面的拉力，产生Ⅰ型应力强度因子，P_Ⅱ和P_Ⅲ在裂纹面内；P_Ⅱ为垂直于裂纹前端的剪力，产生Ⅱ型裂纹应力强度因子；P_Ⅲ为平行于裂纹前端的剪力，产生Ⅲ型裂纹应力强度因子。坐标原点O（0，0，0）取在裂纹面上。在裂纹面上还可取A(-1，0，tan β)和B(0，-1，tan β) 2点，由此可得裂纹面内的矢量OA={-1，0，tan β}和OB={0，-1，tan α}，由于裂纹前端与裂纹的切口线平行，所以P_Ⅲ//OA={-1，0，tan β}。P_Ⅰ垂直于OA和OB，故P_Ⅰ//(OA×OB)={tan β，tan α, 1}，P_Ⅱ//[(OAOB) OA]={tan α tan β，-(1+tan² β), tan ?}。载荷矢量P={0，0，P}。将P在P_Ⅰ，P_Ⅱ和P_Ⅲ 3个方向投影，可求得3个分力的大小：

把力分解成P_Ⅰ，P_Ⅱ和P_Ⅲ 3个分量之后，就可以对每1个分量求解其对应的应力强度因子。K_Ⅰ，K_Ⅱ和K_Ⅲ总的裂纹应力强度因子为：

图3　复合裂纹的位置

Fig.3　Displacement of composite mode cracks

2　应力强度因子的有限元计算

2.1　断裂模型及网格

在建模时，取平板的宽度b=0.2 m，裂纹长度为2a，材料的弹性模量为2.23×10¹¹ Pa，泊松比为0.36。考虑到裂纹的对称性，有限宽板的建模只考虑1/4的模型即可，网格如图4所示。

2.2　Ⅰ型裂纹

Ⅰ型裂纹又称张开型裂纹，其载荷垂直裂纹面，取裂纹的宽度a=0.02 m，单向拉伸应力为100 MPa，在对称边界上加上对称边界，对裂纹尖端进行奇异化处理，然后，应用有限元软件对模型进行结构分析，得到其应力云图。可以看出，裂纹尖端及其附近的应力有明显的应力集中。通过有限元软件ANSYS计算断裂参量的功能，算出其应力强度因子为K_Ⅰ=18.04 MPa?m^1/2，据应力强度因子手册可知，此时K_Ⅰ 通过计算得到解析解为K_Ⅰ=17.84 MPa?m^1/2，采用ANSYS所得计算结果与解析解相对误差为1.12%。

图4　断裂模型

Fig.4　Model of crack

为了研究Ⅰ型裂纹的扩展，计算了裂纹的尺寸和所受载荷变化时其裂纹应力强度因子的变化，并将其与解析解进行对比，求出其相对误差，结果见表1和  表2。

2.3　Ⅱ型裂纹

Ⅱ型裂纹又称滑开型裂纹。在计算Ⅱ型裂纹时，裂纹模型加载和材料的属性与Ⅰ型裂纹的一样，所加剪切应力为100 MPa，可以应用表面效应单元SURF153加载剪力，应用ANSYS有限元软件对模型进行结构分析，得到其应力云图。可以看出，在裂纹尖端同样有明显的应力集中，通过有限元软件ANSYS计算断裂参量的功能，算出其应力强度因子为K_Ⅱ=17.81 MPa?m^1/2，据应力强度因子手册，此时K_Ⅱ=Cσ(πa)^1/2，通过计算得到解析解为K_Ⅱ=17.84 MPa?m^1/2，相对误差为-0.17%。

为了研究Ⅱ型裂纹的扩展，本文计算了Ⅱ型裂纹的尺寸和所受剪切应力变化时其裂纹应力强度因子的变化，并将其与解析解进行对比，求出其相对误差，结果见表3和表4。

从表3和表4可以看出，ANSYS的计算结果很准确，但在裂纹较大和较小时相对误差较大，这是由于划分网格时，裂纹较大时裂纹处的单元不够密，较小时裂纹处的网格产生了变形。通过计算得出应力强度因子的3个分量K_Ⅰ、K_Ⅱ、K_Ⅲ后，就可以利用式(6)求出复合裂纹的应力强度因子，从而判断含裂纹构件是否安全。

表1　Ⅰ型裂纹应力强度因子随裂纹尺寸的变化

Table 1　SIF of mode Ⅰ crack in different sizes                        MPa?m^1/2

表2　Ⅰ型裂纹应力强度因子随拉伸应力的变化

Table 2　SIF of mode Ⅰ crack in different stresses                       MPa?m^1/2

表3　Ⅱ型裂纹应力强度因子随裂纹尺寸的变化

Table 3　SIF of mode Ⅱ crack in different sizes                        MPa?m^1/2

表4　Ⅱ型裂纹应力强度因子随剪切应力的变化

Table 4　SIF of mode Ⅱ crack in different stresses                        MPa?m^1/2

3　结　论

a. 计算复合型裂纹应力强度因子时可以将其简化成Ⅰ，Ⅱ和Ⅲ型应力强度因子分量进行计算，降低了复合裂纹计算的复杂性。

b. 计算裂纹的应力强度因子时，Ⅰ型和Ⅱ型应力强度因子的有限元计算结果和解析解很接近。

c. 裂纹较大和裂纹较小时计算误差相对较大。这是因为在划分网格时，裂纹较大时裂纹处的网格不够密，较小时裂纹处的单元产生了变形。

d. 通过表面效应单元SURF153在模型上加载了剪切力，在计算II型裂纹时，误差很小，表明剪切力的加载可以在ANSYS中实现。

e. 通过ANSYS计算得到的应力强度因子能较好地预测含裂纹构件的安全性，相对实验方法更经济、方便、可靠。

参考文献：

[1] 胡卫华, 吕运冰. 周期性Ⅰ-Ⅱ复合型裂纹的应力强度因子[J]. 武汉理工大学学报, 2002, 26(1): 123-125.

HU Wei-hua, LU Yun-bing. Stress intensity factors of periodic mixed crack Ⅰ-Ⅱ[J]. Journal of Wuhan University of Technology, 2002, 26(1): 123-125.

[2] 朱哲明, 谢和平. 裂纹表面承受均布载荷时的应力强度因子及裂纹张开位移的边界配位法[J]. 计算结构力学与实践, 1997, 14(2): 182-188.

ZHU Zhe-ming, XIE He-ping. Collocation method of the SIF and COD for the crack surface subjected to uniform loading[J]. Journal of Computational Structural Mechanics and Applications, 1997, 14(2): 182-188.

[3] WANG Yuan-han, WU You-lun, YU fei. Fracture calculation of bending plates by boundary collocation method[J]. Applied Mathematics and Mechanics: English Edition, 2003, 24(6): 684-690.

[4] Wang Y H, Tham LG, Lee P K K, Tsui Y. A boundary collocation method for cracked plates[J]. Computers and Structures, 2003, 81: 2621-2630.

[5] 谢和平, 朱哲明, 范天佑. 脆性岩石断裂破坏机理的边界配位法分析[J]. 力学学报, 1998, 30(2): 238-246.

XIE He-ping, ZHU Zhe-ming, FAN Tian-you. The analysis of fracture of brittle rock material by using the boundary collation method[J]. Acta Mechanica Sinica, 1998, 30(2): 238-246.

[6] 尹奇志, 肖金生, 吕运冰, 等. 孔边应力集中和裂纹尖端应力强度因子的有限元分析[J]. 武汉理工大学学报, 2002, 26(1): 47-50.

YIN Qi-zhi, XIAO Jin-sheng, LU Yun-bing, et al. Finite element analysis of stress concent ration a round a hole and stress intensity factor at crack tips[J]. Journal of Wuhan University of Technology, 2002, 26(1): 47-50.

[7] Ligoria S A, Knight G M S， Ramachandra M D S. Three-dimensional finite element analysis of a semi-elliptical circumferential surface crack in a carbon steel pipe subjected to a bending moment[J]. Journal of Strain Analysis for Engineering Design, 2005, 40(6): 525-533.

[8] Dechaumphai P, Phongthanapanich S, Sricharoenchai T. Combined delaunay triangulation and adaptive finite element method for crack growth analysis[J]. Acta Mechanica Sinica: English Series, 2003, 19(2): 162-171.

[9] 黄佩珍, 师俊平, 李中华. 求解界面裂纹应力强度因子的高次权函数法[J]. 固体力学学报, 2000, 21(2): 166-170.

HUANG Pei-zhen, SHI Jun-ping, LI Zhong-hua. Higher order weight function method for analysis of stress intensity factors for interface cracks[J]. Acta Mechanica Solid Sinca, 2000, 21(2): 166-170.

[10] Heo S P, Yang W H. Stress intensity factor analysis of elliptical corner cracks in mechanical joints by weight function method[J]. International Journal of Fracture, 2002, 115(4): 377-399.

[11] MA Kai-ping, LIU Chun-tu. Semi-weight function method on computation of stress intensity factors in dissimilar materials[J]. Applied Mathematics and Mechanics: English Edition, 2004, 25(11): 1241-1248.

[12] Heo S P, Yang W H. Mixed-mode stress intensity factors and critical angles of cracks in bolted joints by weight function method[J]. Archive of Applied Mechanics, 2002, 72(2): 96-106.

[13] CHEN Feng, CAO Ping, RAO Qiu-hua, et al. Mode Ⅱ fracture analysis of double edge cracked circular disk subjected to different diametral compressions[J]. J Cent South Univ Technol, 2004, 11(1): 63-68.

[14] 郭少华, 孙宗颀. 压应力下脆性椭圆型裂纹的断裂规律[J]. 中南工业大学学报: 自然科学版, 2001, 32(5): 457-460.

GUO Shao-hua, Sun Zhong-qi. Fracture law of an elliptical brittle crock under compressive coading[J]. Journal of Central South University of Technology: Science and Technology, 2001, 32(5): 457-460.

[15] 钱士强, 高　桦. 全复合型裂纹疲劳门槛限估算[J]. 机械强度, 1995, 17(1): 17-21.

QIAN Shi-qiang, GAO Hua. Estimation of (Ⅰ+Ⅱ+Ⅲ) mixed- mode fatigue thresholds[J]. Journal of Mechanical Strength, 1995, 17(1): 17-21.

[16] Blackburn W S, Hellen T K. 二维与三维应力强度因子的有限元计算[C]//斯坦利 P. 断裂力学的工程应用. 黄玉珊, 张相周，等译. 北京: 机械工业出版社, 1982.

Blackburn W S, Hellen T K. Computation of 2-D and 3-D stress intensify factors by FEM[C]//Stanley P. Fracture Mechanics in Engineering Practice. HUANG Yu-shan, ZHANG Xiang-zhou, translate. Beijing: Mechanical Industry Press, 1982.

[17] Patterson C. 论有限元素能量法在计算应力强度因子上的应用[C]//斯坦利 P. 断裂力学的工程应用. 黄玉珊, 张相周, 等译. 北京: 机械工业出版社，1982.

Patterson C. Study of the application of finite element energy method in computation of stress intensify factors[C]//Stanley P. Fracture Mechanics in Engineering Practice. HUANG Yu-shan, ZHANG Xiang-zhou, translate. Beijing: Mechanical Industry Press, 1982.

收稿日期：2006-05-16

基金项目：教育部高等学校博士点专项科研基金资助项目(20020533001)

作者简介：徐　慧(1958-)，男，湖南常德人，教授，从事有限元及其计算方法研究

通讯作者：徐　慧，男，教授，博士；电话：0731-8836762；E-mail: xuhui@mail.csu.edu.cn