循环载荷下考虑累积塑性破坏的船体缺口板CTOD理论及数值模拟研究

董 琴 ， 杨 平 ， 余志锋

（1.高性能船舶技术教育部重点实验室（武汉理工大学），武汉430063；2.武汉理工大学 交通学院，武汉 430063）

0 引 言

近年来，随着船舶主尺度的不断加大，为减轻结构重量而采用的高强度钢越来越多，船舶结构变得越来越“软”，船舶结构的应力及变形也越来越大，当应力达到材料的屈服强度时，裂纹尖端附近区域就会产生塑性应变。当塑性应变随着循环次数的增加而累积到一定程度时，必然会加速结构的断裂失效。在弹塑性断裂力学中，裂纹尖端张开位移（CTOD）和J积分是描述韧性结构材料裂纹扩展的主要参数（Hutchinson，1982）[1]；但由于J积分不适应于卸载状态，故而很少应用于循环载荷下的结构断裂评估。因而，研究建立循环载荷下基于累积塑性应变的裂纹尖端张开位移（CTOD）评估方法具有十分重要的意义。

裂纹尖端张开位移（CTOD）指的是裂纹尖端处受载荷后两个表面所张开的相对位移，其值的大小反映出材料抵抗疲劳裂纹的产生和扩展的能力，故CTOD可以作为评价材料抗脆性断裂特性的主要参量。先前，一些仅适用于线弹性状态的基于应力强度因子K的CTOD计算模型相继被提出（Dugdale，1960；Cottrell，1975；Yaowu，1998）[2-4]。Dugdale模型直接建立了外载荷与裂纹参数的对应关系而得到广泛的研究和应用。姜翠香[5]基于Dugdale模型对船舶加筋板CTOD进行了分析，找出了外载荷、刚度比等因素对其的影响规律。 有限元法(Potirniche et al，2003；Wu，2009；Chen，2011，2015)是目前研究弹塑性断裂问题的有效途径之一，因而被广泛应用于各种断裂评估中；它是通过对裂纹尖端应力场、位移场的有限元计算，研究裂纹尖端张口位移CTOD及相关参数。Potirniche等人[6]采用有限元方法计算了稳定的微观结构小裂纹其裂纹尖端的塑性区尺寸及CTOD值，将单晶体塑性理论的二维小应变本构关系引入ANSYS有限元软件中，得到了具有较高精度的CTOD值。Wu等人[7]提出了基于CTOD的数值结果估算试件断裂韧性的有效方法，方法中以轴向切口拉伸试样为研究对象，采用弹塑性有限元方法计算试样的CTOD值等。我国学者陈景杰[8-9]基于裂纹最大张口位移CMOD对CTOD进行了数值分析，消除了外载荷、模型尺寸、材料特性及裂纹长度等因素对CTOD的影响，简化了有限元计算模型。另外，也有学者基于裂纹尖端局部应变相继对CTOD进行了研究，并提出了仅适用于低应变、静载荷条件下的CTOD计算模型（Schwalbe，1994；Linkens et al，2000；Jayadevan et al，2004）[10-12]。 日本学者 Shimanuki Hiroshi[13]基于局部应力与线弹性断裂力学的关系，提出了评估局部应力的新方法，并以焊趾处含有裂纹的十字形接头为研究对象，采用二维有限元方法，讨论了局部应变与CTOD间的对应关系，实现了基于局部应变确定CTOD的思想；对于含裂纹结构在复杂载荷状况下应用这些计算模型引起了广泛讨论，目前，很多学者在复杂载荷下对CTOD进行了相应研究，并得出了相应结论。和提出CTOD计算模型一样，Bjerke等（2011）[14]，Yi等（2014）[15]和 Ostby（2005）[16]对发生大范围塑性变形的含裂纹管道结构进行了数值计算，并提出了一种简单的基于塑性应变的CTOD评估方法。Zhang（2014）[17]讨论了海洋管道结构承受大范围塑性应变时基于应变的CTOD评估模型，考虑了裂纹尺寸、材料硬化指数以及外载荷对CTOD的影响。

