廖三豐

（1. 中國航發(fā)湖南動力機(jī)械研究所，湖南 株洲 412002；2. 航空發(fā)動機(jī)振動技術(shù)航空科技重點(diǎn)實驗室，湖南 株洲 412002）

高轉(zhuǎn)速是渦軸發(fā)動機(jī)轉(zhuǎn)子最典型特征之一，渦軸發(fā)動機(jī)燃?xì)獍l(fā)生器轉(zhuǎn)子采用的是高速兩支點(diǎn)轉(zhuǎn)子，即采用1-0-1支承結(jié)構(gòu)，其工作轉(zhuǎn)速一般超過40000r/min，往往在二階臨界轉(zhuǎn)速以上。轉(zhuǎn)子不平衡微弱變化將產(chǎn)生較大的不平衡激振力，所以前后彈性支承采用的是擠壓油膜阻尼結(jié)構(gòu)（SFD）。SFD的阻尼剛度表現(xiàn)出明顯的非線性特征，在轉(zhuǎn)子不平衡量識別過程中如何處理支承參數(shù)的非線性一直是一個較大的難題。

上個世紀(jì)70年代以來，發(fā)展了許多基于線性理論的無試重動平衡的方法[1-4]，無論是振型平衡法還是影響系數(shù)法都是建立在線性理論基礎(chǔ)之上的，動平衡技術(shù)通常需要多次啟動施加配重，此類方法都以降低轉(zhuǎn)子響應(yīng)為目標(biāo)，不關(guān)注轉(zhuǎn)子實際工作中不平衡量惡化規(guī)律，無法達(dá)到轉(zhuǎn)子全工作轉(zhuǎn)速內(nèi)不平衡量的識別。在旋轉(zhuǎn)機(jī)械故障診斷方面，也無法準(zhǔn)確考慮到支撐結(jié)構(gòu)非線性因素的影響，轉(zhuǎn)子不平衡過大的情況只停留定性的層面上[5-7]。國內(nèi)外已有大量關(guān)于轉(zhuǎn)子不平衡在線識別的研究，饒柱石[8]等基于參數(shù)識別、子結(jié)構(gòu)模態(tài)分析及阻抗匹配原理，通過建立轉(zhuǎn)子——軸承系統(tǒng)頻響函數(shù)實現(xiàn)對轉(zhuǎn)子不平衡量在線識別，該方法在識別過程中因計算規(guī)模小而具有快速的特征。畢士華等[9]提出了一種考慮了油膜力的轉(zhuǎn)子不平衡量時域方法，該方法通過求解線性矩陣方程組，無法給出不平衡量具體分布及相位關(guān)系。張召翠[10]等基于遺傳算法開展非線性軸承參數(shù)及不平衡量參數(shù)識別。劉淑蓮等[11]采用全息譜技術(shù)提取轉(zhuǎn)子共頻分量，從而識別非線性系統(tǒng)的不平衡量。以上識別方法計算復(fù)雜、識別變量多，在渦軸發(fā)動機(jī)實際工程應(yīng)用中難以推廣。

目前，在發(fā)動機(jī)整機(jī)動力學(xué)仿真分析研究中和轉(zhuǎn)子不平衡惡化的在線監(jiān)控的工程應(yīng)用中，缺乏一套高效的不平衡量識別都方法[12-13]。本文以彈性支承應(yīng)力為目標(biāo)函數(shù)，開展非線性支承參數(shù)下高速轉(zhuǎn)子不平衡量參數(shù)識別方法研究，實現(xiàn)了全工作轉(zhuǎn)速下轉(zhuǎn)子不平衡的“可視化”。整機(jī)試車中轉(zhuǎn)子振動測點(diǎn)布局受到發(fā)動機(jī)結(jié)構(gòu)緊湊的制約，彈性支承應(yīng)力值大小是為數(shù)不多較準(zhǔn)確表征整機(jī)環(huán)境下高速轉(zhuǎn)不平衡響應(yīng)狀況的振動信號。

