块体表面裂纹应力强度因子有限元方法研究

1、第29卷 第1期2009年1月北京理工大学学报TransactionsofBeijingInstituteofTechnologyVol.29 No.1Jan.2009块体表面裂纹应力强度因子有限元方法研究肖涛, 左正兴, 廖日东(北京理工大学机械与车辆工程学院,北京 100081)摘 要:研究三维表面裂纹建模及其计算精度问题.采用有限元方法,对奇异单元和J积分两种求解裂纹强度因子方法的计算精度进行研究,并进一步研究了计算精度对网格密度和网格形式的依赖程度.研究发现,对奇异单元法而言,影响计算结果准确度的主要原因是裂纹尖端区的尺寸,网格密度的影响不大;对J积分法而言,网格形式和网格密度的影响都

2、不大,适合处理复杂结构的应力强度因子计算.关键词:应力强度因子;1/4奇异单元法;J积分法;有限元法中图分类号:O346 1 文献标识码: 文章编号:1001 0645(2009)01 0009 05FiniteElementAnalysisMethodfortheStressIntensityFactorofSurfaceCracksXIAOTao, ZUOZheng xing, LIAORi dong(SchoolofMechanicalandVehicularEngineering,BeijingInstituteofTechnology,Beijing100081,China)Abst

3、ract:Investigatesthemodelingmethodsforthreedimensionalsurfacecrackproblems.Basedonfiniteelementmethod,1/4 pointdisplacementmethodandJ integralmethodareusedtocalculatethestressintensityfactorsalongtheentirecrackfront.Severalmeshingcasesarestudied.Itisshownthat,whenusing1/4 pointdisplacementmethod,the

4、stressintensityfactorresultsaresensitivetothecrackfrontregionsizeandnon sensitivetothemeshsize.ForJ integralmethod,theresultsarenon sensitivetoboththemeshsizeandmeshstyle,andsuitabletocalculatethestressintensityfactorofcomplexstructures.Keywords:stressintensityfactor(SIF);1/4 pointdisplacementmethod

5、;J integralmethod;finiteelementmethod三维表面裂纹是内燃机、飞机、压力容器结构中的常见缺陷方式,也是最危险的裂纹之一,对含裂纹构件进行准确的应力分析是保证结构安全可靠的关键.从断裂力学角度出发,一切断裂强度的评估都离不开裂纹尖端应力强度因子的求解,但是由于三维问题本身的复杂性,很难得到有限体三维裂纹问题的准确解,研究者们均借助数值方法来求解,其中最普遍的就是有限元法.Newman Raju利用有限元法详细讨论了拉伸载荷作用下有限体表面裂纹的应力强度因子,并在计算结果的基础上,推导出相应的公式;此后,很1多学者都在探讨利用有限元法对各种各样的裂纹问题进行研究,

