【热点导读】:建筑工程项目管理探析 高层建筑施工技术 浅谈投标报价的分析方法
基于ANSYS的有限弹性裂纹体的应力强度因子计算
关键词:断裂力学;有限单元法;ANSYS;应力强度因子
摘 要:论文在总结断裂力学各种行为研究方法的基础上,采用有限单元法思想,利用ANSYS软件建立裂纹体有限元模型。通过计算得出I型裂纹尖端的应力强度因子,其计算结果与理论值吻合良好。表明模型的选取和网格的划分是合理的,具有可靠的精度,其结果完全可以指导工程设计。
1 前言
随着现代生产技术的高速发展,新材料、新工艺在航空、航天、压力容器、核反应堆、机械、土木等领域得到了广泛应用,结构在高速、高温、高压等环境中使用时,按照传统强度理论设计的结构在应用中却出现大量的断裂事故。目前关于断裂力学的研究,理论上只能求解较简单的模型或做出较强的假设条件;通过实验探求其规律性的成本较高、周期较长。因此对于断裂力学问题,尤其是三维裂纹问题目前大多借助数值方法进行研究。裂纹问题严格来说都是三维问题,并且工程中最后发生事故的裂纹问题的物体(如机械构件,土建结构)总是有限弹性的。因此对有限弹性裂纹问题进行三维分析在实际工程上有重要的现实意义[1]。
迄今在平面断裂力学中已形成了一整套相当成熟的计算方法,但在三维断裂问题中,尤其在表面裂纹方面,还有很多问题有待进一步探讨,本文正是在这方面进行了探索和研究。
通常板表面裂纹应力强度因子可以统一表示为:
(1)
只是不同的解给修正因子F赋予不同的表达式。
由Newman的研究成果可以看出,他主要考虑了拉伸或弯曲载荷下半椭圆表面裂纹的各种裂纹形状(a/c)和板宽、板厚对裂纹前缘应力强度因子的影响,而没有考虑板长的影响,也就是没有考虑平行边界对应力强度因子的影响。
2 裂纹类型与有限元法
2.1 裂纹类型
在断裂力学中,按裂纹受力情况和裂纹面的相对位移方向将裂纹分为三种基本类型[2]:即张开型(I型),滑移型(II型),撕裂型(III型),如图1。两种或三种基本类型的组合称之为复合型裂纹问题。图1中所示三种基本类型的裂纹模型的受力特点如下:I型裂纹受垂直于裂纹面的拉应力作用;II型裂纹受平行于裂纹面而垂直于裂纹前缘的剪应力作用,又称为面内剪切型;III型裂纹受既平行于裂纹面又平行于裂纹前缘的剪应力作用,对应于反平面剪切,又称为面外剪切或纵向剪切。其对应的应力强度因子分别称为I型、II型和III型应力强度因子,记作KI,KII,KIII。
2.2 有限元方法求裂纹应力强度因子
在三种基本裂纹类型中,I型裂纹是工程上最为常见且最为危险的一种裂纹类型[3]。本文只分析I型裂纹,并采用以下方法求解应力强度因子。而作为数值解法之一的有限单元法[4]已经成为工程设计分析领域中一个强有力的计算工具,它能模拟非常复杂的构件,并得到能满足实际工程需要的数值精度要求。
图1 裂纹基本类型
由断裂力学的基本理论可知,裂纹尖端处的应力场和应变场有如下的表达形式:
(2)
计算出各点的位移矢量,最后根据计算的位移得到裂纹尖端附近的应力强度因子。
接下来将通过ANSYS计算裂纹体的应力强度因子,由于裂纹尖端应力的奇异性,我们采用四分之一节点单元(quarter-point element)来模拟裂缝尖端的应力应变场,而在其它区域布置等参单元,从而建立含表面裂纹的三维有限体的有限元模型,得到裂缝尖端的应力应变场的数值解,求出其应力强度因子。
3 有限弹性裂纹体的应力强度因子计算
3.1 算例分析
本文计算模型(由对称性,模型为整体的1/4)如图3所示,厚度t=6.35mm,宽度2b=254mm,含有长度为2a=50.8mm的中心穿透裂纹的平板,在均匀拉应力 =3890Pa作用下,弹性模量E=207GPa,泊松比μ=0.3,求解在拉力作用下,裂纹尖端的应力强度因子。
图2 裂纹前缘法平面示意图 图3 均布拉力作用下含裂纹板计算模型
首先使用三维单元SOLID95和三维单元SOLID45求解,最后用二维平面应变单元PLANE82求解。宏程序FRACT可以将三维SOLID45模型通过添加1/4位置的中间节点转化成SOLID95裂纹尖端单元。
3.2 建模与求解
(1) 定义生成裂纹尖端单元的宏程序。
定义将二维SOLID45模型通过添加1/4位置的中间节点转化成SOLID95裂纹尖端单元的宏程序。
(2) 在前处理模块设置工程选项、分析类型、单元类型和材料参数。
(3) 在柱坐标系下定义节点,并定义单元连接。
用单元定义命令过指定节点定义单元,用单元生成命令按前1个单元的节点连接样式,循环8次生成7个单元,每次单元的节点号增量20生成围绕在裂纹尖端的单元。前面定义过位于1号节点的8个重合的节点就是为了方便这一组单元的定义而设置的。
按前8个单元的节点连接样式,循环9次生成64个单元,每次单元的节点号增量1。这一句命令生成围绕在前面生成的那一组单元外围。过指定节点定义单元。按前1个单元的节点连接样式,循环3次生成2个单元,每次单元的节点号增量为1。过指定节点定义单元,得到的最外侧左面一列单元(如图4所示)。
改变单元为第2类单元类型。将单元从第1种类型改变为第2种类型,重复执行8次改变单元类型的操作,将单元1到单元8改变为第2种单元类型。最后合并重合节点,这样就可以将坐标原点的点合并为一个节点。
图4 含裂纹板1/4模型整体单元布置 图5含裂纹板1/4模型的位移约束和荷载
(4) 定义位移约束。
此时的有限元模型变成如图5所示的样子。在左侧对称截面上有左右对称的位移约束,底部截开的截面上有上下对称的位移约束。所有节点的Z向位移都被固定,上表面有表面拉应力作用。
(5) 在求解模块开始求解。
(6) 在后处理模块中,定义路径后计算应力强度因子。
计算具有对称边界的一般模型的应力强度因子KI,系统提示计算结果给出的应力强度因子结果是0.0369 。提取应力强度因子计算结果,保存到变量KI1。
(7) 定义完成J积分计算的宏程序JIN1,在后处理模块用J积分确定应力强度因子KI。
表达式“CON1=207E3/(1-(0.3*0.3))”是从J积分到应力强度因子KI的转换关系。表达式“KI2=SQRT(CON1*JINT)”从J积分值计算应力强度因子KI,分别显示应力强度因子的计算结果KI1和KI2。其值为KI1=0.0369,KI2=0.0363。
(8) 定义数组,将结果保存到这些数组。然后改用PLANE82单元重新计算。
(9) 定义关键点和线,对称位移条件后将面划分为面单元。
通过对边缘关键点、线的定义和相应的单元划分尺寸定义,位移约束等。最终通过面单元网格划分命令得到如图6所示的单元网格。
图6利用自由网格得到含裂纹板1/4模型
(10) 进入求解模块求解。
(11) 在后处理模块显示计算结果。
在后处理模块中,利用和步骤(6)类似的方法,得到用KCAL命令计算出的应力强度因子。用二维PLANE82单元计算KI1=0.0368。
(12)由J积分确定应力强度因子KI 。
先切换到柱坐标系,用节点选择命令提取到裂纹附近的节点编号NOD4,NOD5,NOD6和NOD7用数据块STINFC并指定路径上的节点。利用和步骤(7)相同的方法,得到显示的计算结果KI1和KI2。用二维PLANE82单元计算的应力强度因子KI1=0.0368,KI2=0.0367。
(13) 定义数组,将结果保存到这些数组。从数据库文件中恢复计算结果后,将其输出到结果文件。分别将采用SOLID95和 SOLID45 (3-D 分析)计算结果、PLANE82 (2-D 分析)计算结果以及文献[6](本文将该文献结果作为理论计算结果即精确值,其值为0.0356 )给出的应力强度因子列于表1。
表1 用不同的二维单元和三维单元计算含裂纹板的裂尖应力强度因子KI
4 结论
由前述计算可见,采用ANSYS软件建立的裂纹体应力强度因子有