熊 茜,馬 輝,2,官 宏
(1.東北大學(xué) 機(jī)械工程與自動化學(xué)院,沈陽 110819;2.航空動力裝備振動及控制教育部重點試驗室(東北大學(xué)),沈陽 110819)
葉片是航空發(fā)動機(jī)的核心部件之一,在航空工業(yè)中得到了廣泛的應(yīng)用。由于葉片長期受到交變載荷激勵的影響,極易誘發(fā)葉片振動并引發(fā)高周疲勞,進(jìn)而導(dǎo)致葉片出現(xiàn)裂紋,改變結(jié)構(gòu)的性能及工作穩(wěn)定性。因此,研究旋轉(zhuǎn)裂紋葉片的振動特性對提高葉片的可靠性及裂紋故障診斷具有重要的工程意義。
在過去的幾十年里,葉片振動得到了大量的研究[1-2],提出了集中質(zhì)量模型[3]、半解析模型[4-6]與有限元模型[7-8]等。為了模擬裂紋面在振動過程中出現(xiàn)的呼吸現(xiàn)象,部分學(xué)者提出了雙線性裂紋模型[9-11]。由于離心應(yīng)力和彎曲應(yīng)力的共同作用,裂紋除了存在打開與閉合兩種狀態(tài)外,還存在部分開啟狀態(tài)。因此,Zeng等[12]采用接觸有限元對升速過程的壓氣機(jī)旋轉(zhuǎn)裂紋葉片動力學(xué)特性進(jìn)行了研究。Xie等[13]考慮拉應(yīng)力與彎曲應(yīng)力的影響,提出了一種基于應(yīng)力來模擬旋轉(zhuǎn)裂紋葉片裂紋面間呼吸效應(yīng)的建模方法。Long等[14]采用傅里葉級數(shù)模擬裂紋面間的呼吸效應(yīng),對裂紋葉片的非線性特性進(jìn)行了研究。Zhang等[15]基于開閉裂紋面的面積提出了含呼吸效應(yīng)的穿透裂紋葉片動力學(xué)模型。上述研究中,大多是基于彈簧、接觸單元等對裂紋葉片的呼吸效應(yīng)進(jìn)行模擬的,這些呼吸模型難以表征裂紋面接觸區(qū)域的變化情況。
在實際工作環(huán)境中,由于葉片經(jīng)常受到氣流激振力、離心力等綜合影響,裂紋形式不僅僅表現(xiàn)為直裂紋,還存在斜裂紋的情況,如圖1所示。目前大部分學(xué)者只針對矩形截面懸臂梁以及軸受扭轉(zhuǎn)產(chǎn)生的斜裂紋進(jìn)行了研究[16-17]。Ma等[18]基于平面單元與梁單元建立了斜裂紋懸臂梁的有限元模型,分析了不同裂紋參數(shù)及載荷參數(shù)對斜裂紋懸臂梁動力學(xué)特性的影響。Liu等[19]對含斜裂紋轉(zhuǎn)軸的動力學(xué)特性進(jìn)行了研究。Ma等[20]采用三角插值函數(shù)模擬斜裂紋懸臂梁的剛度變化,并指出斜裂紋角度也是影響裂紋懸臂梁固有特性的因素之一。綜上所述可知,目前很少有研究考慮裂紋角度對葉片固有特性與響應(yīng)特性的影響。
圖1 斜裂紋示意Fig.1 Schematic of slant crack
上述文獻(xiàn)調(diào)查表明,很少有學(xué)者建立含呼吸效應(yīng)的斜裂紋扭形旋轉(zhuǎn)葉片的動力學(xué)模型。因此,本文首先推導(dǎo)了新的含斜裂紋梁單元的應(yīng)力強(qiáng)度因子與柔度矩陣,提出了一種根據(jù)斜裂紋面接觸區(qū)域來計算斜裂紋梁單元時變剛度的方法,建立了含呼吸效應(yīng)的旋轉(zhuǎn)斜裂紋扭形葉片的動力學(xué)模型;其次,通過與ANSYS有限元模型對比驗證了本文模型的有效性;最后,討論了裂紋葉片在超諧共振狀態(tài)與第1階共振狀態(tài)下,不同裂紋角度對裂紋葉片振動特性的影響。
含斜裂紋的旋轉(zhuǎn)扭形葉片模型如圖2所示,其中,OXYZ為固定坐標(biāo)系,OrXrYrZr為旋轉(zhuǎn)坐標(biāo)系,oxyz為葉片任意一個截面的局部坐標(biāo)系。L為裂紋葉片的長度,Lc為葉片根部與裂紋尖端的距離,B為葉根處葉片寬度,H為葉根處葉片厚度。扭形斜裂紋葉片在旋轉(zhuǎn)坐標(biāo)系OrXrYrZr中繞OrZr軸旋轉(zhuǎn),Rd為輪盤半徑,ω為葉片旋轉(zhuǎn)角速度。h、b分別為葉片任意橫截面處的厚度和寬度,其中,h=(1-τhx/L)H,b=(1-τbx/L)B。
圖2 旋轉(zhuǎn)扭形斜裂紋葉片結(jié)構(gòu)圖Fig.2 Structural diagram of rotating twisted slant crack blade
在本文中,將扭形葉片分為健康與含裂紋區(qū)域兩部分,其中,健康部分扭形葉片采用Timoshenko梁單元來建立,Timoshenko梁單元第k個單元的質(zhì)量矩陣、剛度矩陣、陀螺矩陣、旋轉(zhuǎn)軟化矩陣、離心剛化矩陣等詳見文獻(xiàn)[23]。
如圖3所示為含斜裂紋單元的扭形葉片示意圖。圖4所示為葉片斜裂紋梁單元示意圖,oexeyeze為斜裂紋單元坐標(biāo)系,ocxcyczc為斜裂紋面坐標(biāo)系,其中,斜裂紋面坐標(biāo)系為斜裂紋單元坐標(biāo)系繞oeze軸旋轉(zhuǎn)θ角而成。bc為裂紋單元的寬度,hc為裂紋單元的高度,lm為裂紋單元的長度,a為斜裂紋上的裂紋深度,斜裂紋面的高度為hsc=hc/cosθ,x為裂紋尖端距左節(jié)點i的距離。對于斜裂紋梁單元,其受軸向拉力P1與P7,剪切力P2、P3、P8與P9,扭矩P4與P10,彎矩P5、P6、P11與P12。
圖3 含斜裂紋的旋轉(zhuǎn)扭形葉片示意Fig.3 Schematic of rotating twisted blade with slant crack
圖4 斜裂紋梁單元示意Fig.4 Schematic of slant crack beam element
在斷裂力學(xué)中,按照力學(xué)特性可將裂紋形式劃分為:I型(張開型)裂紋、II型(滑開型)裂紋、III型(撕開型)裂紋。針對這3種不同的裂紋形態(tài),將其對應(yīng)的應(yīng)力強(qiáng)度因子分別記為:KI、KII與KIII。對于實際的裂紋體,由于裂紋產(chǎn)生的應(yīng)力集中使得裂紋尖端的應(yīng)力強(qiáng)度因子比理論值大,因此,考慮應(yīng)力集中的影響,將Ki(i=I、II、III)的表達(dá)式記為
(1)
式中:σi、Fi(γ)(i=1,2,3)分別為無裂紋體的應(yīng)力和應(yīng)力強(qiáng)度因子修正系數(shù)[16],其中γ為裂紋深度比,γ=a/hsc。
由Castigliano定理[16]可知,裂紋深度為a所引起的附加廣義位移μi為
(2)
式中:Pi、U分別為相應(yīng)的載荷和總的應(yīng)變能,總應(yīng)變能為無裂紋單元的應(yīng)變能與裂紋單元應(yīng)變能的總和,即[16]
(3)
式中:Un、Uc分別為無裂紋梁單元的應(yīng)變能和由結(jié)構(gòu)裂紋擴(kuò)展產(chǎn)生的附加應(yīng)變能。
進(jìn)一步可得到柔度矩陣系數(shù)為[24]
(4)
當(dāng)結(jié)構(gòu)中的裂紋存在一個角度時,傳統(tǒng)的直裂紋應(yīng)力強(qiáng)度因子已不再適用。因此,考慮裂紋角度的影響,結(jié)合圖4與式(1),本文重新推導(dǎo)了斜裂紋尖端的應(yīng)力強(qiáng)度因子:
(5)
式中:θ為裂紋面與oeye軸的夾角,當(dāng)θ=0°時斜裂紋退化為直裂紋,Ac為裂紋面面積,κ為剪切形狀系數(shù),對于矩形截面,本文取κ=5/6。
令:
(6)
結(jié)合式(3)與式(5)求得由裂紋引起的附加應(yīng)變能之后,再根據(jù)式(4)可以推導(dǎo)出由斜裂紋引起的附加柔度:
(7)
式中E′=E/(1-υ2),m=1+υ,υ為泊松比,其余項均為0。
結(jié)合式(4)與式(7)可求得裂紋單元柔度矩陣G為
(8)
進(jìn)一步得到開裂紋單元的剛度矩陣:
(9)
需要注意的是,為了方便起見,本文做出如下假設(shè):假設(shè)裂紋的存在對葉片的質(zhì)量矩陣、科氏力矩陣、離心剛化矩陣、旋轉(zhuǎn)軟化矩陣等沒有影響,均與健康梁單元一致,只考慮裂紋引起的結(jié)構(gòu)剛度矩陣改變。
裂紋的呼吸效應(yīng)會導(dǎo)致旋轉(zhuǎn)扭形斜裂紋葉片的剛度出現(xiàn)周期性變化,進(jìn)而使旋轉(zhuǎn)裂紋葉片的振動呈現(xiàn)非線性現(xiàn)象。目前,對于直裂紋的呼吸模型,大多數(shù)是方波模型、余弦模型等,此類呼吸模型體現(xiàn)的剛度變化過程與實際的裂紋開閉模式仍然有較大的差別,無法真實體現(xiàn)出裂紋的時變剛度特性;對于斜裂紋的呼吸模型,許多學(xué)者[24]基于應(yīng)力強(qiáng)度因子為零法對斜裂紋轉(zhuǎn)子的呼吸特性進(jìn)行了研究,然而該模型無法準(zhǔn)確表示在旋轉(zhuǎn)過程中裂紋面由于振動而導(dǎo)致的裂紋面接觸面積的變化。因此,本文根據(jù)一階剪切變形梁理論,改進(jìn)了一種根據(jù)斜裂紋面開閉區(qū)域面積來計算斜裂紋梁單元時變結(jié)構(gòu)剛度的方法[26]。根據(jù)裂紋面處裂紋的接觸特性可分為全開裂紋、全閉裂紋以及部分開裂紋,如圖5所示。其中,斜裂紋面上的彎矩可表示為
圖5 斜裂紋面呼吸狀態(tài)示意Fig.5 Schematic of breathing state of crack surface
(10)
式中:Iyc、Izc分別為裂紋處關(guān)于oyc軸、ozc軸的截面慣性矩,χ為斜裂紋面在單元坐標(biāo)系中oeyeze平面的投影與葉根處截面的夾角。
根據(jù)裂紋面積與接觸面積的比值可得圖5中3種裂紋形式下的裂紋時變剛度。圖5中E、F兩點為動點,與裂紋尖端分別相交于G、H兩點,其中:
(11)
根據(jù)圖5(b)~5(e)及式(7)與式(9)可得,部分開裂紋單元時變剛度為
(12)
閉裂紋單元裂紋時變剛度為
(13)
(14)
在旋轉(zhuǎn)坐標(biāo)系中將各個單元矩陣進(jìn)行組集,得到旋轉(zhuǎn)扭形斜裂紋葉片的動力學(xué)方程為
(15)
C=α1M+β1K
(16)
式中:α1=4πfcfc2ξ1/(fc+fc2)和β1=ξ2/(π(fc+fc2)),fc、fc2分別為葉片的前兩階固有頻率,ξ1、ξ2分別為與固有頻率相對應(yīng)的模態(tài)阻尼比,本文取ξ1=ξ2=0.02。
采用ANSYS接觸有限元模型來驗證本文所建模型的正確性。圖6所示為ANSYS Solid186單元結(jié)合Conta174與Targe170單元建立的旋轉(zhuǎn)扭形斜裂紋葉片示意圖,其中,葉片葉根處采用固定約束,即約束葉根所有節(jié)點的自由度,在ANSYS軟件中,兩個斜裂紋面間的摩擦影響忽略不計[27]。葉片的幾何參數(shù)見表1。表1中材料參數(shù)為:彈性模量E=2.068×1011Pa,密度ρ=7 861 kg/m3,泊松比υ=0.3。
表1 葉片幾何參數(shù)Tab.1 Geometric parameters of the blade
首先,選取斜裂紋角度為10°時,本文模型與ANSYS有限元模型得到的不同裂紋位置與裂紋深度下扭形斜裂紋葉片的第1階固有頻率及誤差如圖7所示。其中,當(dāng)裂紋深度比為0時表示健康葉片,本文模型與ANSYS模型的一階固有頻率分別為676.25,683.23 Hz,誤差為-1.02%。由圖7(b)誤差圖可知,本文模型與ANSYS有限元模型的最大誤差的絕對值為2.20%,證明本文所建模型的正確性。從圖7(a)中還可以看出,隨著裂紋深度比的增加,固有頻率逐漸減小;而隨著裂紋位置比的增加,固有頻率在逐漸增大。
圖7 本文模型與ANSYS有限元模型第1階固有頻率對比Fig.7 Comparison of first natural frequency between proposed model and finite element model in ANSYS
首先,對一個周期內(nèi)的裂紋呼吸狀態(tài)進(jìn)行了對比。考慮呼吸效應(yīng)的旋轉(zhuǎn)裂紋葉片第1階固有頻率可由下式近似計算:
fbreathing=2fofh/(fo+fh)
(17)
式中:fo、fh分別為開裂紋葉片第1階固有頻率和健康葉片第1階固有頻率。取裂紋位置比λ=Lc/L=0.2;裂紋深度比ζ=a/hsc=0.25;氣動力幅值A(chǔ)=0.1 MPa;斜裂紋角度θ=10°,計算100個周期且其中每個周期計算128個載荷步。在一個周期內(nèi)取9個點進(jìn)行靜力求解,如圖8、9所示。
圖8 一個周期內(nèi)的采樣點Fig.8 Sampling points in one period
圖9中百分?jǐn)?shù)表示閉裂紋區(qū)域占整個裂紋面的面積比,從圖中可以看出,當(dāng)氣動力為正時,隨著氣動力幅值的增大,裂紋面的接觸面積也相應(yīng)增大;氣動力幅值減小,接觸面積也相應(yīng)減小。當(dāng)氣動力幅值為負(fù)時表現(xiàn)為開裂紋,兩個裂紋面無接觸,因此圖9中下半周期只給出了3個采樣點的呼吸狀態(tài)。本文模型與ANSYS接觸模型均表現(xiàn)出了同樣的規(guī)律,因此本文模型可以判斷裂紋面的呼吸狀態(tài)及接觸區(qū)域變化情況。
圖9 呼吸狀態(tài)對比Fig.9 Comparison of breathing states
在非線性系統(tǒng)中,當(dāng)外激勵頻率的倍頻接近系統(tǒng)固有頻率時產(chǎn)生的共振稱為超諧共振;相反地,當(dāng)外激勵頻率的分頻接近系統(tǒng)固有頻率時產(chǎn)生的共振稱為亞諧共振。由于超諧共振的產(chǎn)生所需要的激振頻率相對較小,而第1階共振時產(chǎn)生的振動幅值較大,在實際工程應(yīng)用中,超諧共振和第1階共振往往會引起較大的危險。因此本文將激勵頻率fe分別取0.5倍和1.0倍的葉片第1階呼吸頻率來驗證旋轉(zhuǎn)扭形斜裂紋葉片的響應(yīng)特性。
取斜裂紋角度為10°,裂紋位置比λ=0.2;裂紋深度比ζ=0.25;氣動力幅值A(chǔ)=0.1 MPa;激勵頻率fe=0.5fbreathing與fe=1.0fbreathing;轉(zhuǎn)速n=60fe/ke,ke為靜子葉片數(shù),本文取10。裂紋葉片葉尖點彎曲方向的振動響應(yīng)結(jié)果如圖10所示。為了凸顯幅值較低的頻率成分,本文給出對數(shù)頻譜圖,由圖10可知,在超諧共振與第1階共振狀態(tài)下均出現(xiàn)了倍頻成分,說明裂紋的出現(xiàn)引起了非線性現(xiàn)象。在FFT譜圖中,0fe為常值分量,表示葉片裂紋引起的平衡位置的偏移量。從圖10中還可以得出:超諧共振狀態(tài)下的1.0fe幅值與2.0fe幅值大小更加接近,這是因為超諧共振狀態(tài)下2.0fe對應(yīng)的頻率接近裂紋葉片的第1階固有頻率。
圖10 θ=10°時葉尖點Yr向的振動響應(yīng)Fig.10 Vibration response of blade tip in Yr-direction at θ=10°
將本文模型與ANSYS模型得到的超諧共振與第1階共振狀態(tài)下的振動響應(yīng)進(jìn)行對比,計算時間與倍數(shù)關(guān)系見表2。其中計算機(jī)硬件條件為:Intel(R),i7-9700,CPU 3.00 GHz;本文模型自由度為600,采用MATLAB R2020b進(jìn)行計算,而ANSYS有限元模型的自由度則達(dá)到了31 239,采用ANSYS R19.2計算。根據(jù)表2中倍數(shù)關(guān)系可知,本文模型的計算效率遠(yuǎn)遠(yuǎn)大于ANSYS軟件計算效率,其中,tA表示ANSYS有限元模型計算時間,tp表示本文模型計算時間。
表2 ANSYS模型與本文模型計算效率對比Tab.2 Comparison of calculation efficiency between ANSYS model and proposed model
本文取裂紋位置比λ=0.2、裂紋深度比ζ=0.50,研究斜裂紋角度變化對裂紋葉片第1階固有頻率的影響,不同斜裂紋角度的固有頻率對比如圖11所示。從圖11中可以看出,當(dāng)裂紋角度從0°變化到80°時,根據(jù)下式可知,固有頻率增加了約3%,即隨著裂紋角度的增加,裂紋接觸面積發(fā)生改變,間接改變實際裂紋位置,使得等效裂紋位置離葉根越遠(yuǎn),系統(tǒng)結(jié)構(gòu)剛度越大,固有頻率增加。
圖11 不同裂紋角度下的固有頻率Fig.11 Natural frequencies under different slant crack angles
(18)
式中f0、f80為裂紋角度為0°和80°時的固有頻率。
本文研究了在定轉(zhuǎn)速下不同的裂紋角度對旋轉(zhuǎn)裂紋葉片響應(yīng)特性的影響。取轉(zhuǎn)速n=1 000 r/min、λ=0.2、ζ=0.50,分別研究裂紋葉片在第1階共振狀態(tài)和超諧共振狀態(tài)下,不同裂紋角度對旋轉(zhuǎn)裂紋葉片非線性特性的影響,其三維對數(shù)譜圖如圖12所示,倍頻成分始終存在于頻譜圖中,說明隨著裂紋角度的增加,由裂紋呼吸效應(yīng)誘發(fā)的非線性振動也始終存在。
圖12 不同裂紋角度下的三維譜圖Fig.12 Spectrum cascades under different slant crack angles
由于對數(shù)圖無法直觀的體現(xiàn)頻率成分對應(yīng)幅值的變化規(guī)律,本文將超諧共振和第1階共振狀態(tài)下的0fe、1.0fe、2.0fe幅值變化對比于圖13中。由圖13可知,隨著裂紋角度的增加,各個頻率成分所對應(yīng)的幅值均減小,這是因為裂紋角度越大,系統(tǒng)剛度越大,位移響應(yīng)的幅值則越小。同時,在第1階共振狀態(tài)下,斜裂紋角度從0°變化到80°時,采用式(18)相同的計算方法可知,其1.0fe幅值降低了約40%。
1)考慮裂紋角度對葉片動力學(xué)建模的影響,本文推導(dǎo)了新的含裂紋角度的應(yīng)力強(qiáng)度因子表達(dá)式與柔度矩陣;考慮葉片呼吸效應(yīng),建立了含斜裂紋旋轉(zhuǎn)扭形葉片的動力學(xué)模型,該模型能精確模擬斜裂紋葉片振動。
2)隨著裂紋角度的增加(即裂紋的等效位置離葉根越遠(yuǎn)),旋轉(zhuǎn)裂紋葉片的固有頻率越大,同時其倍頻處所對應(yīng)的幅值越來越小;并且在第1階共振狀態(tài)時1.0fe幅值變化最大,降低了約40%。
3)通過與ANSYS有限元模型對比驗證了本文模型的有效性,本文模型極大地提高了計算速度,在計算動力學(xué)響應(yīng)特性時約為ANSYS有限元軟件的22倍。