氣動(dòng)力降階模型的優(yōu)化識(shí)別方法

張 珺，李立州，羅 驍，張新燕

(1 太原學(xué)院 數(shù)學(xué)系， 太原 030001; 2 中北大學(xué) 機(jī)電工程學(xué)院， 太原 030051)

基于CFD和CSD耦合的數(shù)值求解方法是研究葉片和機(jī)翼氣動(dòng)彈性振動(dòng)的主要方法[1-5]，但其計(jì)算效率極低[6-16]，工程應(yīng)用不便。氣動(dòng)力降階模型(aerodynamic reduced order model,aerodynamic ROM)是描述擾動(dòng)對(duì)葉片和機(jī)翼氣動(dòng)力影響的簡(jiǎn)化數(shù)學(xué)模型[6]。通過氣動(dòng)力降階模型研究葉片和機(jī)翼的氣動(dòng)彈性振動(dòng)已成為重要方向。Silva[6]提出了線性和非線性氣動(dòng)力Volterra級(jí)數(shù)降階模型；Su等[7]用Volterra級(jí)數(shù)降階模型研究了葉柵的顫振，認(rèn)為降階模型的計(jì)算效率極高；Liou[8]用Volterra級(jí)數(shù)法研究了跨音速條件下葉片的顫振；He[9]用諧波平衡法(harmonic balance)研究了含分離流動(dòng)的葉片流場(chǎng)； Ekici等[10]用諧波平衡法研究了葉片顫振，認(rèn)為該方法與勢(shì)流理論符合較好；Ashcroft等[11]用諧波平衡法研究了亞音速和跨音速二維壓氣機(jī)葉柵的顫振特性，發(fā)現(xiàn)該方法可以準(zhǔn)確地預(yù)測(cè)顫振；Balajewicz[12]研究了多輸入多輸出Volterra級(jí)數(shù)的氣動(dòng)力降階模型；Skujins[13]用降階模型研究了高超聲速飛行器的動(dòng)態(tài)氣動(dòng)力；Jirasek[14]用降階模型研究了X-31的動(dòng)態(tài)氣動(dòng)力；Wang[15]建立了基于Tchebyshev多項(xiàng)式的二階Volterra氣動(dòng)力降階模型；Milansese[16]用連續(xù)的脈沖識(shí)別了Volterra降階模型的核函數(shù)。

目前，氣動(dòng)力降階模型核函數(shù)都采用脈沖法或階躍法識(shí)別。對(duì)于多輸入多輸出的問題，需要多次的CFD計(jì)算或者試驗(yàn)。為減少氣動(dòng)力降階模型識(shí)別時(shí)CFD模擬或氣動(dòng)力試驗(yàn)的次數(shù)，本文提出為每一個(gè)核函數(shù)建立一個(gè)擬合曲線，通過優(yōu)化曲線的參數(shù)，使氣動(dòng)力降階模型的氣動(dòng)力輸出和CFD模型試驗(yàn)結(jié)果一致，以此辨識(shí)多輸入多輸出氣動(dòng)力降階模型的核函數(shù)。

1 氣動(dòng)力降階模型的優(yōu)化識(shí)別方法

1.1 時(shí)變系統(tǒng)的Volterra級(jí)數(shù)降階模型

對(duì)于均勻采樣的時(shí)變系統(tǒng)，其離散的Volterra 級(jí)數(shù)模型可以寫為[10]：

(1)

1.2 基于Volterra級(jí)數(shù)的葉片氣動(dòng)力降階模型

本文以二維葉片位移擾動(dòng)引起氣動(dòng)力變化為例，介紹葉片氣動(dòng)力降階模型核函數(shù)的優(yōu)化識(shí)別方法。葉片模型如圖1所示，葉片的扭轉(zhuǎn)和彎曲引起葉片氣動(dòng)升力和升力矩的變化，該系統(tǒng)模型可以用{CL(t),CM(t)}=Ψ{v(t),α(t)}表示。α和v分別為葉片的扭轉(zhuǎn)和彎曲位移，CL和CM為升力和力矩。

該系統(tǒng)的一階Volterra級(jí)數(shù)氣動(dòng)力降階模型可表示為：

(2)

1.3 核函數(shù)的優(yōu)化識(shí)別方法

