正交各向异性圆柱壳裂纹形态求解及其自由振动特性①

谢宗伟 谭大鹏 王 彤 李 霖 殷梓超 吴佳峰

(浙江工业大学机械工程学院 杭州 310014)

0 引言

薄壁圆柱壳振动响应是一个复杂的非线性动力学问题,对其裂纹形态进行求解是该领域的一个难点。目前针对此类问题的研究对象主要为各向同性材料结构[1-4]。近年来,正交各向异性圆柱壳开始在航空航天、核电、石油等行业得到广泛应用[5-6],其结构特性引起国内外学者的关注。Liu 等人[7]基于Donnell-Mushtari 壳体理论,给出了具有经典边界条件的正交异性圆柱壳的解析过程,并得到了具有确定系数的封闭振动解。在此基础上,张爱国等人[8]为了解决边界条件的局限性问题,利用改进的Fourier 级数法,得到了正交各向异性圆柱壳自由振动的控制方程,研究了不同边界条件对正交各向异性圆柱壳结构振动行为的影响。此外,波传播法将圆柱壳的边界条件简化成对应的梁边界,可以有效地研究一般边界条件下正交各向异性圆柱形壳的自由振动特性。汪志强等人[9]基于Flügge 壳体理论和波传播方法讨论了正交各向异性圆柱壳的自由振动问题。Pang 等人[10]基于Reissner-Naghid 壳体理论,建立了控制运动方程,并利用波传播法研究了一般边界条件下,正交各向异性圆柱壳的自由振动特性。

在实际工作中,正交各向异性圆柱壳经常承受复杂的冲击作用,因而产生一定幅度的振动。在某些特殊情况下,冲击振动会致使薄壳产生疲劳裂纹,从而引发断裂,造成严重的生产事故。因此,为了确保其工作的安全,含有裂纹的正交各向异性圆柱壳的振动响应特性研究具有重要意义。关于结构中的裂纹分析,Joshi 等人[11]提出了一种含局部裂纹的正交各向异性薄板模型,分析了3 种边界条件对第一振型频率的影响。Lai 等人[12]研究了热环境条件下具有裂纹缺陷的正交各向异性矩形板的振动,考虑了裂纹矩形板的均匀加热载荷,研究了有裂纹或无裂纹板的临界屈曲温度。在这些研究中,裂纹被假设是平行于板的一侧边缘。Googarchin 等人[13]利用线弹簧模型建立壳体动力学模型,得出了含裂纹正交各向异性圆柱壳自由振动的解析解,但其只能分析轴向和周向2 种特定角度的裂纹,存在一定的局限性。

综合以上研究发现,目前关于含裂纹缺陷的正交各向异性圆柱壳振动特性的研究较少,建立的裂纹模型存在不能充分模拟裂纹的形态、扩展角度和传播方向等问题。针对上述问题,本文对含有斜裂纹的正交各向异性圆柱壳的自由振动进行了研究,基于Kirchhoff-Love 壳体理论,得到壳体运动的齐次偏微分方程组;利用线弹簧模型来表示壳体上的表面裂纹,根据壳体的正交各向异性,分析裂纹区域附近的附加应力,从而将裂纹元素引入到系统微分方程组中;求解方程组,得到含表面裂纹条件下正交各向异性圆柱壳的振动响应特性,在此基础上讨论壳体参数以及裂纹形态等对自由振动固有频率的影响,并基于不同阶数的固有频率下降比等值线图对裂纹形态进行识别,得到裂纹的几何空间分布。本文提出的方法可以简单有效地求解含表面斜裂纹的正交各向异性圆柱壳的自由振动问题,可为正交各向异性圆柱壳裂纹扩展规律及其裂纹识别方法的研究提供新的思路。

1 建模与求解

1.1 正交各向异性圆柱壳力学模型

建立如图1 所示的正交各向异性圆柱壳力学模型,为了方便计算,采用圆柱坐标系进行研究。r、θ和x分别代表径向坐标、周向坐标与轴向坐标。圆柱壳的中面半径为R、厚度为h、长度为L。

图1 正交各向异性圆柱壳模型及坐标系示意图

正交各向异性圆柱壳在坐标方向上的弹性模量为Ex、Eθ,泊松比为νx、νθ,根据Betti 原理[14],弹性模量和泊松比之间有以下关系:

选取壳体微元柱面为研究对象,设轴向、周向和径向的位移分别为u、v、w,根据Kirchhoff-Love 壳体理论,壳体的应变表达式为[15]