1 渦軸發(fā)動機(jī)高速轉(zhuǎn)子不平衡特征

1.1 轉(zhuǎn)子非線性因素

渦軸發(fā)動機(jī)燃?xì)獍l(fā)生器轉(zhuǎn)子為典型的高速兩支點(diǎn)轉(zhuǎn)子，見圖1。在裝配中，中心拉桿上的一、二級壓緊螺母給整個轉(zhuǎn)子施加軸向預(yù)緊力，但在工作過程中會產(chǎn)生損失，主要表現(xiàn)在兩方面原因：一是發(fā)動機(jī)轉(zhuǎn)子在工作過程中，軸流壓氣機(jī)部分主要受到逆氣流方向的軸向力，渦輪部分受到氣流方向的軸向力，抵消了一部分預(yù)緊力；二是高速轉(zhuǎn)子的超高工作轉(zhuǎn)速，巨大的離心力會引起轉(zhuǎn)子變形，由于Poisson效應(yīng)會造成零部件發(fā)生軸向收縮。變小的端齒壓緊力造成“壓不住”的情況發(fā)生，連接結(jié)構(gòu)的不穩(wěn)定，轉(zhuǎn)子不平衡量隨之發(fā)生非線性變化。

高速轉(zhuǎn)子前、后支承一般是帶有鼠籠式彈性支承的定心擠壓油膜阻尼器結(jié)構(gòu)，其中結(jié)構(gòu)形式見圖2，滑油通過供油孔進(jìn)入供油槽，滑油在彈性支承與阻尼器之間的間隙中形成一圈油膜，轉(zhuǎn)子運(yùn)轉(zhuǎn)產(chǎn)生偏心量，受擾動油膜產(chǎn)生的油膜力阻礙轉(zhuǎn)子發(fā)生偏心，流動的滑油把它 “吸收”的轉(zhuǎn)子振動動能“帶”出去。擠壓油膜阻尼器剛度阻尼參數(shù)理論推導(dǎo)方法常采用“短”“長”軸承假設(shè)，阻尼器長度與直徑之比小于0.25，兩端不封嚴(yán)時，采用短軸承假設(shè)，密封較好則采用長軸承假設(shè)，具體推導(dǎo)過程及計算公式在相關(guān)設(shè)計手冊已有推導(dǎo)，本文不再贅述。

圖1 燃發(fā)轉(zhuǎn)子結(jié)構(gòu)簡圖

圖2 擠壓油膜阻尼器結(jié)構(gòu)圖

帶有擠壓油膜阻尼器的高速轉(zhuǎn)子模型見圖3，K1、K2分別是前后支承剛度，Cs和Ks分別是擠壓油膜阻尼器阻尼與剛度。基于轉(zhuǎn)子試驗臺中振動數(shù)據(jù)的不平衡量識別過程中K1、K2采用彈性支承剛度，而基于整機(jī)試車中振動數(shù)據(jù)的轉(zhuǎn)子不平衡量識別中的K1、K2應(yīng)包含機(jī)匣的剛度，是一個隨工頻變化的動剛度。剛度、阻尼的非線性會直接影響發(fā)動機(jī)轉(zhuǎn)子不平衡量識別的準(zhǔn)確性。

圖3 轉(zhuǎn)子支承剛度阻尼模型

1.2 不平衡量識別中參數(shù)描述

考慮動剛度及非線性油膜剛度阻尼的轉(zhuǎn)子運(yùn)動方程可以寫成以下形式