現(xiàn)有氣動(dòng)力降階模型的核函數(shù)采用階躍和脈沖信號(hào)識(shí)別。用階躍信號(hào)辨識(shí)出的核函數(shù)比較穩(wěn)定且能部分地反映系統(tǒng)的非線性特性，因此在實(shí)踐中應(yīng)用較為廣泛。用階躍信號(hào)識(shí)別核函數(shù)的步驟是：首先獲得葉片的穩(wěn)態(tài)流場(chǎng)和氣動(dòng)力，然后在穩(wěn)態(tài)流場(chǎng)上增加單位階躍的擾動(dòng)，獲得在該階躍擾動(dòng)下葉片氣動(dòng)力的響應(yīng)，最后用如下公式計(jì)算氣動(dòng)力降階模型的核函數(shù)：

(3)

① 用CFD或者試驗(yàn)獲得已知擾動(dòng)下葉片的氣動(dòng)力響應(yīng)；

② 確定氣動(dòng)力降階模型的輸入輸出變量和核函數(shù)的階次；

③ 為每一個(gè)核函數(shù)建立一條擬合曲線；

④ 用核函數(shù)的擬合曲線代替核函數(shù)，構(gòu)造氣動(dòng)力降階模型；

⑤ 優(yōu)化所有的核函數(shù)擬合曲線的參數(shù)，使氣動(dòng)力降階模型的輸出與CFD或試驗(yàn)得到的氣動(dòng)力一致。

2 算例驗(yàn)證

為驗(yàn)證本文提出的方法，以圖2葉片為對(duì)象，建立氣動(dòng)力降階模型。取中間葉片的彎曲位移和俯仰角為輸入，取中間葉片的氣動(dòng)力為輸出。采用CFD模型獲得輸入輸出之間的關(guān)系。CFD模擬采用fluent，空氣，理想氣體，Spallart-Allmaras粘性模型，無滑移壁面，葉片流場(chǎng)進(jìn)、出口見圖2。出口壓力101 325 Pa，進(jìn)口總壓120 500 Pa，溫度為300 K，進(jìn)氣角30°。求解時(shí)，確定先求得穩(wěn)態(tài)流場(chǎng)，求得葉片的升力h0為-797 N，力矩u0為 -12.22 N·m。在穩(wěn)態(tài)流場(chǎng)的基礎(chǔ)上，用udf程序調(diào)整不同時(shí)刻中間葉片的彎曲位移和俯仰角，進(jìn)而獲得位移擾動(dòng)條件下葉片的氣動(dòng)力，用于核函數(shù)的辨識(shí)。

根據(jù)模型識(shí)別理論，系統(tǒng)的激勵(lì)信號(hào)必須充分激勵(lì)系統(tǒng)所有模態(tài)。因此，在本文中噪聲信號(hào)被用做葉片的擾動(dòng)。葉片的彎曲位移擾動(dòng)采用v=0.000 1×(random(t)-0.5)，葉片的俯仰角擾動(dòng)采用α=0.1×(random(t)-0.5)。用udf將該信號(hào)加載到葉片fluent模型上，得到葉片的氣動(dòng)力，見圖3和圖4。

為了驗(yàn)證核函數(shù)的精度，圖5也給出了直接用階躍信號(hào)識(shí)別的核函數(shù)。從結(jié)果來看：用優(yōu)化方法識(shí)別的核函數(shù)和用階躍方法識(shí)別的核函數(shù)完全相同。這說明了本文提出方法是可行的。另外，考慮到本文的方法可以一次性識(shí)別多輸入多輸出系統(tǒng)所有的核函數(shù)并且具有平均和容錯(cuò)能力，因此本文的方法具有更大的優(yōu)勢(shì)。

將優(yōu)化得到的核函數(shù)和階躍信號(hào)得到的核函數(shù)代入葉片氣動(dòng)力降階模型，計(jì)算噪聲信號(hào)下葉片的氣動(dòng)力，結(jié)果見圖3和圖4。為方便比較，用CFD模型計(jì)算的葉片氣動(dòng)力也顯示在圖3和圖4中。從圖3和圖4 可以看出，降階模型和CFD模型得到的葉片氣動(dòng)力結(jié)果完全一致。這進(jìn)一步說明了本文提出方法是可行的。