其中,t为时间,ρ为圆柱壳的密度。对于正交各向异性材料,Nx、Nθx、Nxθ和Nθ是单元上的薄膜应力,Mx、Mxθ和Mθ是单元内的弯矩,它们与位移u、v、w的关系为

式(5)、(6)中,Aij(i,j=1,2,6) 是拉伸刚度,Dij是弯曲刚度,且Dij=Aijh2/12,i,j=1,2,6,Aij的表达式为

其中,Gxθ为壳体材料的剪切模量。

式(4)所示为完整正交各向异性圆柱壳自由振动方程,为了模拟裂纹的存在,还需应用合适的裂纹模型。本文利用线弹簧模型来模拟表面裂纹,将该模型引入正交各向异性圆柱壳模型的重点和难点主要在于,需要根据壳体的正交各向异性计算裂纹区域的附加应力与附加弯矩。

1.2 裂纹形态求解方法

线弹簧模型[16]是基于Kirchhof-Love 壳体理论提出的解决板壳结构表面裂纹问题的简化计算模型,其主要优点是将表面裂纹这种三维问题降为二维问题进行处理,从而降低了分析和计算过程的复杂性。假设裂纹的尺寸与壳体尺寸相比很小,当曲率的影响可以忽略不计时,用板问题的裂纹柔度系数来模拟壳体裂纹是一个简单有效的方法。因此,本文利用线弹簧模型来解决正交各向异性圆柱壳的裂纹形态求解问题。

本文假设裂纹位于壳体的中心,图2 所示为含裂纹壳体微元柱面作用力和力矩分布,图中Qxr和Qθr为微元上的横向剪切力,和为裂纹区域的附加薄膜应力,和为裂纹区域的附加弯矩。

图2 含裂纹正交各向异性圆柱壳微元的力和力矩分布

Moazzez 等人[2]计算了各向同性壳体裂纹区域和远离裂纹区域的力和弯曲之间的关系。根据壳体的正交各向异性,这些关系可以表示为

式(8)～(11)中,a=Rθc是表面裂纹长度的一半,θc是裂纹张开角的一半,αc为表面裂纹与x轴之间的夹角,αtt、αtb和αbb是对称加载条件(模式1)下的裂纹柔度系数;ctt、ctb和cbb是反对称加载条件(模式2 和模式3)下的裂纹柔度系数。

关于对称加载条件下的裂纹柔度系数,模式1的应力强度因子可表示为

其中,σ和m分别为名义拉应力和名义弯曲应力。gt和gb分别为结构的薄膜和弯曲应力强度因子,分别可以表示为

其中,ξ是裂纹深度与壳体厚度之比。裂纹柔度系数可表示为

其中,i,j=t,b,s为裂纹的深度。由上式对裂纹深度进行积分,可得:

在非对称加载条件下,cij同样为ξ的函数[17]:

其中,ft和fb可以表示为

将式(8)～(11)中的裂纹项带入式(4)中可得到含表面裂纹正交各向异性圆柱壳的自由振动方程:

1.3 含裂纹圆柱壳自由振动方程的求解

将式(5)、(6)带入式(23)中,可得到3 个位移方向上的运动方程:

以简支边界条件为例,该边界条件为

采用波传播方法进行研究,考虑周向和轴向模态间的耦合作用,可将自由振动方程位移解的形式设为[18]

其中,Umn、Vmn、Wmn分别对应x、r和θ方向上的位移幅值,λ为x方向的特征波数。当边界条件为两端简支时,λ=mπ/L。n为周向的模态阶数,m为轴向半波数,ω为轴向的圆频率。将所设的3 个方向的位移解代入式(24)中的圆柱壳自由振动位移方程组中,约去公共因子化简得到系统特征方程的矩阵形式:

经求解,式(27)中系统特征方程的系数矩阵各元素为

若式(27)有非零解,则行列式的值det([Lij])=0。若已知壳体的特征波数,将特征波数代入到行列式中,按行列式规则进行展开,可得到关于自由振动固有频率ω的六次方程:

式中,ai为方程的系数,通过求解该方程,可得到3对复数根,其中最小的正实根即为对应耦合模态下正交异性圆柱壳的自由振动固有频率。本文讨论的裂纹形态因素包括裂纹的长度、角度和深度,通过求解式(37)可以得到含有不同形态裂纹的圆柱壳的固有频率,进而分析固有频率随着裂纹形态的变化规律,由此探寻基于固有频率的裂纹形态识别方法。

2 数值分析

2.1 方法验证