式中，q為轉(zhuǎn)子廣義坐標(biāo)系矩陣，m為高速轉(zhuǎn)子系統(tǒng)質(zhì)量矩陣，d為包含系統(tǒng)內(nèi)阻尼和陀螺力矩矩陣，se為包含系統(tǒng)動剛度的矩陣。為減小篇幅，以下參數(shù)介紹不帶下標(biāo)，α為待識別不平衡量位置矩陣（一般選取轉(zhuǎn)子動平衡去材料界面位置），F(xiàn)表示不平衡力，β為擠壓油膜阻尼器位置矩陣，f為油膜力。為減少參數(shù)識別變量，se通過前期的計算得到，而擠壓油膜阻尼器剛度阻尼通過試驗實測或理論公式計算得到，至此方程中待識別參數(shù)為F。轉(zhuǎn)子工作狀態(tài)下的不平衡力主要表現(xiàn)為工頻激勵。將不平衡力Fi沿水平和豎直方向分解，所得到的不平衡力表達(dá)式為

其中，F(xiàn)i為不平衡量大小，ω為轉(zhuǎn)子工作轉(zhuǎn)速，φi表示不平衡量角度。本文轉(zhuǎn)子不平衡識別主要用于整機(jī)動力學(xué)研究的變量輸入，不涉及動平衡試驗，則各位置等效不平衡量周向角度關(guān)系以相對角度表示。則(1)式可描述為

(3)式中待識別不平衡待識別參數(shù)為Fi與φi，共計 2n-1 個。

渦軸發(fā)動機(jī)整機(jī)環(huán)境中的鼠籠彈性支承肋條應(yīng)力水平是最能直接表征轉(zhuǎn)子不平衡狀態(tài)的試驗值。圖4為典型渦軸發(fā)動機(jī)鼠籠彈性支承結(jié)構(gòu)，軸承旋轉(zhuǎn)載荷作用下，肋條應(yīng)力σ矩陣是轉(zhuǎn)子軸承位置位移響應(yīng)q矩陣及位移角度ψ的函數(shù)，描述為

圖4 彈性支承結(jié)構(gòu)

求解在一個轉(zhuǎn)子旋轉(zhuǎn)周期內(nèi)彈性支承應(yīng)力矩陣的有效值得到σeq，同樣提取軸承位置位移響應(yīng)矩陣有效值qeq，則式(4)式可變?yōu)橐粋€與位移角度無關(guān)的函數(shù)表達(dá)式

矩陣σeq為參數(shù)識別的目標(biāo)函數(shù)變量。

至此，參數(shù)識別中的自變量及目標(biāo)函數(shù)已經(jīng)明確。

2 轉(zhuǎn)子系統(tǒng)參數(shù)識別

采用建立起的高速轉(zhuǎn)子ANSYS有限元模型，基于ANSYS開展支承參數(shù)識別，并以ISIGHT為平臺，開展發(fā)動機(jī)轉(zhuǎn)子不平衡量識別參數(shù)識別。

2.1 支承參數(shù)識別

當(dāng)擠壓油膜阻尼器結(jié)構(gòu)確定后，其Cs和Ks是供油壓力、滑油粘度、轉(zhuǎn)子工作轉(zhuǎn)速、和偏心率函數(shù)。一般情況下，供油壓力和滑油粘度是確定的，工作狀態(tài)下的Cs和Ks僅是轉(zhuǎn)子工作轉(zhuǎn)速及偏心率的函數(shù)。基于ANSYS非線性支承參數(shù)轉(zhuǎn)子不平衡響應(yīng)計算流程見圖5，通過不斷“修正”支承剛度、阻尼使轉(zhuǎn)子不平衡響應(yīng)達(dá)到一個穩(wěn)定值。由于試驗數(shù)據(jù)的離散性，識別過程中基于實測擠壓油膜阻尼器Cs和Ks數(shù)據(jù)采用線性插值獲取。每次迭代得到的擠壓油膜阻尼器剛度阻尼帶入到動力學(xué)有限元模型中即可計算得到支撐位置的油膜力，代入到（1）式即可實現(xiàn)對動力學(xué)方程的求解。

圖5 支承參數(shù)識別流程

2.2 不平衡量識別