从以上分析可知，目前对CTOD的研究要么局限于线弹性状态，要么仅在静载荷条件下结合塑性应变理论对CTOD进行有限元分析；对于承受低周循环载荷下的含裂纹船体板的断裂评估目前很少有触及。本文结合船舶结构在航行过程中承受的低循环高应力外载荷条件，以含中心穿透裂纹船体板为研究对象，建立一种基于裂纹尖端局部区域累积塑性应变的CTOD计算模型；并对循环载荷的应力比、应力幅的影响进行数值分析。通过本文的分析为评估船舶构件在实际承载状态下的断裂失效提供参考。

1 理论分析

式中：σy为材料屈服极限，dn取决于材料的硬化行为。对于幂硬化材料，其裂纹尖端的应力应变场是Hutchinson，Rice 和 Rosengren（HRR）[19-20]奇异场。 在此种情况下，dn是 Ramberg-Osgood 硬化指数 n 的

对于理想弹塑性材料，Shih[18]对承受单向载荷作用下的含中心穿透裂纹的无限大板进行研究并给出了J积分和无限大板裂纹尖端张口位移（CTOD）之间的计算公式：函数并且由Shih[21]在平面应力和平面应变下的数值做了系列表格。对于平面应力的情况，dn由以下三次多项式给定：

循环载荷下，对于表现Masing行为的材料，循环J积分ΔJ有可靠的数学基础。公式（1）中所有的参量可由对应的循环参量表示，由Kumar等[22]和Heitmann等[23]的研究可知，循环J积分ΔJ可近似由弹性（小范围屈服）和塑性部分叠加而成。假设在平面应力的情况下，对于I型裂纹，我们可得到：

式中：n′为循环硬化指数，σcy为循环屈服强度，E为弹性模量，a为裂纹半长，Δσ为应力幅值，Δεpl为塑性应变幅值。在循环J积分ΔJ弹性部分，采用有效应力幅值考虑裂纹闭合效应，Δσeff=σmax-σop，σmax和σop分别表示为最大应力和裂纹张口应力。裂纹张口应力σop可通过经验公式估算或者由数值分析结果获取。本文中，裂纹张口应力通过Newman[24]得到，在低周疲劳作用下较小的裂纹张口应力对应较大的最大应力：

力和单轴极限拉伸应力的平均值。

由文献[25]Chaboche模型中塑性应变与累积塑性应变间关系微分可得：

式中：Δpn+1是第n+1次载荷循环后在裂纹尖端产生的等效塑性应变增量，通过Newton-Raphson迭代即可求得；故第n+1次循环后对应的累积塑性应变增量Δεn+1可得；依次更新相应参量，可求得每次对应循环下的塑性应变增量。

在循环载荷下，在第n+1次循环状态下裂纹尖端塑性应力应变关系可用Ramberg-Osgood表达为[26-27]：

将（6）式代入文献[25]求解Δpn+1的Newton-Raphson迭代公式中，即可求得第n+1次循环后裂纹尖端塑性应变增量Δpn+1，然后利用公式（5）即可获得第n+1次循环后裂纹尖端累积塑性应变增量Δεn+1。

本文中心穿透裂纹船体板仅承受单轴循环载荷作用，故在单轴循环应力作用下第n+1次循环后裂纹尖端累积塑性应变增量与该次循环下塑性应变增量关系如下：

由（5）～（7）式可得经过n+1次循环后累积递增塑性应变为：

故（8）式即为经过n+1次循环后中心穿透裂纹船体板在裂纹尖端位置处的累积递增塑性应变值。

由（3）式和（8）式可知，本文结合船舶结构在航行过程中承受的低循环高应力外载荷条件，以中心穿透裂纹船体板为研究对象，建立了基于裂纹尖端局部区域累积递增塑性应变的CTOD计算模型。

