萬 周 劉璟澤 張大海 陳 強(qiáng) 唐振寰 費(fèi)慶國(guó)
(1東南大學(xué)機(jī)械工程學(xué)院, 南京 211189)(2東南大學(xué)江蘇省空天機(jī)械裝備工程研究中心, 南京 211189)(3中國(guó)航空發(fā)動(dòng)機(jī)集團(tuán)有限公司湖南動(dòng)力機(jī)械研究所, 株洲 412002)
航空發(fā)動(dòng)機(jī)的高保真建模是其研制過程中的重要問題[1].支承參數(shù)是航空發(fā)動(dòng)機(jī)動(dòng)態(tài)特性的關(guān)鍵影響因素,對(duì)于航空發(fā)動(dòng)機(jī)整機(jī)動(dòng)力學(xué)建模至關(guān)重要,而支承結(jié)構(gòu)中的非線性因素對(duì)支承參數(shù)識(shí)別存在一定影響,因此有必要開展考慮非線性的支承參數(shù)識(shí)別研究.
準(zhǔn)確有效的航空發(fā)動(dòng)機(jī)部件模型是精確識(shí)別支承參數(shù)的基礎(chǔ).近年來,航空發(fā)動(dòng)機(jī)部件建模技術(shù)不斷發(fā)展,基于梁?jiǎn)卧暮娇瞻l(fā)動(dòng)機(jī)建模方法利用非線性力模型表征支承結(jié)構(gòu)中的非線性因素,目前已成為十分有效的建模方法[2-5].
獲得準(zhǔn)確的轉(zhuǎn)靜子模型后,一些學(xué)者針對(duì)航空發(fā)動(dòng)機(jī)支承參數(shù)識(shí)別問題開展了相關(guān)研究.屈美嬌等[6]運(yùn)用支持向量機(jī)(SVM)與遺傳算法實(shí)現(xiàn)了航空發(fā)動(dòng)機(jī)靜態(tài)支承剛度識(shí)別.費(fèi)慶國(guó)等[7]提出了一種基于徑向基神經(jīng)網(wǎng)絡(luò)(RBFNN)的有限元模型修正方法.以上研究雖未能實(shí)現(xiàn)復(fù)雜非線性結(jié)構(gòu)的參數(shù)識(shí)別,但其提出的機(jī)器學(xué)習(xí)參數(shù)識(shí)別方法為解決復(fù)雜非線性結(jié)構(gòu)的動(dòng)力學(xué)反問題提供了新的研究思路.
隨著機(jī)器學(xué)習(xí)的發(fā)展,性能強(qiáng)大的神經(jīng)網(wǎng)絡(luò)層出不窮[8-10].長(zhǎng)短期記憶(LSTM)神經(jīng)網(wǎng)絡(luò)在處理時(shí)間序列時(shí)表現(xiàn)出優(yōu)秀的長(zhǎng)期預(yù)測(cè)能力,使得求解復(fù)雜非線性結(jié)構(gòu)的動(dòng)力學(xué)反問題成為可能.
本文提出一種基于LSTM的航空發(fā)動(dòng)機(jī)整機(jī)支承剛度識(shí)別方法.首先,基于整機(jī)模型獲得不同支承剛度對(duì)應(yīng)的位移響應(yīng);然后,使用深度學(xué)習(xí)網(wǎng)絡(luò)構(gòu)建支承剛度與位移響應(yīng)之間的非線性關(guān)系;最后,利用訓(xùn)練好的深度學(xué)習(xí)網(wǎng)絡(luò)實(shí)現(xiàn)支承剛度識(shí)別.本文通過與真實(shí)值的對(duì)比驗(yàn)證了該方法的有效性,通過與基于RBFNN和SVM構(gòu)建的識(shí)別模型進(jìn)行對(duì)比,驗(yàn)證了該方法的優(yōu)越性.
由軸承外圈、軸承內(nèi)圈和數(shù)個(gè)滾珠組成的滾動(dòng)軸承模型如圖1所示,其中,O1為軸承外圈圓心,O2為軸承內(nèi)圈圓心,r為軸承滾道內(nèi)徑,R為滾道外徑,Nb為滾珠個(gè)數(shù),ωb為軸承保持架的轉(zhuǎn)速.
對(duì)于本文所研究的滾動(dòng)軸承,其滾珠在滾道內(nèi)等距排列,在滾道內(nèi)做純滾動(dòng)運(yùn)動(dòng),軸承內(nèi)外圈與滾珠的接觸符合非線性赫茲接觸理論.設(shè)軸承游隙為δ0,對(duì)于第j個(gè)滾珠,只有當(dāng)其與滾道的法向接觸變形量,δj=xcosθj+ysinθj-δ0>0時(shí),該滾珠與滾道之間才會(huì)產(chǎn)生作用力Fj,基于非線性赫茲接觸理論,軸承產(chǎn)生的軸承力為[11]
圖1 滾動(dòng)軸承示意圖
(1)
(2)
式中,F(xiàn)j=Cbδj3/2H(δj),Cb為赫茲接觸剛度,H為亥維賽函數(shù);θj為第j個(gè)滾珠轉(zhuǎn)過的角度,θj=ωbt+2π(j-1)/Nb,ωb=ωrr/(R+r),ωr為轉(zhuǎn)子轉(zhuǎn)速.滾動(dòng)軸承的變剛度振動(dòng)(VC)頻率為轉(zhuǎn)子旋轉(zhuǎn)頻率的BN倍,即ωVC=ωrBN=ωrNbr/(R+r)[12].
圖2為擠壓油膜阻尼器動(dòng)力學(xué)模型.其中,Os為阻尼器內(nèi)環(huán)中心,xs與ys分別為Os的橫坐標(biāo)和縱坐標(biāo),Oc為阻尼器外環(huán)中心,xc與yc分別為Oc的橫坐標(biāo)和縱坐標(biāo),Rs與Rc分別為阻尼器內(nèi)外環(huán)半徑,F(xiàn)t與Fr分別對(duì)應(yīng)油膜的切向與徑向作用力.
圖2 擠壓油膜阻尼器動(dòng)力學(xué)模型
基于短軸承半油膜模型,可知切向油膜力與徑向油膜力分別為[13]
(3)
(4)
其中
(5)
(6)
式中,μ為滑油黏度;Ro為油膜半徑,Ro=(Rs+Rc)/2;Lo為油膜長(zhǎng)度;Co為油膜間隙,Co=Rc-Rs.
作用于水平、垂直方向的油膜力分量為
Fox=-Frcosφ+Ftsinφ
(7)
Foy=-Frsinφ-Ftcosφ
(8)
本文的研究對(duì)象為如圖3所示的航空發(fā)動(dòng)機(jī)整機(jī)模型.圖中,klx、kly和clx、cly分別為左側(cè)鼠籠彈性支承的橫、縱向支承剛度和阻尼,krx、kry和crx、cry分別為右側(cè)鼠籠彈性支承的橫、縱向支承剛度和阻尼,cblx、cbly和cbrx、cbry分別為左右側(cè)滾動(dòng)軸承橫、縱向阻尼,kflx、kfly和kfrx、kfry分別為前后安裝節(jié)橫、縱向支承剛度,cflx、cfly和cfrx、cfry分別為前后安裝節(jié)橫、縱向阻尼.Folx、Foly和Forx、Fory分別為左右側(cè)擠壓油膜阻尼器產(chǎn)生的橫、縱向油膜力,F(xiàn)blx、Fbly和Fbrx、Fbry分別為左右側(cè)滾動(dòng)軸承產(chǎn)生的橫、縱向軸承力.
圖3 航空發(fā)動(dòng)機(jī)整機(jī)模型
基于Timoshenko梁?jiǎn)卧c剛性圓盤單元對(duì)轉(zhuǎn)子進(jìn)行建模,將機(jī)匣視為不旋轉(zhuǎn)的軸,基于Timoshenko梁?jiǎn)卧⑵溆邢拊P?各單元節(jié)點(diǎn)編號(hào)如圖3所示,則航空發(fā)動(dòng)機(jī)整機(jī)動(dòng)力學(xué)方程為
(9)
結(jié)構(gòu)阻尼矩陣Cs采用瑞利阻尼的形式,即Cs=α1M+α2K,其中α1、α2為常數(shù),可通過結(jié)構(gòu)固有頻率與阻尼比求得[5].
長(zhǎng)短期記憶神經(jīng)網(wǎng)絡(luò)是一種帶有“門”結(jié)構(gòu)的循環(huán)神經(jīng)網(wǎng)絡(luò).其在隱層中引入了能夠判斷舊有信息是否有用的單元(見圖4),解決了處理長(zhǎng)時(shí)間序列時(shí)易陷入梯度消失或梯度爆炸的問題,具有更強(qiáng)的長(zhǎng)期預(yù)測(cè)能力[14].
圖4 LSTM的基本單元
長(zhǎng)短期記憶神經(jīng)網(wǎng)絡(luò)基本單元中參量的計(jì)算公式如下[14]:
ft=sigmoid(Wfxxt+Wfhht-1+bf)
(10)
it=sigmoid(Wixxt+Wihht-1+bi)
(11)
gt=tanh(Wgxxt+Wghht-1+bg)
(12)
ot=sigmoid(Woxxt+Wohht-1+bo)
(13)
St=gt?it+St-1?ft
(14)
ht=τ(St)?ot
(15)
式中,ft、it、gt、ot、ht和St分別為遺忘門、輸入門、輸入節(jié)點(diǎn)、輸出門、中間輸出和狀態(tài)單元的狀態(tài);Wfx、Wfh、Wix、Wih、Wgx、Wgh、Wox、Woh分別為相應(yīng)門與輸入xt和中間輸出ht-1相乘的矩陣權(quán)重;bf、bi、bg、bo分別為相應(yīng)門的偏置;?表示向量中元素按位相乘.
本文利用航空發(fā)動(dòng)機(jī)整機(jī)模型獲得各類樣本,然后以樣本中的位移響應(yīng)為輸入、支承剛度為輸出訓(xùn)練網(wǎng)絡(luò)模型,最后將實(shí)測(cè)位移響應(yīng)輸入到事先訓(xùn)練的網(wǎng)絡(luò)模型中,直接識(shí)別出支承剛度的實(shí)際值.該方法主要步驟如下:
①建立航空發(fā)動(dòng)機(jī)整機(jī)模型,分析其結(jié)構(gòu)特點(diǎn),確定各處支承剛度的取值范圍.
②在目標(biāo)轉(zhuǎn)速下分析支承剛度對(duì)輪盤橫向位移響應(yīng)的影響,選擇影響顯著的支承剛度作為仿真變量.
③在仿真變量的取值范圍內(nèi)選取典型的剛度值,使用數(shù)值積分方法求解整機(jī)模型的不平衡響應(yīng),獲得目標(biāo)轉(zhuǎn)速下對(duì)應(yīng)于不同支承剛度的橫向位移響應(yīng).
④預(yù)處理橫向位移響應(yīng),將支承剛度值轉(zhuǎn)化為序列,以預(yù)處理后的位移響應(yīng)為輸入,支承剛度序列為輸出構(gòu)建訓(xùn)練樣本與微調(diào)樣本.
⑤在支承剛度取值范圍內(nèi)隨機(jī)選取剛度值,計(jì)算對(duì)應(yīng)的位移響應(yīng),并以步驟④中的方式將其構(gòu)造為測(cè)試樣本.
⑥建立以LSTM為核心層的深度學(xué)習(xí)網(wǎng)絡(luò)模型,使用訓(xùn)練樣本預(yù)訓(xùn)練該網(wǎng)絡(luò).
⑦設(shè)定預(yù)訓(xùn)練最大容許相對(duì)誤差δpt,使用微調(diào)樣本測(cè)試已預(yù)訓(xùn)練網(wǎng)絡(luò)模型,當(dāng)識(shí)別相對(duì)誤差小于δpt時(shí),使用微調(diào)樣本微調(diào)網(wǎng)絡(luò)模型,當(dāng)識(shí)別相對(duì)誤差大于δpt時(shí),重復(fù)步驟⑥.
⑧設(shè)定微調(diào)最大容許相對(duì)誤差δft,使用測(cè)試樣本測(cè)試已微調(diào)網(wǎng)絡(luò)模型,當(dāng)識(shí)別相對(duì)誤差小于δft時(shí),認(rèn)為該模型即為最終的支承剛度識(shí)別模型,當(dāng)識(shí)別相對(duì)誤差大于δft時(shí),重復(fù)步驟⑦.
⑨將實(shí)驗(yàn)測(cè)得的不平衡位移響應(yīng)輸入識(shí)別模型,輸出序列的平均值即為目標(biāo)轉(zhuǎn)速下目標(biāo)支承剛度的識(shí)別值.
3.1.1 模型參數(shù)
航空發(fā)動(dòng)機(jī)整機(jī)模型的材料參數(shù)如表1所示.用D表示梁?jiǎn)卧睆?,L表示梁?jiǎn)卧L(zhǎng)度,下標(biāo)數(shù)字表示節(jié)點(diǎn)編號(hào),其幾何尺寸如下:D1-2=D2-3=D4-5=10 mm,D3-4=12 mm,L1-2=L6-7=25 mm,L2-3=L7-8=165 mm,L3-4=L8-9=190 mm,L4-5=L9-10=180 mm,機(jī)匣內(nèi)徑為100 mm,外徑為110 mm,轉(zhuǎn)盤直徑為80 mm,轉(zhuǎn)盤厚度為10 mm.
整機(jī)模型支承參數(shù)如下:clx、cly、crx、cry、cblx、cbly、cbrx、cbry設(shè)為20 kN/(m·s-1),cflx、cfly、cfrx、cfry設(shè)為2 kN/(m·s-1).滾動(dòng)軸承參數(shù)如表2所示,擠壓油膜阻尼器參數(shù)如表3所示.
對(duì)于結(jié)構(gòu)阻尼矩陣Cs,取第1階阻尼比為 0.02,第2階阻尼比為0.04,求得轉(zhuǎn)子結(jié)構(gòu)阻尼矩陣的α1為5.62,α2為7.1×10-5,機(jī)匣結(jié)構(gòu)阻尼矩陣的α1為64.48,α2為3.6×10-6.
表1 部件材料參數(shù)
表2 滾動(dòng)軸承參數(shù)
表3 擠壓油膜阻尼器參數(shù)
3.1.2 模型驗(yàn)證
軸承力和油膜力作為支承結(jié)構(gòu)中的非線性因素,其準(zhǔn)確性對(duì)整機(jī)振動(dòng)特性影響較大,因而本文將分別對(duì)軸承力模型和油膜力模型進(jìn)行驗(yàn)證.為了減小轉(zhuǎn)子與機(jī)匣之間的耦合振動(dòng)對(duì)轉(zhuǎn)子非線性響應(yīng)的影響,在驗(yàn)證過程中將klx、kly、krx、kry設(shè)為15 MN/m,kflx、kfrx設(shè)為5 GN/m,kfly、kfry設(shè)為7 GN/m.
本文通過對(duì)比分析滾動(dòng)軸承的VC振動(dòng)規(guī)律驗(yàn)證軸承力模型的有效性[12].圖5(a)給出了轉(zhuǎn)速為100 r/min、軸承游隙為20 μm時(shí),基于整機(jī)模型計(jì)算獲得的節(jié)點(diǎn)3橫向位移響應(yīng);圖5(b)為橫向位移響應(yīng)對(duì)應(yīng)的幅值頻譜.
(a) 時(shí)域波形
(b) 幅值頻譜
由圖5(a)可知,轉(zhuǎn)子運(yùn)動(dòng)呈周期性,并表現(xiàn)出滾珠在滾道中運(yùn)動(dòng)的過程;由圖5(b)可知,轉(zhuǎn)子振動(dòng)表現(xiàn)出VC頻率(滾珠通過頻率)及其諧波.本文計(jì)算結(jié)果與文獻(xiàn)[12]的計(jì)算結(jié)果相符,因而驗(yàn)證了本文軸承力模型的有效性.
本文通過驗(yàn)證擠壓油膜阻尼器的減振作用來驗(yàn)證油膜力模型的有效性[13].圖6為不同轉(zhuǎn)速下轉(zhuǎn)子軸頸處(節(jié)點(diǎn)11)的穩(wěn)態(tài)位移響應(yīng)振幅.
圖6 振動(dòng)幅值對(duì)比
由圖6可知,擠壓油膜阻尼器在臨界轉(zhuǎn)速附近具有明顯的減振作用,油膜力對(duì)轉(zhuǎn)子振動(dòng)的抑制使得11號(hào)節(jié)點(diǎn)的振幅一直小于油膜間隙Co,擠壓油膜阻尼器的減振作用得以驗(yàn)證,從而驗(yàn)證了本文油膜力模型的有效性.
3.1.3 變量選擇
考慮到支承剛度與位移響應(yīng)之間的非線性映射關(guān)系可能并不單調(diào),本文采用文獻(xiàn)[6]中的變量選擇方法選取仿真變量.
根據(jù)工程經(jīng)驗(yàn),klx=kly,krx=kry,考慮其取值范圍為10~20 MN/m;kflx、kfrx取值范圍為45~55 MN/m;安裝節(jié)縱向剛度很強(qiáng),不考慮其對(duì)位移響應(yīng)的影響,將kfly、kfry設(shè)為70 MN/m.
本文的目標(biāo)轉(zhuǎn)速為3 000 r/min,通過比較該轉(zhuǎn)速下各剛度對(duì)節(jié)點(diǎn)3橫向位移響應(yīng)的影響,選擇影響顯著的支承剛度作為仿真變量.
分別設(shè)klx、kly、krx、kry為15 MN/m,kflx、kfrx為50 MN/m,單獨(dú)改變其中一個(gè)支承剛度,求解對(duì)應(yīng)的位移響應(yīng).以一段位移響應(yīng)的信號(hào)能量為比較對(duì)象,信號(hào)能量等于數(shù)字信號(hào)序列所有元素的平方和,分析各剛度對(duì)信號(hào)能量的影響,影響規(guī)律如圖7所示.
由圖7(c)和(d)可知,鼠籠彈性支承的剛度變化對(duì)信號(hào)能量的影響最為顯著.因此,將兩側(cè)鼠籠彈性支承的剛度klx、kly和krx、kry作為仿真變量,并將kflx、kfrx設(shè)為50 MN/m.
3.1.4 樣本構(gòu)造
在鼠籠彈性支承的取值范圍10~20 MN/m內(nèi)每隔1 MN/m取一個(gè)剛度值,則每個(gè)支承剛度都有11種取值,將序號(hào)為奇數(shù)的剛度取值進(jìn)行組合,獲得36種不同的剛度組合以構(gòu)造訓(xùn)練樣本.將序號(hào)為偶數(shù)的5種取值組合成另外25種剛度組合,用來構(gòu)造微調(diào)樣本.為了驗(yàn)證網(wǎng)絡(luò)模型的準(zhǔn)確性,在剛度取值范圍內(nèi)隨機(jī)生成5種剛度值,用來構(gòu)造25個(gè)測(cè)試樣本.
(a) 前安裝節(jié)橫向剛度kflx
(c) 左側(cè)鼠籠彈性支承剛度klx、kly
基于Newmark-β法獲得這些剛度組合對(duì)應(yīng)的位移響應(yīng),以節(jié)點(diǎn)1與節(jié)點(diǎn)3的橫向位移響應(yīng)為原始數(shù)據(jù),基于z-score方法將其標(biāo)準(zhǔn)化.同時(shí)將剛度值轉(zhuǎn)化為與位移響應(yīng)長(zhǎng)度相同的等值序列,并將其縮小107倍以統(tǒng)一輸入與輸出的量級(jí).最后以位移響應(yīng)為輸入數(shù)據(jù),對(duì)應(yīng)的剛度值序列為輸出數(shù)據(jù)構(gòu)造各個(gè)樣本.
建立如圖8所示的深度學(xué)習(xí)網(wǎng)絡(luò)模型.其中,輸入層節(jié)點(diǎn)數(shù)為2,長(zhǎng)短期記憶網(wǎng)絡(luò)層隱含層單元數(shù)為256,全連接層節(jié)點(diǎn)數(shù)分別為128和2,輸出層節(jié)點(diǎn)數(shù)為2,丟棄層歸零概率為50%,網(wǎng)絡(luò)模型的損失函數(shù)為均方誤差函數(shù).
圖8 深度學(xué)習(xí)網(wǎng)絡(luò)模型框架
取δpt=10%,設(shè)初始學(xué)習(xí)率為0.01,基于Adam優(yōu)化算法預(yù)訓(xùn)練網(wǎng)絡(luò),數(shù)次預(yù)訓(xùn)練后,微調(diào)樣本的支承剛度識(shí)別相對(duì)誤差分別降為5.28%與 3.11%,開始對(duì)網(wǎng)絡(luò)模型進(jìn)行微調(diào).
取δft=1%,設(shè)初始學(xué)習(xí)率為1.0×10-4,使用微調(diào)樣本對(duì)網(wǎng)絡(luò)模型進(jìn)行全局權(quán)重微調(diào).微調(diào)網(wǎng)絡(luò)權(quán)重后,25個(gè)測(cè)試樣本的識(shí)別相對(duì)誤差絕對(duì)值的均值已小于1%,將此時(shí)的深度學(xué)習(xí)網(wǎng)絡(luò)作為支承剛度識(shí)別模型.
本文為位移響應(yīng)附加信噪比為5%的高斯白噪聲模擬試驗(yàn)數(shù)據(jù),為了減小噪聲隨機(jī)性的影響,基于25個(gè)仿真樣本構(gòu)造了500個(gè)含噪聲的新樣本.為了驗(yàn)證本文方法的優(yōu)越性,在其他條件不變的情況下,分別基于RBFNN與SVM構(gòu)建支承剛度識(shí)別模型,并對(duì)含噪聲樣本進(jìn)行識(shí)別,3種網(wǎng)絡(luò)的識(shí)別結(jié)果如表4所示.
表4 3種網(wǎng)絡(luò)識(shí)別相對(duì)誤差絕對(duì)值的平均值 %
由表4可知:無噪聲時(shí),3種網(wǎng)絡(luò)的識(shí)別精度較高,LSTM具有最高的識(shí)別精度;5%噪聲時(shí),3種網(wǎng)絡(luò)均受到一定影響,但LSTM受到的影響最小,仍保持較高的識(shí)別精度.
以本文方法識(shí)別相對(duì)誤差最大的樣本為分析對(duì)象,其剛度識(shí)別值對(duì)應(yīng)的位移響應(yīng)與剛度真實(shí)值對(duì)應(yīng)的位移響應(yīng)如圖9所示,圖中實(shí)、虛線分別表示剛度真實(shí)值、剛度識(shí)別值對(duì)應(yīng)的3號(hào)節(jié)點(diǎn)橫向位移響應(yīng).由圖可知,剛度識(shí)別值對(duì)應(yīng)的位移響應(yīng)與真實(shí)值對(duì)應(yīng)的位移響應(yīng)高度一致,利用剛度識(shí)別后的航空發(fā)動(dòng)機(jī)整機(jī)模型能夠獲得更加準(zhǔn)確的位移響應(yīng),進(jìn)一步驗(yàn)證了本文方法的有效性.
圖9 橫向位移響應(yīng)對(duì)比
1) 建立了包含擠壓油膜阻尼器的航空發(fā)動(dòng)機(jī)轉(zhuǎn)子-滾動(dòng)軸承-支承-機(jī)匣耦合動(dòng)力學(xué)模型,通過仿真研究了滾動(dòng)軸承的VC振動(dòng)規(guī)律與擠壓油膜阻尼器的減振作用,并與現(xiàn)有文獻(xiàn)進(jìn)行對(duì)比,驗(yàn)證了本文模型的有效性.
2) 以航空發(fā)動(dòng)機(jī)整機(jī)模型為識(shí)別對(duì)象,基于LSTM構(gòu)建了支承剛度識(shí)別模型,該模型在無噪聲情況下的識(shí)別誤差小于1%,在5%噪聲情況下的識(shí)別誤差小于2%.
3) 相比于RBFNN與SVM,本文構(gòu)建的深度學(xué)習(xí)網(wǎng)絡(luò)具有更高的識(shí)別精度,表明了本文方法在解決支承剛度識(shí)別問題時(shí)的優(yōu)越性.
4) 本文方法利用神經(jīng)網(wǎng)絡(luò)的泛化特性對(duì)支承剛度進(jìn)行識(shí)別,無需迭代計(jì)算,避開了動(dòng)力學(xué)反問題中復(fù)雜的尋優(yōu)過程,大大減少了單次運(yùn)算的計(jì)算量,降低了計(jì)算機(jī)的運(yùn)算負(fù)荷.