轉(zhuǎn)子不平衡量識別過程可歸結(jié)為優(yōu)化問題。基于公式（4），獲取支承位置位移響應(yīng)有效值，化彈支應(yīng)力實驗數(shù)據(jù)為支承位置位移響應(yīng)數(shù)據(jù)。參數(shù)識別過程的數(shù)學(xué)模型可描述為

數(shù)學(xué)模型中：Re表示仿真獲取的彈性支承位移響應(yīng)有效值向量；Ex表示試驗中支承位置位移響應(yīng)有效值向量；V表示仿真值向量與計算值向量之差的2范數(shù)；F表示待識別的不平衡量；f表示變量上下邊界；φ為不平衡量相對相角；?表示相對相角上下邊界，一般情況下取0～2π。

由于轉(zhuǎn)子支承參數(shù)的高度非線性，本文選用Pointer優(yōu)化器，該求解器能夠高度非線性和非連續(xù)性設(shè)計空間問題，通過自動捕捉設(shè)計空間的信息，然后靈活的組合算法形成一個參數(shù)識別策略。Pointer優(yōu)化器包括的算法有：是線性單純形法（linear simplex），序列二次規(guī)劃法（sequential quadratic programming），最速下降法（downhill simple）和遺傳算法（genetic algorithms）。

3 算例研究

以某渦軸發(fā)動機(jī)高速轉(zhuǎn)子為不平衡量識別對象，該轉(zhuǎn)子為3A+1C壓氣機(jī)結(jié)構(gòu)，一級燃?xì)鉁u輪，轉(zhuǎn)子支點(diǎn)跨距為495mm，轉(zhuǎn)子總質(zhì)量為20.3030 kg，軸慣性矩Ixc=62848.30 kg·mm2，直徑慣性矩 Iyc=429375.08 kg·mm2，Izc=429375.08 kg ·mm2?；贏NSYS建立該轉(zhuǎn)子梁單元有限元模型，見圖6。

圖6 高速轉(zhuǎn)子梁單元模型

通過仿真分析發(fā)現(xiàn)，工作過程中，彈性支承處于彈性變形階段，其應(yīng)力有效值σeq和軸承位置位移響應(yīng)qeq呈線性關(guān)系。建立起前后支承彈性支承有限元分析，計算得到σeq=h(qeq)函數(shù)關(guān)系

高速轉(zhuǎn)子前后彈性支承剛度分別為1.72×107N/m和2.15×107N/m，擠壓油膜阻尼器參數(shù)見表1，其中滑油滑油密度ρ=904 kg/m3，運(yùn)動粘度υ=4.2×10-6m2/s。

表1 轉(zhuǎn)子前后支承參數(shù)

以轉(zhuǎn)子一級壓氣機(jī)盤及渦輪盤為等效不平衡量識別位置，式（2）中轉(zhuǎn)子基頻分量，前后平衡截面等效不平衡力分別為：。本算例中待識別參數(shù)為一級壓氣機(jī)盤處等效不平衡量渦輪盤處等效不平衡量和相對相角

以該轉(zhuǎn)子在轉(zhuǎn)子試驗臺上試驗數(shù)據(jù)為目標(biāo)，彈性支承水平與垂直肋條上應(yīng)力值相當(dāng)，提取轉(zhuǎn)速 32900 r/min、36000 r/min、38000 r/min、41600 r/min和42500 r/min轉(zhuǎn)速下轉(zhuǎn)子彈性支承水平肋條上應(yīng)力有效值見表2?；谑剑?）轉(zhuǎn)換得到圖7，該轉(zhuǎn)子彈性支承應(yīng)力隨轉(zhuǎn)速先增加，而后呈減小趨勢。高速轉(zhuǎn)子為典型的剛性轉(zhuǎn)子，轉(zhuǎn)子工作轉(zhuǎn)速位于二階擺動到三階彎曲轉(zhuǎn)速之間，轉(zhuǎn)子工作轉(zhuǎn)速、不平衡量和油膜力是轉(zhuǎn)子臺上彈性支承應(yīng)力發(fā)生變化的主要影響因素。