2 数值分析及讨论

采用大型有限元软件ANSYS对船体板进行有限元计算。在进行有限元模拟时，采用8节点PLANE82号单元进行网格划分，材料本构关系模型选用Chaboche非线性强化准则。船体板几何模型及其有限元计算网格模型如图1所示，2L为缺口板长度，2W为缺口板宽度，2a为缺口板裂纹长度，裂纹尖端局部区域采用奇异网格单元并在局部区域进行网格细化处理。有限元模型中选取船体板常用304号钢材，其相应材料参数选取如表1所示。

表1 304钢的应变疲劳参数Tab.1 Material parameters of 304 stainless steel

图1 船体板的几何模型及1/4有限元模型Fig.1 Geometry model and finite element model of ship plate with central through crack

2.1 循环载荷下裂纹尖端应力应变曲线

通过数值分析，基于Chaboche模型得到船体板在循环载荷作用下其裂纹尖端应力应变迟滞回线如图2所示。从图2中可以看出，在循环载荷作用下，随着循环次数的增加，船体板裂纹尖端附近区域产生了明显的塑性应变的累积增加。

2.2 船体板累积塑性与循环次数关系曲线

对图1给出的船体板的有限元模型，分析应力比和应力幅的变化研究，其对裂纹尖端累积塑性应变的影响。在计算中选取 L=150 mm，W=130 mm，a=10 mm。 变更应力比 R=0.1，R=0.2和 R=0.3，应力幅 σa=240 MPa，σa=260 MPa和 σa=280 MPa。 经过系列有限元计算，得出不同应力比和不同应力幅下裂纹尖端累积塑性应变随循环次数的变化关系曲线如图3和图4所示。

图2 基于Chaboche循环塑性模型的裂纹尖端应力应变曲线Fig.2 The stress-strain curve near the crack tip of Chaboche cyclic plastic model

图3 不同应力比下累积塑性应变随循环次数变化曲线Fig.3 The curve of accumulative plastic strain vs cyclic number under different stress ratio

图4 不同应力幅下累积塑性应变随循环次数变化曲线Fig.4 The curve of accumulative plastic strain vs cyclic number under different stress amplitude

从图3和图4可以看出，当应力比或应力幅一定时，裂纹尖端累积塑性应变随循环次数的增大而逐渐增加并趋于稳定。随应力比或应力幅的增大，裂纹尖端累积塑性应变也增加，表明循环载荷下应力比和应力幅对裂纹尖端累积塑性变形影响较大。

2.3 船体板裂纹尖端张口位移与循环次数关系曲线

对图1给出的船体板的有限元模型，分析应力比和应力幅的变化，研究其对裂纹尖端张口位移的影响。在计算中选取L=150 mm，W=130 mm，a=10 mm。变更应力比R=0.1，R=0.2和R=0.3，应力幅σa=240 MPa，σa=260 MPa和σa=280 MPa。经过系列有限元计算，得出不同应力比和不同应力幅下裂纹尖端张口位移随循环次数的变化关系曲线如图5和图6所示。

从图5和图6可以看出，当应力比或应力幅一定时，循环载荷下中心穿透裂纹板CTOD随循环次数的增加而逐渐增大并趋于稳定；随应力比或应力幅的增大，裂纹尖端CTOD也随之增加，说明循环载荷下应力比和应力幅对裂纹尖端CTOD影响较大。

图5 不同应力比下CTOD随循环次数变化曲线Fig.5 The curve of CTOD vs cyclic number under different stress ratio

图6 不同应力幅下CTOD随循环次数变化曲线Fig.6 The curve of CTOD vs cyclic number under different stress amplitude

3 循环载荷下CTOD与裂纹尖端累积塑性应变间的关系曲线

由（3）式可知，在平面应力状态下，船体板在循环载荷作用下其裂纹尖端张口位移变化量ΔCTOD与累积递增塑性应变Δεp存在一定的函数关系。本文对图1给出的模型在循环载荷作用下，研究应力比、应力幅的变化对ΔCTOD和Δεp关系曲线的影响。