为了验证所提出计算方法的准确性,将本文计算方法所得结果与已发表文献进行对比。通过将式(8)～(11)中的裂纹长度a设置为0,可以得到现有模型的无裂纹结果。选取文献[19]中给出的4种材料情况进行验证,其具体材料属性如表1 所示,将正交各向异性圆柱壳的边界条件设置为两端简支。文献结果和采用本文方法计算所得结果如表2所示(壳体参数:相对长度L/R=2,相对厚度h/R=0.01,轴向半波数m=l。无量纲固有频率Ω可以表示为Ω=ωR/

表1 圆柱壳的材料属性

表2 无量纲固有频率理论计算结果与文献中结果对照

文献[19]和[20]在建立正交各向异性圆柱壳模型时,分别使用了Novozhilov 理论和Flügge 理论。从表2 中的比较结果中可以看出,本文方法计算所得结果与2 组文献中的计算结果相比,整体误差很小,其最大值均不超过0.9%。此外,本文方法在分析裂纹问题时,利用线弹簧模型将复杂的三维问题转化为二维问题,可使用Matlab 编写程序进行求解,具有计算简单、求解效率高的优点。

2.2 壳体参数对固有频率的影响

为了分析壳体参数对自由振动固有频率的影响,结合参考文献采用如表3 所示的圆柱壳参数进行计算分析。

表3 圆柱壳材料参数

取L/R=5,m=1,无量纲频率Ω随周向模态数n的变化曲线如图3(a)所示,当n=0、1、2 时,各组h/R对应的Ω几乎重合;当n＞2 时,h/R=0.05、0.1对应的Ω随着n的上升迅速增加,且h/R值越大,增加速度越快,而h/R=0.001、0.01对应的Ω值则持续下降;当n＞4 后趋于平稳或呈缓慢增加趋势。取h/R=0.1,其Ω随n的变化曲线如图3(b)所示,当n=0、1 时,壳体的L/R值越小则Ω值越高;当n＞1 时,L/R=10、20、30、40 时对应的Ω变化曲线几乎重合,其变化趋势总体上与图3(a)中曲线h/R=0.1 基本相同。由图3 可知,当n＞2 时,Ω的值只受相对厚度h/R的影响,对相对长度L/R的变化不敏感。

图3 无量纲频率随周向模态数变化曲线(m=1)

取L/R=20,n=1,其Ω随m的变化曲线如图4(a)所示,取不同的h/R值的壳体,Ω均随m的上升而增加,且各条曲线几乎重合,说明随着m的变化,Ω基本不受h/R值的影响。取壳体h/R=0.1,其Ω随m的变化曲线如图4(b)所示,随着m的上升,Ω均呈增加趋势,且壳体L/R值越小,Ω增加得越快。由图4 可知:随着m的变化,Ω只受L/R的影响,且L/R值越小,Ω增加得越快。

图4 无量纲频率随轴向模态数变化曲线(n=1)

取L=5 m、h/R=0.01、m=1、Eθ=10 GPa,具有不同Ex/Eθ比的圆柱壳体的无量纲固有频率Ω随周向波数n的变化曲线如图5(a)所示。当n＜4时,各种壳体对应的Ω逐渐减小;当n＞4 时,Ω的值逐渐增大。同时,当Ex/Eθ比越大时,无量纲固有频率Ω的值越小。

图5 不同Ex/Eθ 值的圆柱壳无量纲频率变化曲线

取L=5 m、h/R=0.01、n=1、Eθ=10 GPa,具有不同Ex/Eθ比的圆柱壳体的无量纲固有频率Ω随轴向波数m的变化曲线如图5(b)所示,随着m的增加,各种壳体对应的Ω逐渐增大,当m＞4 时,Ω的值趋于平稳。同样地,当Ex/Eθ比越大时,无量纲固有频率Ω的值越小。

2.3 裂纹形态对圆柱壳固有频率的影响

正交各向异性圆柱壳上的裂纹形态如图6 所示。为了方便表示裂纹的尺寸,定义一个无量纲参数β=θc/θmax(β∈[0,1]),其中θc是裂纹的半开角,θmax是裂纹的最大半开角。由于线弹簧模型仅在裂纹尺寸远小于壳体尺寸(特别是其半径)的情况下有效,因此,假设θmax=2°。

图6 裂纹形态示意图

为了表示裂纹的角度,引入了一个无量纲参数η=αc/αmax(η∈[0,1]),αc为表面裂纹与x轴之间的夹角,αmax=90°。另外,ξc为裂纹的深度,则裂纹深度与壳体厚度之比ξ=ξc/h(本文限制ξ∈[0,0.9])。定义第一振型(m=1,n=1) 固有频率下降比τ=[(ωi-ωc)/ωi]×100%,其中ωi和ωc分别为完整壳体和含裂纹壳体的第一振型固有频率。

在ABAQUS 软件中建立含裂纹正交各向异性圆柱壳体的三维有限元模型,并进行模态分析,以验证裂纹模型的准确性。在建立有限元模型时,根据式(25)中的约束条件对简支边界进行模拟(如图7所示)。根据在裂纹附近对有限元模型的网格进行加密处理(如图8 所示),在裂纹尖端进行进一步加密(如图9 所示),以减少裂纹尖端的应力奇异性,提高模型的计算精度。

图7 有限元模型及边界条件设置

图8 有限元模型的网格划分及裂纹区域加密

图9 裂纹尖端的网格加密细节

采用不同的网格划分策略对有限元计算结果可能有较大的影响,因此需要进行网格收敛性验证。从图10 中可以得出结论,对于全局种子大小小于0.1 的网格配置,第一振型固有频率几乎保持不变,因此在建立有限元模型时,将全局种子大小设置为0.1,局部种子为0.01～0.05。

图10 网格收敛性验证

表4 中给出了理论计算得出的简支边界条件下,含裂纹壳体固有频率与有限元仿真结果的对比(圆柱壳的材料参数与表3 中的相同,L=5 m,R=1 m,h=0.01 m,β=1,η=1,m=1)。结果表明,理论计算的结果与有限元仿真结果吻合性良好,因此建立的模型可以用来分析裂纹对正交各向异性圆柱壳固有频率的影响。

表4 含裂纹圆柱壳固有频率的理论计算结果与有限元仿真计算结果的对照(m=1)

图11 研究了不同裂纹尺寸对壳体第一振型(m=1,n=1)固有频率下降比τ的影响。图11(a)所示为η=1、L=20 m、R=1 m 时,不同相对厚度圆柱壳的τ随β变化而变化的曲线,可见固有频率随着裂纹长度增大而减小,这是由于裂纹的存在使壳体刚度减小所导致的,当h/R越小时,固有频率下降越明显。图9(b)所示为η=1、h=0.01 m、R=1 m时,不同相对长度圆柱壳的τ随β变化而变化的曲线,从图中可以看出,对于具有恒定相对厚度的壳体,随着裂纹长度的增大,固有频率逐渐减小,同时,对于恒定的裂纹尺寸,在L/R值较大的壳体中,固有频率的变化更明显。

图11 不同裂纹尺寸的频率下降比

图12 研究了不同裂纹深度对壳体第一振型固有频率下降比τ的影响。图12(a)所示为β=1、η=1、L=20 m、R=1 m 时,不同相对厚度圆柱壳的τ随裂纹相对深度ξ变化而变化的曲线,可见固有频率随着裂纹深度增大而减小,当h/R越小时,固有频率下降越明显。

图12 不同裂纹深度的频率下降比

图12(b)所示为β=1、η=1、h=0.01 m、R=1 m 时,不同相对长度圆柱壳固有频率随ξ变化而变化的曲线。从图中可以看出,对于具有恒定相对厚度的壳体,随着裂纹深度的增大,固有频率逐渐减小,同时,对于恒定的裂纹尺寸,在较长的壳体中,固有频率的变化更明显。

图13 研究了裂纹角度对裂纹壳体第一振型固有频率下降比τ的影响,图13(a)中所示为β=1、R=2 m、L/R=40 时,不同相对厚度圆柱壳的τ随η变化而变化的曲线。从图中可以看出,当η ＜0.5时,随着η的增大,各组固有频率逐渐减小;当η=0.5,即αc=90 °时,固有频率减小到最小;当η ＞0.5 时,固有频率逐渐增大。另外,固有频率的变化在R/h较大时更明显。

图13 不同裂纹角度的频率下降比

图13(b)中(β=1、R=2 m、R/h=40),不同相对长度圆柱壳的τ随η的变化与图13(a)有相似的趋势,但固有频率的变化在L/R值较小时更明显。

2.4 基于固有频率下降比的裂纹形态识别

对于含有斜裂纹的正交各向异性圆柱壳,裂纹形态的不同,会产生不同的受力状态,引起不同的裂纹扩展行为,因此对裂纹形态的识别是十分重要的。

取圆柱壳的L=40 m、R=0.5 m、h=0.01 m,定义m=2、n=1 时的固有频率下降比为τ’。分别绘制不同阶固有频率下降比τ和τ’ 随裂纹长度和角度变化的三维曲面图,如图14 所示。

由图14 可知,相同几何参数的正交各向异性圆柱壳中含有不同角度或长度的裂纹时,会得到不同的固有频率下降比,因此可以绘制固有频率下降比在不同裂纹长度β和相对角度η下的等值线图。由于裂纹点既在第一条等值线也在第二条等值线上,则这2 条等值线必会有交点,交点信息即为裂纹信息。基于上述思路,假设正交各向异性圆柱壳上有一裂纹,裂纹的长度β=0.8、角度为η=0.4,通过之前的计算得一阶固有频率下降比τ为0.48%,二阶固有频率下降比τ’ 为0.91%。根据这两阶固有频率下降比,得到其对应的固有频率下降比等值线图,如图15 所示。从图中可以看到这2 条等值线相交于β=0.8,η=0.4,即识别出了裂纹的长度和角度。该裂纹形态如图16 所示,裂纹的长度为28.2 mm,角度为36 °,深度为4 mm。裂纹形态的影响因素包括裂纹的长度、角度和深度,这些因素的改变均会导致固有频率下降比的变化。绘制不同阶固有频率下降比的等值线图,根据不同等值线图的交点可以确定裂纹长度、角度和深度,从而有效地识别了正交各向异性圆柱壳上的裂纹形态。

图15 基于等值线的裂纹形态识别

图16 裂纹形态示意图

3 实验验证

根据第2 节中建立的有限元仿真模型仿真来搭建相关实验平台(如图17 所示),采用LMSTest.Lab振动测试平台对正交各向异性圆柱壳进行模态频率测试。选择玻璃/环氧树脂材料管道来作为正交各向异性圆柱壳模型,其几何尺寸与材料参数如表5所示。

图17 实验测试平台

为了验证正交各向异性圆柱壳裂纹求解方法的准确性,确定实验的基本流程:利用冲击力锤在正交各向异性圆柱壳上施加一个锤击激励,并在其上布置加速度传感器,以获得圆柱壳的频率响应函数,在此基础上利用LMS Test.Lab 对模态进行提取,以得到圆柱壳的各阶自由振动固有频率。实验原理如图18所示。

图18 实验原理示意图

对圆柱壳进行模态频率测试,获得其各阶自由振动固有频率结果(m=1,边界条件为简支边界),并与理论和仿真值进行对比,其结果如表6 所示。从表中可以看出理论值与实验结果存在一定误差(这与外界环境干扰和实验边界条件的设置等因素有关),但误差不超过5%,且理论值与仿真结果较吻合。

表6 完整圆柱壳固有频率对比(m=1)

表7～9 分别给出了通过实验方法测得的裂纹长度、深度和角度改变时的正交各向异性圆柱壳第一阶(m=1,n=1)固有频率下降比。从各表中可以看出,理论所得的固有频率下降比与实验结果的误差均不超过3.5%,可以认为所提方法具有较好的准确性。

表7 不同裂纹长度的固有频率下降比(ξ=0.7,η=1)

表8 不同裂纹深度的固有频率下降比(β=0.7,η=1)

表9 不同裂纹角度的固有频率下降比(β=0.7,ξ=0.7)

4 结论

基于Kirchhoff-Love 壳体理论推导了正交各向异性圆柱壳偏微分形式的自由振动方程,结合线弹簧模型,得到了含斜裂纹正交各向异性圆柱壳的自由振动响应特性,分析了壳体的几何参数与裂纹的尺寸、角度等因素对自由振动固有频率的影响,并进行了仿真和实验验证,主要研究结果可归纳如下。

(1) 当轴向模态数m一定且周向模态数n＞2时,完整正交各向异性圆柱壳的无量纲固有频率Ω的值只受相对厚度h/R的影响,对相对长度L/R的变化不敏感;当周向模态数n一定时,各阶轴向模态数对应的无量纲固有频率Ω只受相对厚度L/R的影响,对相对长度h/R的变化不敏感。

(2) 表面裂纹的存在会导致正交各向异性圆柱壳体自由振动各阶固有频率的减小。对于给定的L/R值,表面裂纹对h/R较小的壳体的固有频率的影响较大,对于具有特定裂纹和h/R值的壳体,L/R值较大的壳体会有较大的频率变化。

(3) 裂纹的角度对壳体的固有频率有一定程度的影响,当η＜0.5 时,随着η的增大,各组固有频率逐渐减小;当η=0.5,即αc=90°时,固有频率减小到最小;当η＞0.5 时,固有频率逐渐增大,这一变化在较短和较薄的壳体中更为明显。

(4) 通过对裂纹长度和角度的分析,得到了不同阶固有频率及其变化规律,并利用不同阶固有频率下降比等值线图实现了对裂纹形态的识别,可为正交各向异性圆柱壳裂纹扩展规律及其裂纹识别方法的研究提供新的思路。