表2 彈性支承水平肋條應(yīng)力

圖7 軸承位移響應(yīng)有效值

本算例采用的仿真計算機(jī)為Intel i5四核處理器，主頻3.4 GHz，8 GB內(nèi)存容量。

在識別程序中輸入轉(zhuǎn)速、前后軸承位移響應(yīng)有效值，以F1、F2和φ為待識別參量，開展不平衡量識別研究。為減小篇幅，本文只列出轉(zhuǎn)速32900 r/min下轉(zhuǎn)子不平衡量識別過程，在參數(shù)識別程序中分別輸入前、后軸承位移響應(yīng)。以前、后支承位移響應(yīng)為目標(biāo)的識別過程見圖8～圖9，每個迭代步中嵌套有圖5中支承參數(shù)識別流程，不同平衡量下識別的油膜阻尼、剛度也存在差異性，整個不平衡量及支承參數(shù)識別過程約7 min，其他轉(zhuǎn)速下識別時間為5～10 min。

圖8 前支承位移響

圖9 后支承位移響應(yīng)

表3 彈不平衡量識別結(jié)果

?

輸入其他轉(zhuǎn)速及對應(yīng)的目標(biāo)值，識別結(jié)果見表3，繪制工作轉(zhuǎn)速內(nèi)轉(zhuǎn)子不平衡量變化曲線圖，圖10表明，轉(zhuǎn)子不平衡量隨高速轉(zhuǎn)子轉(zhuǎn)速上升呈現(xiàn)放大的趨勢，高速轉(zhuǎn)子前、后不平衡量相對相角約180度。

圖10 全工作轉(zhuǎn)速內(nèi)不平衡量變化

通常情況下，該類轉(zhuǎn)子裝機(jī)前不平衡量會控制在一定標(biāo)準(zhǔn)范圍內(nèi)，經(jīng)過一輪試車之后，對運(yùn)轉(zhuǎn)后的轉(zhuǎn)子做不平衡量測試，不平衡量往往會出現(xiàn)不同程度的放大，根據(jù)筆者的工程經(jīng)驗，轉(zhuǎn)子不平衡量最大可放大數(shù)十倍。通過本算例的研究發(fā)現(xiàn)，工作中的轉(zhuǎn)子不平衡量比裝機(jī)前需用不平衡量放大了10倍到20倍左右。

4 結(jié)論

本文建立起渦軸發(fā)動機(jī)高速兩支點(diǎn)轉(zhuǎn)子有限元動力學(xué)模型，選取轉(zhuǎn)子動平衡位置處等效不平衡量及相對相角為待識別變量，面向轉(zhuǎn)子彈性支承肋條應(yīng)力的工頻信號，基于ANSY編寫擠壓油膜阻尼器阻尼、剛度識別程序，基于ISIGHT軟件采用Pointer求優(yōu)化器完成渦軸發(fā)動機(jī)高速兩支點(diǎn)轉(zhuǎn)子不平衡量的識別，研究發(fā)現(xiàn)：

1)本文所提出的不平衡量識別方法識別效率較高，在文中的算例中，單個工作轉(zhuǎn)速下不平衡量識別的時間為5～10min，該方法在工程設(shè)計中具有較強(qiáng)的可操作性；

2)識別結(jié)果表明工作中的轉(zhuǎn)子不平衡量會發(fā)生不同程度的放大，本文的算例研究發(fā)現(xiàn)，該方法識別的結(jié)果能反映轉(zhuǎn)子在工作轉(zhuǎn)速內(nèi)不平衡量變化趨勢。

本文所識別得到的結(jié)果可作為整機(jī)振動仿真分析的輸入條件，為發(fā)動機(jī)整機(jī)不平衡響應(yīng)研究中的不平衡量的施加提供參考。