摘要 描述了两种基于有限单元计算面形裂纹应力强度因子的方法,建议了一种创造三维有限单元网格的途径。计算方法的精度通过和其它解析解或数值解的比较得到了说明。
关键词 应力强度因子 有限元 损伤容限设计 断裂评定
无论在损伤容限设计还是在缺陷评定阶段,工程师们需要知道正在分析的构件中裂纹的应力强度因子,因为判断含裂纹构件的断裂,或者计算剩余疲劳寿命大多依赖于这一参量。因此,在断裂力学发展中,如何求取应力强度因子一直是一个重要的课题。当前已有许多方法可用来计算应力强度因子,较为典型的有解析法、边界配位法、有限单元法、边界元素法、体力法、权函数法和线弹簧模型。利用这些方法,大量的应力强度因子解已经获得,已出版的应力强度手册[1]中收编了许多典型的解。尽管如此,工程师们仍然会感到自己所需要的应力强度因子解很难找到,这是因为要解决的工程问题往往是一些受复杂载荷的构件,包含的裂纹也往往是一些不规则裂纹。
本文简单介绍了两种基于三维有限单元法计算面形裂纹应力强度因子的方法。有限单元法已经成为工程设计分析领域中一个强有力的计算工具,它能模拟非常复杂的构件。基于有限元的应力强度因子计算方法,自然也将具有卓越的工程能力。除了计算方法的介绍以外,还将简单描述一种简化网格的生成方法。最后提供了一些所得到的典型应力强度因子解,并和大家熟知的解进行了比较,以说明本文所描述的方法的可靠性。
1 计算方法简介
(图片)
图1 裂纹前沿单元网格
图1表示了笔者[2,3]建议的裂纹尖端网格形式。网格由3个半环共12个20节点三维等参单元组成,每个半环有4个单元,其中第一个半环单元的中节点被移至1/4点位置,以模拟裂纹尖端应力和应变场的奇异性[4,5]。根据线弹性断裂理论[6~8],裂纹尖端的位移场可表示为
(图片)
(图片)
ut=0 (3)
式中,n和t分别为裂纹前沿的法线方向和切线方向;z为垂直于裂纹平面方向(见图2);μ为材料剪切模量。(图片)
图2 曲线裂纹前沿坐标系 1.1 1/4点位移法
根据式(2),若裂纹表面上(θ=180°)某一点的垂直于裂纹平面的位移已知,那么应力强度因子(图片) 用从有限单元法求出的1/4点位移uz(1/4)(见图1)求应力强度因子的方法即为1/4点位移法,即(图片) 显然这一方法很简单,但要求裂纹尖端附近的应力应变场能较好地被模拟。
1.2 三维J积分法
J积分法实际上是一种能量方法,近来被广泛用来计算应力强度因子,因为两者可按下式转换
(图片) 用有限元计算J积分的方法通常有2种:回路积分法与虚拟裂纹扩展法。后者最早由Parks[9]和Hellen[10]独立提出,主要是通过移动有限元模型的节点位置来模拟裂纹扩展。最近deLorenzi[11,12]已经根据连续介质力学成功地推导了面形裂纹的能量释放计算公式(图片) 式中符号请参见文献[11]。商用有限元分析软件ABAQUS已经收编了上述基于“虚拟裂纹扩展原理”的三维J积分方法,可以直接获得面形裂纹前沿各角节点和中节点的J积分值。
2 网格生成
为了方便三维网格的生成,笔者建议了图3所示的分块方法。整个网格由裂纹块和非裂纹块组成,裂纹块包含了裂纹,其网格需要进行精心划分,因为它对应力强度因子有相当大的影响。相对于裂纹块,非裂纹块显得不是很重要,一般网格可以划分得粗一些,但必须使载荷和约束有一个比较合理的分配。裂纹块和非裂纹块网格的组合可采用“多点位移约束法”,使其位移在其接触面上基本保持连续,从而使组合后的网格保持变形协调。对于复杂含裂纹构件,只要合理地设计裂纹块和非裂纹块即可大大简化网格的生成。这种网格形式也特别适合于疲劳裂纹扩展的模拟计算,疲劳裂纹扩展后的新网格重新生成可以只对裂纹进行,而保持非裂纹块网格不变[2,3]。(图片)
图3 结构分块示意图 图4表示了一个笔者建议的三维裂纹块网格的生成过程。首先生成一个8节点等参单元二维网格,在这个网格里定义开裂表面和裂纹前沿,然后将其扩展成一个三维网格,单元的扩展有2种,见图4。一个非裂纹前沿单元扩展成一个20节点单元,而一个裂纹前沿单元扩展成6个20节点单元。围绕裂纹前沿,网格呈放射形。最后组合成的裂纹尖端网格具有图1所示的结构。显然,这种网格扩展方法进一步简化了三维网格的建立。(图片)
图4 三维裂纹块网格的形成 3 计算程序
笔者最近发展了如图5所示的计算应力强度因子的计算机程序。对于一些常见的工程裂纹,如平板中的半椭圆表面裂纹、压力容器中的半椭圆表面裂纹、平板孔边的角裂纹、圆棒中的表面裂纹等,只要输入一些关键的几何参数和载荷参数,程序将自动输出应力强度因子结果。该程序目前主要和商用软件ABAQUS相连接,在UNIX操作系统下运行。(图片)
图5 应力强度因子计算流程图 4 结果比较
4.1 无限大体中的椭圆形埋藏裂纹
图6表示了沿椭圆形埋藏裂纹前沿的应力强度因子结果。椭圆裂纹短长轴之比为a/c=0.5。结果表明,根据位移法和J积分法获得的结果非常相近。两者的误差约为1%。图中的虚线为Irwin[6]的解析解,本文结果略小于Irwin的解,最大误差小于3%。(图片)
图6 埋藏椭圆裂纹前沿应力强度因子变化 4.2 有限厚度平板中的半椭圆表面裂纹
近20年来,这一裂纹得到了深入的研究,大量的应力强度因子结果已经发表[13],其中,Newman_Raju[14]根据他们自己的有限单元分析结果提出的应力强度因子方程,被认为具有比较可靠的精度,并且易编程。这一方程已经在工程中得到了广泛的应用。图7比较了笔者和Newman_Raju的结果,笔者的结果由1/4点位移法获得。裂纹形状为半圆形a/c=1。结果表明,无论在裂纹深度点还是在表面点,笔者的结果和Newman—Raju的解误差约在5%以内,全面的比较见文献[2]。(图片)
图7 平板中半椭圆表面裂纹前沿应力强度因子比较
(a)拉伸载荷 (b)弯曲载荷 4.3 压力容器中的半椭圆形内表面裂纹
图8所示为压力容器筒体的局部网格及应力强度因子结果。筒体的壁厚和内径之比为t/Ri=0.1,载荷为均匀内压力。Raju—Newman[15]的结果也示于图8中。结果表明,笔者的结果和Raju—Newman的结果基本一致,只是当裂纹深度达到4/5的壁厚时,Raju—Newman的结果相对低一些。从图8a中可以看出,筒体网格由裂纹块和非裂纹块组成。进一步的比较见文献[16]。(图片)
图8 受内压容器内表面裂纹
(a)网格 (b)应力强度因子比较 4.4 孔边的半椭圆表面裂纹
图9比较了笔者的结果[17]和Shivarkumar—Newman[18]及Zhao等人[19,20]的结果。Shivarkumar—Newman的结果基于有限元分析,而Zhao等人用的是权函数方法。可见,结果同样符合得较好。(图片)
图9 孔边表面裂纹应力强度因子比较 5 结束语
本文描述了两种基于有限单元法计算面形应力强度因子的方法,建议了一种裂纹尖端的网格形式,提供了一种能方便创造包含裂纹构件三维网格的思路,最后比较了几个工程中常见裂纹的应力强度因子结果。结果表明,笔者的结果和一些已被广泛接受的应力强度因子解符合得相当好,这说明本文所描述的方法具有可靠的计算精度,可用来计算复杂工程裂纹的应力强度因子。
9/14/2004
|