3.1 应力比影响

考虑循环载荷下应力比R对船体板ΔCTOD和Δεp间函数关系的影响时，选取应力比分别为R=0.1，R=0.2和R=0.3，应力幅σa=240 MPa。在材料参数保持不变的情况下分别计算裂纹尖端张开位移ΔCTOD和裂纹尖端累积塑性应变Δεp，经过系列有限元计算得到的数据值如图7所示，并与（3）式的结果进行对比。

图7表明，在应力比和应力幅保持不变时，在循环载荷作用下船体板其裂纹尖端张开位移CTOD随着累积塑性应变增加而增大。在应力幅一定，应力比增加时，其趋势亦是如此。理论结果与有限元预测结果基本重合，分散性较小，说明其具有一定的可行性和较高的预测精度。

图7 不同应力比下CTOD随累积塑性应变变化曲线Fig.7 The curve of CTOD vs accumulative plastic strain under different stress ratio

图8 不同应力幅下CTOD随累积塑性应变变化曲线Fig.8 The curve of CTOD vs accumulative plastic strain under different stress amplitude

3.2 应力幅影响

考虑循环载荷下应力幅σa对船体板ΔCTOD和Δεp间函数关系的影响时，选取应力幅值分别为σa=240 MPa，σa=260 MPa，σa=280 MPa，应力比R=0.1。在材料参数保持不变的情况下分别计算裂纹尖端张开位移ΔCTOD和裂纹尖端累积塑性应变Δεp，经过系列有限元计算得到的数据值如图8所示，并与（3）式的结果进行对比。

图8表明，在应力幅和应力比保持不变时，在循环载荷作用下船体板其裂纹尖端张开位移CTOD随着累积塑性应变增加而增大。在应力比一定，应力幅增加时，其趋势亦是如此。理论结果与有限元预测结果基本重合，分散性较小，说明其具有一定的可行性和较高的预测精度。

4 结 论

本文基于弹塑性断裂力学理论，从循环J积分着手，以裂纹尖端累积塑性应变为重要参量，建立循环载荷下船体板CTOD理论模型，并在有限元模拟中分析了应力比、应力幅等相关因素影响，得到了如下结论：

（1）在循环载荷作用下，裂纹尖端局部区域将产生塑性应变的累积，随着循环次数的增加，裂纹尖端的累积塑性应变将逐渐增加并最终趋于稳定。

（2）研究表明：应力比、应力幅值对裂纹尖端张开位移-累积塑性应变关系曲线影响较为明显。任其中一者保持不变，另外一项增加时，船体板裂纹尖端张开位移随裂纹尖端累积塑性应变的增加而增大。这一趋势与理论结果吻合较好。

参 考 文 献：

[1] Hutchinson J W.Fundamentals of the phenomenological theory of nonlinear fracture mechanics[J].Journal of Applied Mechanics,1982,49：103-197.

[2]Dugdale D S.Yielding of steel sheets containing slits[J].Journal of the Mechanics and Physics of Solids,1960,8：100-108.

[3]Cottrell A H.Mechanisms of fracture,the 1963 Tewksbury lecture[R].Tewksbury：Tewksbury 146 SymPo.On Fracture,1963：1-27.

[4]Yaowu S,Siying S,Hidekazu M.Finite element analysis on relationships between the J-integral and CTOD for stationary cracks in welded tensile specimens[J].International Journal of Pressure Vessels and Piping,1998,75：197-202.

[5]Jiang Cuixiang.Research on fracture and crack arrest in ship structures[D].Wuhan：Ph.D.Dissertation of Huazhong University of Science and Technology,2005.

[6]Potirniehe G P,Daniewiez S R.Analysis of crack tip plasticity for microstructurally small cracks using crystal plasticity theory[J].Engineering Fracture Mechanics,2003,70：1623-1643.