6、其中包括裂纹尖端应力应变的模拟2、裂纹扩展过程的模拟,以及裂纹闭合现象的模拟3,所有这些都是建立在精确的网格处理的基础之上.为此,很多学者探讨了裂纹尖端网格的划分方法.Barry等对二维裂纹的建模方法进行了研究,得出了比较理想的结果,Lin5等介绍了利用有限元法进行有限大体中三维椭圆裂纹的建模方法,Courtin等指出了利用J积分方法计算SIF的准确性和优越性,但是以上所有的研究都只是针对单元的类型,而没有对网格划分方法进行研究.64收稿日期:2008 03 27基金项目:国家部委预研项目(40402010105)( ),女,E michelle.xtcom;左正兴( 男,教授,.10北京理工

7、大学学报第29卷实际上,在对裂纹尖端进行处理时,网格的划分方法对计算结果的精度影响很大,处理不当将会导致计算结果的分散性很大,可信度降低.奇异单元法和J积分方法是应力强度因子计算最常用的两种方法,作者利用商业软件ABAQUS,详细讨论奇异单元法和J积分方法的网格划分方法,并对网格密度对计算精度的影响进行分析,从而提高断裂建模问题的一次成功率.2 半椭圆表面裂纹计算以有限大体中的表面裂纹为对象,分别采用奇异单元法和J积分计算SIF,并对计算过程和计算结果进行比较.模型如图1所示,有限块体两端受到拉伸应力的作用,在垂直于力的平面上存在一个半椭圆表面裂纹,裂纹面示意图如图1(b)所示.图中,a为裂纹

8、前缘最深点A处的深度;a,b为半椭圆裂纹的两轴;t为块体横截面的宽度;w为块体横截面的长度.块体的几何尺寸可以用量纲-参数#=t/w(块体横截面长宽比)表示,裂纹的几何尺寸可以用量纲-参数=a/b(半椭圆裂纹的形状比)和%=a/t(裂纹前缘A点的相对深度)表示.1 数值计算方法1 1 奇异单元法裂纹尖端的应力场具有奇异性,靠近裂纹尖端的各应力分量都与r-1/2成正比,r为所求应力点与裂纹尖端的距离.当r 0时,应力急剧增长,在常规的有限元法中,用多项式表示单元内部应力和位移,在奇异点附近不能很好地反映应力的变化.Barsoum7证明,将20,节点等参单元的中间节点移到裂纹尖端的1/4处,可以模

9、拟裂纹尖端的应力应变奇异性.这样,在求得裂纹尖端应力应变解后,取裂尖处1/4角点处裂纹张开位移,采用位移计算公式直接求得应力强度因子为K!=z(1/4)4(1- )(1/4).(1)式中:K!为!型裂纹应力强度因子(SIF);r(1/4)为1/4角点到裂尖的距离;sz(1/4)为1/4角点处裂纹的张开位移;E为弹性模量; 为材料的泊松比.Lin和Smith5指出,应用奇异单元法时,裂纹的附近一定区域的网格需要细化.但是在实际的处理过程中,却很难把握.本文中将对此种方法的裂纹划分方法和网格密度要求进行分析研究.1 2 J积分法J积分法是一种能量方法,是Rice8于1968年首次提出,积分路径 从

10、裂纹下表面任一点出发,沿任意路径绕过裂纹尖端,最后终止于上表面任意一点,由此可以得到裂纹扩展的应变能为Wdy-!ds.J=lim(2) 0 式中: 为积分路径;W为应变能密度;!为积分路图1 有限大体椭圆裂纹示意图Fig.1 Surfacecrackinafinitethicknessplate计算采用的有限块体模型尺寸为=0 50,#=0 25,裂纹的相对深度为%=0 25,0 50,分别采用奇异单元法和J积分方法计算.计算结果以量纲-SIF的形式进行比较,即K*!=!F!(4)式中:K*!为!型量纲-SIF;!F为名义拉伸应力.2 1 奇异单元法的网格划分方法利用奇异单元法对SIF进行计算

11、时,采用文献中常用的网格划分方法,将裂纹构件划分为包含裂纹的裂纹区和非裂纹区.划分网格时,只对裂纹区进行特殊的网格处理即可,如图2所示.图2(a)中的方体为从有限大体中取出的裂纹,径边界上的应力矢量;ds为积分路径上的位移;x,y为裂纹的局部坐标轴;u为位移矢量.在弹性情况下,J积分与SIF有如下关系KI=.-(3)第1期肖涛等:块体表面裂纹应力强度因子有限元方法研究11端区(区域3)和裂纹附近区(区域1和区域2).裂纹尖端区半圆形横截面划分形式如图2(b)所示.裂纹尖端的网格为楔形单元,周围为六面体单元,这样就在裂纹前缘构建了多圈单元.在利用1/4单元法计算SIF时,常常因为建模精度达不到要

12、求而使计算的SIF远远偏离实际的SIF.所以,需要讨论影响SIF计算精度的关键因素.本文中主要讨论以下3种因素对计算结果的影响:裂纹尖端区半径、裂纹尖端区网格密度(即裂纹尖端网格的圈数)及裂纹附近区网格密度.man Raju计算结果最接近,最大误差为2 2%.其余各种情况下,误差都很大,最大达到10 2%.图3 不同裂纹尖端区尺寸下SIF的计算结果Fig.3 SIFdistributionalongthecrackfrontfordifferentsizesintheregion3图2 奇异单元法的网格划分示意图Fig.2 Atypicalfiniteelementmeshforthecrac

13、kedplate由以上分析可知,裂纹尖端区的尺寸对计算结果的影响很大.随着裂纹尖端区半径的增大,SIF值逐渐收敛.为了进一步描述收敛趋势,图4给出了%=0 05,0 25,0 50三种情况下,随着裂纹尖端区尺寸的增大,裂纹应力强度因子的收敛曲线.2.1 1 裂纹尖端区尺寸对SIF的影响讨论裂纹尖端区的尺寸对SIF计算值的影响.为了消除单位的影响,取规格化裂纹区尺寸&=R/a进行讨论,其中,R为裂纹尖端区的半径.针对%=0 25,0 50的模型,采用4种裂纹尖端区尺寸,即&分别为0 05,0 15,0 20,0 30,计算中,保证裂纹附近区的尺寸和网格数量不变.计算结果如图3所示

14、,图中横坐标依次为裂纹前缘从表面点B到最深点A的网格节点.由图3(a)可知,随着裂纹尖端区尺寸的增大,计算结果是逐渐收敛的,在&=0 05时计算结果与其它计算值的差值比较大,为11 3%,其余3种情 况的计算结果相差不大,最大计算误差为0 5%.图中也列出了Newman Raju公式的计算结果,两者的最大误差为3 8%.由图3(b)可知,随着裂纹尖端区尺寸的增大,1图4 裂纹尖端区尺寸对计算结果的影响Fig.4 Effectofdifferentsizesoftheregion3ontheresults由图4可以看出,当&>0 30时,应力强度因子的计算值趋于一条直线.2

15、.1 2 裂纹尖端区的网格密度对SIF的影响裂纹尖端区域的网格是由包含裂纹尖端的许,可以的层改12北京理工大学学报第29卷变裂纹的疏密.分别取4,8,16层单元分析裂纹尖端区的网格密度对SIF计算值的影响.计算中,保证裂纹尖端区的尺寸大小不变.计算结果对比见图5.裂纹附近区的网格疏密将会对结果产生较大影响.2 2 J积分的网格相关性采用J积分方法计算SIF时,无须将中间节点移动至1/4单元处.在ABAQUS中,软件自动搜索裂纹周围的每圈单元计算J积分值.J积分具有与路径无关的特点,因此由各圈积分路径所得的积分值应该是相同的.2.2 1 网格形式对J积分的影响为了检验J积分计算结果与网格形式的相

16、关性,对裂纹面采用两种网格划分方法,如图7(a),7(b)所示,其中,方法1来源于文献8.图中标出了两种网格划分方法的前4圈积分路径,在计算结果的处理上,选择第24圈积分值取平均.计算结果与奇异单元法的比较见图7(c).图5 不同裂纹尖端区网格密度下SIF计算值的比较Fig.5 ComparisonofSIFdistributionalongthecrackfrontfordifferentelementsizes由图5可知,裂纹尖端区内网格密度对计算结果的影响很小,对于同一尺寸,随着层数的增多,SIF的计算值有增大的趋势,但增大的幅度很小,例如从4层增加到16层时,计算值的最大误差只有0 4

17、%.因此,裂纹尖端区内的网格层数不是关键的考虑因素.2.1 3 裂纹附近区网格密度对SIF的影响通过改变裂纹附近区单元的长度讨论裂纹附近区域的网格密度对计算结果的影响.采用量纲-参数=l/a(单元尺寸l与裂纹最深点A的深度之比)描述网格的疏密.取%=0 25,0 50,也考虑3种网格长度即=0 05,0 10,0 20,计算结果对比如图6所示.图6 裂纹附近区不同网格密度下量纲-SIF比较Fig.6 ComparisonofSIFdistributionalongthecrackfrontfordifferentelementsizeintheregion1and2图7 两种裂纹面的划分方法及

18、计算结果比较Fig.7 Twomethodsofthemeshandtheirresults由图6可见,整体上裂纹附近区网格疏密对SIF的影响不大,裂纹密度增大,计算值有收敛的趋势.需要指出的是,本次计算采用的模型是裂纹尖由图7(c)可知,两种网格划分方法得出的结果一致,网格形式对计算结果的影响不大,图7(c)中也列出了采用奇异单元法的计算结果,3种方法计K*!0 %.第1期肖涛等:块体表面裂纹应力强度因子有限元方法研究参考文献:132.2 2 网格密度对J积分的影响同样采用量纲-参数考察网格密度对于积分的影响.%=0 25时,取=0 2,0 5,1 0;%=0 50时,取=0 1,0 3,0

19、 5,计算结果如图8所示.1NewmanJJ,RajuIS.AnempiricalstressintensityfactorequationforthesurfacecrackJ.EngineeringFractureMechanics,1981,15:185 192.2ZhangB,GuoWL.Three dimensionalstressstatearoundquarter ellipticalcornercracksinelasticplatessubjectedtouniformtensionloadingJ.EngineeringFractureMechanics,2007,74:3

20、86 398.3MatosDP,NowellD.Ontheaccurateassessmentofcrackopeningandclosingstressesinplasticity inducedfatiguecrackclosureproblemsJ.EngineeringFrac tureMechanics,2007,74:1579 1601.4BarryDF,KevinZT.Anevaluationoffracturemechanicsquarter pointdisplacementtechniquesusedforcomputingstressintensityfactorsJ.S

21、tructures,1999,21:406 415.5LinXB,SmithRA.Finiteelementmodelingoffatiguecrackgrowthofsurfacecrackedplates,partI:thenu mericaltechniqueJ.EngineeringFractureMechanics,1999,63:503 522.6CourtinS,GardinC.AdvantagesoftheJ integralapproachforcalculatingstressintensityfactorswhenu singthecommercialfiniteelementsoftwareABAQUSJ.EngineeringFractureMechan