热交换管内壁子午面半椭圆裂纹应力强度因子数值计算*

李东方，杨海波，毛朝晖，张成龙

(1.衢州职业技术学院 机电工程学院，浙江 衢州 324000；2.北京科技大学 机械工程学院，北京 100083；3.中国一重集团大连工程技术有限公司 核电装备技术部，辽宁 大连 116600；4.浙江逊龙机器人科技有限公司，浙江 金华 321000)

0 引 言

承受高温高压等复杂载荷下的核电站一回路的镍基合金或奥氏体不锈钢构件，其应力腐蚀开裂是影响核电安全的重要因素之一[1-4]。要对这类薄壁圆筒构件进行疲劳寿命预判或可靠性分析，就需研究其裂纹扩展情况，因而，开展其应力腐蚀机理研究具有重要意义。

应力强度因子(K)是反映裂纹尖端弹性应力场强弱的物理量之一。目前，国内外学者对K值的求解进行了大量的研究。黄士振等[5]用边界元法计算了厚壁圆筒内壁轴向表面浅裂纹受内压的应力强度因子，比较了形状相同而深度不同以及同一深度的不同形状裂纹的解，及表面裂纹与贯穿裂纹的解；白杨[6]用有限元法对带有环向内裂纹的薄壁钢管结构进行了断裂力学计算分析；鲁丽君等[7]利用ANSYS有限元软进行实体建模法，建立了孔边角裂纹和三维穿透裂纹的模型；徐呈祥等[8]采用ANSYS有限元软件建立了受内压作用的压力管道，及其表面具有不同周向夹角斜裂纹的三维模型，并对裂纹尖端的应力强度因子进行了有限元计算及分析；张丽屏等[9]采用含真实裂纹的三维有限元法，研究了温度与压力作用下应力强度因子；Hongjun Yu等[10]提出了一种用于提取具有复杂界面的正交各向异性材料的应力强度因子的建模方法，并结合扩展有限元法和DII积分法，研究了板中的直裂纹和复合材料中的弯曲裂纹应力强度因子分布情况；Pavol Dlhy等[11]对注射成型聚合物部件由于材料收缩而产生的裂纹进行了研究，建立了求解应力强度因子的有限元模型，获得了不同裂纹长度中心裂纹的应力强度因子；Matías Braun等[12]采用格点法和线性弹性本构模型，进行了裂纹扩展速度对应力强度因子的影响的研究。但是这些方法非常复杂，计算量大，而且模型建立时还需要使用多个不易获得的参数。综合来看，目前，热-机耦合的复杂载荷下的内壁子午面含半椭圆裂纹的研究文献较少。

本文将基于断裂力学理论，采用有限元法，对I型裂纹的应力强度因子分布及影响因素进行分析与研究，研究结果可为热交换管管道安全评估提供理论和数据依据。

1 三维有限元模型建立

笔者的研究对象为核电站热交换管内壁子午面上的I型半椭圆裂纹，该裂纹所在位置及相关参数示意图如图1所示。

图1 裂纹所在位置及相关参数示意图A-裂纹尖端最深点；B-裂纹尖端表面点；a-椭圆裂纹短半轴长度；b-椭圆裂纹长半轴长度；δ-裂纹前沿离散点相对于短轴间夹角，即裂纹尖端角度；Ro-热交换管外径，为25 mm；Ri-热交换管内径，为22 mm；L0-热交换管总长度，取为290 mm；L-裂纹面中心距管底部的长度，取为145 mm

1.1 应力强度因子KI外推法理论

在断裂力学中，将裂纹分为I型(张开型)、II型(滑开型)和III型(撕开型)3种。在极坐标下，当裂纹尖端半径r→0的应力场、应变场和位移场都称作近场[13]，即裂纹尖端附近存在应力集中。目前，常用Irwin理论来确定裂纹尖端的应力强度因子。

在线弹性断裂力学中，当受到I、II或III型中的一种或以上载荷作用，裂纹尖端附近的应力场和位移场用下式来表达，即：

(1)

(2)

式中：E′，μ′—材料当量弹性模量和泊松比；E，μ—材料弹性模量和泊松比；x—裂纹正前方方向；y—裂纹面法线方向；z—裂纹前沿曲线上某点的切线方向；r，θ—裂纹尖端极坐标下的极径和极角；KI，KII，KIII—I，II和III型应力强度因子；σx，σy，σz，τxy，τyz，τzx—裂纹尖端应力场的6个分量。

当为平面应力状态问题时，E′=E、μ′=μ；当为平面应变问题时，E′=E/(1-μ2)、μ′=μ/(1-μ)。对平面应变问题，取裂纹前沿位移场表达式的泰勒展开式的奇异项，则I型应力强度因子KI为：

(3)

式中：uy(r,π)—已知裂纹表面上某一点垂直于裂纹平面的位移。