[7]Wu F W,Ibrahim R N,Das R.Fracture toughness for CNT specimens from numerically obtained critical CTOD values[J].Theoretical and Applied Fracture Mechanics,2009,52：50-54.

[8]Chen Jingjie.Strength analysis method research of cracked ship structure[D].Dalian：Dalian University of Technology,2011.

[9]Chen Jingjie,Huang Yi.A study on evaluation method of crack tip reverse plastic zone size for the center cracked steel plate model under tension-compression cyclic loading[J].Engineering Fracture Mechanics,2015,133：138-151.

[10]Schwalbe K H.The crack tip opening displacement and J integral under strain control and fully plastic conditions estimated by the engineering treatment model for planes tress tension[J].Fracture Mechanics,1994,24：636-651.

[11]Linkens D,Formby C L.Astrain-based approach to fracture assessment-example applications[C]//Fifth International Conference on Engineering Structural Integrity Assessment.Cambridge,2000：45-52.

[12]Jayadevan K R,Østby E,Thaulow C.Strain-based fracture mechanics analysis of pipelines[C]//International Conference on Advances in Structural Integrity.Indian Institute of Science,Bangalore,India,2004：1352-1359.

[13]Shimanuki Hiroshi,INOUE,Takehiro.Study on the CTOD estimation method of the crack in stress concentrated area[C].Japan Welding Association Conference Proceedings,2007,25(l)：230-237.

[14]Bjerke S L,Scultori M.DNV’s strain-based fracture assessment approach for pipeline girth welds[C].International Offshore and Polar Engineering Conference,19-24,2011,21：690-697.

[15]Yi D K,Xiao Z M.One elastic-plastic fracture behavior of a bi-layered composite plate with a sub-interface crack under mixed mode loading[J].Composites Part B：Engineering,2014,60：60-73.

[16]Østby E,Jayadevan K R,Thaulow C.Fracture response of pipelines subject to large plastic deformation under bending[J].International Journal of Pressure Vessels and Piping,2005,82(3)：201-215.

[17]Zhang Y M,Xiao Z M.Strain-based CTOD estimation formulations for fracture assessment of offshore pipelines subjected to large plastic deformation[J].Ocean Engineering,2014,91：9164-72.

[18]Shih C F.Relationship between the J-integral and the crack opening displacement for stationary and extending cracks[J].Journal of the Mechanics and Physics of Solids,1981,29(4)：305-326.

[19]Hutchinson J W.Singular behavior at the end of a tensile crack in a hardening material[J].Journal of the Mechanics and Physics of Solids,1968,16：13-31.

[20]Rice J R,Rosengren G F.Plane strain deformation near a crack tip in a power-law hardening material[J].Journal of the Mechanics and Physics of Solids,1968,16：1-12.

[21]Shih C F.Tables of Hutchinson-Rice-Rosengren singular field quantities[R].MRL E-147,Materials Research Laborato ry,Brown University,1983.

[22]Kumar V,German M D,Shih C F.An engineering approach for elastic-plastic fracture analysis[R].Electric Power Research Institute(Report),EPRI NP,var paging,1981.

[23]Heitmann H H,Vehoff H,Neumann P.Life prediction for random load fatigue based on the growth behavior of microcracks[J].Fracture,1984,12：3599-3606.

[24]Newman J C.A crack opening stress equation for fatigue crack growth[J].International Journal of Fatigue,1984,24：131-135.

[25]Hu Guijuan.Plastic behavior of metals under tension-torsion loading-experimental and numerical research on yield surface evolution[D].Doctor of philosophy degree thesis,Guangxi：Guangxi University,2012.

[26]Landgraf R W,Morrow J.Determination of the cyclic stress-strain curve[J].Journal of Materials,1969,4(1)：176.

[27]Noroozi A H,Glinka G,Lambert S.A two parameter driving force for fatigue crack growth analysis[J].International Journal of Fatigue,2005,27：1277-1296.