3 結(jié)論

現(xiàn)有研究通過階躍響應(yīng)和脈沖響應(yīng)得到氣動(dòng)力降階模型的核函數(shù)。對(duì)于多輸入多輸出問題，需要多次進(jìn)行CFD模擬或試驗(yàn)，計(jì)算或試驗(yàn)的工作量巨大。為解決這一問題，本文提出采用優(yōu)化方法識(shí)別氣動(dòng)力降階模型核函數(shù)，其基本思想是：為每一個(gè)核函數(shù)建立一個(gè)擬合曲線；用優(yōu)化方法調(diào)整擬合曲線的參數(shù)，使氣動(dòng)力降階模型的結(jié)果和CFD或試驗(yàn)結(jié)果一致，一次性獲得所有的核函數(shù)。

二維葉片氣動(dòng)力降階模型核函數(shù)的辨識(shí)結(jié)果表明：優(yōu)化得到的核函數(shù)與用階躍方法識(shí)別的結(jié)果完全相同，說明本文提出方法是可行的。另外，考慮到本文的方法可以一次性識(shí)別多輸入多輸出系統(tǒng)所有的核函數(shù)，并且具有一定的容錯(cuò)能力，因此本文提出方法具有更大優(yōu)勢(shì)。

參考文獻(xiàn)：

[1] AMOO L M.On the design and structural analysis of jet engine fan blade structures[J].Progress in Aerospace Sciences,2013,60:1-11.

[2] VIKAS R,VIPAN K,LOVELEEN K B.Dynamic modeling of underplateform damper used in turbomachinery[J].World Academy of Science,Engineering and Technology,2012,61:536-545.

[3] 陳予恕,張華彪.航空發(fā)動(dòng)機(jī)整機(jī)動(dòng)力學(xué)研究進(jìn)展與展望[J].航空學(xué)報(bào),2011,32(8):1371-1391.

[4] HUANG X Q,HE L,BELL D L.Influence of upstream stator on rotor flutter stability in a low pressure steam turbine stage[J].Proceeding of the Institution of Mechanical Engineering,2006,220(1):25-35.

[5] HSU C N.A study on fluid self-excited flutter and forced response of turbomachinery rotor blade[J].Mathematical Problems in Engineering,2014:1-20.

[6] SILVA W A.Reduced order models based on linear and nonlinear aerodynamic impulse responses[R].AIAA 99-1262.

[7] SU D,ZHANG W W,MA M S,et al.An efficient coupled method of cascade flutter investigation based on reduced order model[R].AIAA 2013-3618.

[8] LIOU M S,YAO W G.Flutter analysis for turbomachinery using Volterra series[C].Proceedings of ASME Turbo Expo 2014.June 16-20,2014,Dusseldorf,Germany.GT2014-25474.

[9] HE L.Harmonic solution of unsteady flow around blades with separation[J].AIAA Journal,2008,46(6),1299-1307.

[10] EKICI K,KIELB R E,HALL K C.The effect of aerodynamic asymmetries on turbomachinery flutter[J].Journal of Fluids and Structures,2013,36:1-17.

[11] ASHCROFT G,FREY C,KERSHEN H P.On the development of a harmonic balance method for aeroelastic analysis[C].11thWorld Congress on Computational Mechanics(WCCM XI),Barcelona,Spain,2014.

[12] BALAJEWICZ M,NITZSCHE F,FESZTY D.Application of multi-input Volterra theory to nonlinear multi-degree of freedom aerodynamic systems[J].AIAA Journal,2010,48(1):56-62.

[13] SKUJINS T,CESNIK C E S.Reduced-order modeling of hypersonic vehicle unsteady aerodynamics[R].AIAA 2010-8217.

[14] JIRASEK A,CUMMINGS R M.Application of Volterra functions to X-31 aircraft model motion[R].AIAA 2009-3629.

[15] WANG Y H,HAN J L.Approach to identification of a second-order Volterra kernel of nonlinear systems by Tchebyshev polynomials method[J].Research Journal of Applied Science,Engineering and Technology,2013,5(20):4950-4955.

[16] MILANSESE A,MARZOCCA P.Volterra kernels identification using continuous time impulses applied to a nonlinear aeroelastic problem[R].AIAA 2009-2321.