1.2 裂纹求解有限元模型建立

基于前节的应力强度因子外推法理论，笔者在ANSYS软件中编写宏命令，用来求解(3)式，即可求得某一点裂纹尖端KI值。实际工作中，热交换管承受着复杂载荷，轴向拉伸载荷，Ta，内外壁承受压力，Pi与Po，内外壁温度，Hi与Ho。有限元模型边界条件为底端固定约束。

笔者采用参数化语言[14]直接建立裂纹三维模型[15]，奇异裂纹前缘单元最终是通过宏命令赋值语句对SOLID45(或70) 8-NODE单元进行添加、移动节点，而形成1/4 SOLID95(或90) 20-NODE等参奇异单元，最终建立得到了三维裂纹有限元模型及裂纹前沿图，如图2所示。

图2 三维裂纹有限元模型及裂纹前沿图

图2中，模型材料为00Cr17Ni14Mo2[16]。00Cr17Ni14Mo2材料特性参数如表1所示。

表1 00Cr17Ni14Mo2材料特性参数

2 实验验证和结果分析

2.1 内壁子午面半椭圆裂纹实验验证

为了验证有限元模型合理性，一方面要对相关建模参数进行调试，另一方面通过实验进行验证。

2.1.1 裂纹尖端建模尺寸r对KI影响

2016年及2017年水稻季不同水肥处理下稻田水量平衡要素如表3所示。从表3可知2016年水稻季降雨量几乎是2017年的2倍，这也是2016年水稻季灌水量较少排水量较多的主要原因。2016年和2017年水稻季各处理平均灌水量分别为223.0 mm和342.5 mm，平均排水量分别为826.9 mm和215.4 mm，2016年水稻季约是2017年的4倍。从降雨量和排水量的对比可以发现，2016水稻季的单次降雨量较大。

为分析裂纹尖端建模尺寸r对KI值的影响，笔者列出KIA和KIB及其对应的相对误差。

裂纹尖端建模尺寸r对KI影响如表2所示。

表2 裂纹尖端建模尺寸r对KI影响

当r分别为0.025 mm、0.030 mm以及0.035 mm时，KIA和KIB值稳定在某一数值。如，对应的KIA分别为52.309 MPa·mm0.5、52.604 MPa·mm0.5与52.832 MPa·mm0.5，其相对误差在1%之内。因此，r可选范围0.025 mm～0.035 mm，能够保证模型稳定性。

2.1.2 KI数值解与实验解、文献解的对比

参照文献[17]中的光弹性法的实验方法，以及文献[18]中的裂纹求解结果，笔者进行实验验证和分析。取内壁半径Ri为60 mm，外壁半径Ro为125 mm，裂纹深度a为3.25 mm，裂纹形状比a/b为0.15，内壁及裂纹面上同时承受1 MPa压力的验证模型，采用机加工方法在试件上制作出裂纹，把管试件两端进行封装，然后施加上述载荷，利用光弹法测得KI值，将KI数值解与实验解、文献解进行对比。

将数据进行处理，可得：KIA数值解与实验解、文献解分别为4.640 MPa·mm0.5、4.935 MPa·mm0.5与4.667 MPa·mm0.5；KIB数值解与实验解、文献解分别为7.321 MPa·mm0.5、8.023 MPa·mm0.5、7.961 MPa·mm0.5。则A点与B点的KI有限元解相对于实验解和文献解的相对误差分别为：-5.98%、-0.58%和-6.87%、-8.84%。结果表明：有限元解与文献解的相对误差在均不超过-8.84%，能满足工程实际需要，上述误差是可以接受的。因而，所建模型是合理的。

2.2 裂纹尺寸与形状对裂纹尖端KI分布规律影响

本研究对裂纹深度a分别为1 mm、1.25 mm、1.5 mm与1.75 mm，裂纹形状比a/b分别为0.1、0.2、0.3 … 0.8和0.9的36种基本模型计算结果进行分析与讨论。

2.2.1 裂纹形状比a/b对裂纹尖端KI分布规律影响

图3 裂纹形状比a/b对裂纹尖端KI分布规律影响

分析图3(a)，在相同a条件下，随着裂纹的扩展，当a/b较小时，KI值从B点到A点呈现先逐渐减小而后又增大到某一极值的变化趋势，且在裂纹尖端A点达到最大值，总体上，呈“W”形分布，KI值在δ约为±60°处出现明显转折；随着a/b逐渐变大，这种趋势越发不明显，KI值从B点到A点逐渐减小到某一最小值(A点)，而在B点到达最大值，总体上，呈“U”形分布；当a/b分别为0.6、0.7、0.8和0.9的KIB值分别为100.080 MPa·mm0.5、100.624 MPa·mm0.5、100.612 MPa·mm0.5和100.198 MPa·mm0.5。KI分布曲线出现转折，说明KI值沿裂纹尖端B点到A点变化过程中，最小值KImin有可能出现在A点或是B点，还有可能出现在A点到B点的某个位置。

分析图3(b～d)，有上述相似结论。此外，在裂纹形状按照从扁圆扩展到近似圆，即a/b从0.1到0.9的变化过程中，a越大，KIB值越小，在a/b较小时出现稳定的分布值。如：a为1.25 mm，a/b分别为0.5、0.6、0.7、0.8和0.9的KIB值分别为110.514 MPa·mm0.5、111.369 MPa·mm0.5、111.475 MPa·mm0.5、111.077 MPa·mm0.5和110.337 MPa·mm0.5；而a为1.5 mm，a/b分别为0.4、0.5、0.6、0.7、0.8和0.9的KIB值分别为120.215 MPa·mm0.5、121.616 MPa·mm0.5、121.881 MPa·mm0.5、121.441 MPa·mm0.5、120.563 MPa·mm0.5与119.416 MPa·mm0.5；当a为1.75 mm，a/b分别为0.3、0.4、0.5、0.6、0.7、0.8和0.9的KIB值分别为128.928 MPa·mm0.5、131.423 MPa·mm0.5、132.092 MPa·mm0.5、131.663 MPa·mm0.5、130.601 MPa·mm0.5、129.189 MPa·mm0.5与127.591 MPa·mm0.5。

2.2.2 裂纹深度a对裂纹尖端KI分布规律影响

裂纹深度a对裂纹尖端KI分布规律影响如图4所示。

在a/b较小时，相同a/b条件下，随着a增大，KI分布趋近于圆滑，KIA值逐渐减小，在δ约为±60°处，不同a的KI值大致相等；当a/b较大时，KImin出现在A点；KIA值随着a增加而急剧减小。当a/b趋近于1时，即裂纹形状趋近于圆形时，KI分布曲线也越接近椭圆形。

2.3 不同复杂载荷对裂纹尖端KI分布规律的影响

不同复杂载荷的组成如表3所示。

表3 不同复杂载荷的组成

通过对不同复杂载荷下的模型进行KI的求解，对KI值进行数据处理，得到了不同复杂载荷下的KI分布，如图5所示。

图5 不同复杂载荷下的KI分布

图5中，KI分布呈现“W”型，而随着载荷中的Pi的下降，KI值下降；同时，分布曲线逐渐呈现“U”字型。

从图5可看出，从KI分布曲线形状和数值上，耦合载荷中的Pi对KI分布规律影响较大。随着载荷中的Po的下降，KI值也下降，同时分布曲线逐渐也呈现“U”字型。说明从形状和数值上，耦合载荷中的Po对KI分布规律影响也较大。但是，Pi对KI值分布影响大于Po。在L1、L6及L9条件下，KIA分别为74.821 49 MPa·mm0.5、126.407 MPa·mm0.5以及252.724 MPa·mm0.5，L9的KImax比L6下的KImax要小49.98%。由此可知，裂纹对内壁所受压力载荷变化比较“敏感”。复杂载荷中的内外壁温度差ΔT均匀增大时，KI分布趋势完全相同，均呈“W”字型，KI值呈均匀增大，这一点，从图5(b)中可清晰地看出。

从图5中还可看出：KIA、KIB及KImin值随ΔT的变化曲线均为线性变化。在L1、L2与L3下，KI分布曲线完全重合在一起，其KImax出现在B点，为98.761 MPa·mm0.5。由此可知，轴向载荷对KI分布影响较小，这也和裂纹模型的受力情况吻合。

3 结束语

本文进行了内壁子午面上含半椭圆裂纹的热交换管在复杂载荷下的应力强度因子KI分布规律的分析与研究，分析了复杂载荷对裂纹扩展的影响，并通过实验解和文献解，验证了数值解模型的可靠性；可出以下结论：

(1)分析了不同裂纹尖端尺寸的KI，表明r可选定在0.025 mm～0.035 mm范围内，其结果相对误差不超过1%。取Ri为60 mm，Ro为125 mm，a/b为0.15，内壁及裂纹面上同时承受1 MPa压力的模型，将数值解和实验解、文献解进行对比，其相对误差最大为-8.84%，3种解相互之间吻合情况良好，能满足工程实际计算需求；

(2)系统地分析了复杂载荷作用下的裂纹扩展趋势。其中，裂纹深度和形状比的增大将会使KI分布从W字型趋近于U字型。得出了内外壁温度差ΔT很大程度上决定了KI值大小，几乎呈线性增长趋势；同时，裂纹对内壁所受载荷变化更加敏感。因而，在热交换管工作过程中，应极力避免载荷峰值的出现，尤其是要关注管内外壁间压力和温度波动，做好监测工作。