張 珺， 李立州， 王沐晨， 路 寬， 原梅妮

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

1 引 言

在航空發(fā)動機(jī)流道中存在多排的葉型、支柱和凸肩等結(jié)構(gòu)。這些結(jié)構(gòu)相互影響，導(dǎo)致某一個流域的流動特性受到其他區(qū)域影響，形成上下游干涉問題。

通常上下游干涉問題的分析方法是對整個流道進(jìn)行數(shù)值模擬[1]。然而，這種方法的代價通常十分高昂[2]，尤其是在需要反復(fù)迭代的某一個葉片優(yōu)化設(shè)計和可靠性分析過程中，這種數(shù)值模擬的計算代價甚至讓人對優(yōu)化設(shè)計和可靠性分析失去興趣。

氣動力降階模型方法ROM(Reduced Order Model)通過簡單的數(shù)學(xué)模型表征線性和非線性氣動力系統(tǒng)的主要特性[3]，提高機(jī)翼和葉片氣動彈性分析[4-12]和形狀優(yōu)化[13-18]的速度，因而近年來得到了廣泛的研究。文獻(xiàn)[12-14]回顧了氣動力降階模型的相關(guān)工作。文獻(xiàn)[16,17,19-21]建議用氣動力降階模型部分或完全地代替氣動模擬。文獻(xiàn)[19]用正交分解法POD討論了翼型反設(shè)計問題。文獻(xiàn)[16,17,20]提出了仿真模型與POD降階模型混合的方法，用仿真模型精確地模擬近壁面的非線性流動，用降階模型模擬平緩的遠(yuǎn)場。

本文針對航空發(fā)動機(jī)上下游干涉的問題，提出了降階模型與全階模型耦合的邊界替代方法。該方法建立了研究流域與相鄰流域交界面的流動特性降階模型，并將該降階模型與研究流域的仿真模型耦合，由此將存在上下游干涉的流動模擬轉(zhuǎn)化為單個流域的流動模擬，同時考慮上下游的耦合效應(yīng)。該方法非常適合航空發(fā)動機(jī)流道中葉片的優(yōu)化設(shè)計和可靠性分析。

2 上下游邊界替代降階模型

圖1給出了上下游邊界替代降階模型方法的主要思想。圖1的流道由三排結(jié)構(gòu)組成。對第二排葉片的流域進(jìn)行分析時，該區(qū)域與其上游流域(第一排葉片)的交界面是上游的出口，也是其進(jìn)口；其與下游流域(第三排葉片)的交界面是下游的進(jìn)口，也是出口。在對第二排葉片流域進(jìn)行模擬時，可以建立該流域與上游流域和下游流域交界面的流動特性降階模型(第一排和第三排對第二排影響的邊界代替降階模型)。根據(jù)壓力和流量在界面守恒的條件，將該流域的仿真模型和交界面降階模型耦合，可以建立僅包含該流域的仿真模型，同時可以考慮上下游影響的第二排的流場仿真模型。

用到的界面降階模型從上下游流域已知的流動分析結(jié)果中辨識得到，包含了上下游流域內(nèi)結(jié)構(gòu)特征對界面流動特性的影響，具體來說，就是界面降階模型能夠在給定的流量下計算出和上下游流域仿真模型相同的界面壓力，或者在給定的界面壓力下計算出和上下游流域仿真模型相同的界面流量。通過研究流域的仿真模型和降階模型耦合，就可以表征整個流道的結(jié)構(gòu)特征對研究流域的影響。

需要指出的是，本文的上下游區(qū)域是相對于研究對象的。圖1中，如果選定第二排葉片所在的流域?yàn)檠芯繉ο?，則第一排葉片的流域是上游，第三排葉片的流域?yàn)橄掠危蝗绻x定第一排葉片所在的流域?yàn)檠芯繉ο?，第二排和第三排葉片的流域一起組成下游。

圖1 航空發(fā)動機(jī)流道模擬邊界替代降階模型Fig.1 Boundary ROM for aeroengine flow passage

采用Volterra級數(shù)法建立上述邊界替代降階模型。當(dāng)然該降階模型也可以采用POD和神經(jīng)網(wǎng)絡(luò)等方法建立。為簡單起見，采用二維流動問題詳述Volterra級數(shù)方法，三維問題可參照施行。

對發(fā)動機(jī)流道中的任意一個區(qū)域進(jìn)行研究，可以用交界面將該流域與其他區(qū)域分開。該研究流域與任意一個相鄰流域交界面的壓力和流量可以表示為

{p(y,t)}=Ψ{m(y,t)}

(1)

式中Ψ為交界面的流動模型，m(y,t)為界面流量分布，p(y,t)為界面壓力分布，y為坐標(biāo)，t為時間。

根據(jù)多輸入多輸出Volterra級數(shù)理論[22，23]，式(1)可以用離散的一階卷積Volterra級數(shù)表示為

(i,j=1,2,…,r) (2)

Bj i的本質(zhì)是i節(jié)點(diǎn)單位流量變化引起j節(jié)點(diǎn)壓力的脈沖響應(yīng)，通常從仿真或?qū)嶒?yàn)數(shù)據(jù)中辨識得到[8,10-12]。本文從包含界面Ψ的相鄰流域的氣動分析結(jié)果中辨識得到。由于這些氣動分析結(jié)果包含了相鄰流域內(nèi)部結(jié)構(gòu)對界面的影響，由此辨識出的Bj i也包含了相應(yīng)的信息。當(dāng)降階模型和研究區(qū)域的仿真模型耦合時界面流動特性就包含了所有區(qū)域?qū)缑娴挠绊憽?/p>

如果只考慮靜態(tài)響應(yīng)，式(2)可簡化為

(i,j=1,2,…,r)(3)

由于式(2,3)只保留了一階項(xiàng)，忽略了高階項(xiàng)，因此是一種降階模型。

對于選定的研究區(qū)域，建立仿真模型。用上述方法建立其上下游界面的邊界替代降階模型。將仿真模型與降階模型耦合，可以建立選定區(qū)域的混合模型。這樣，在選定區(qū)域進(jìn)行流動分析時，可以大幅減少計算量，非常適合發(fā)動機(jī)葉片的優(yōu)化設(shè)計和可靠性分析。

3 邊界替代降階模型的算例驗(yàn)證

3.1 算例說明

用圖2的算例驗(yàn)證本文方法。整個流域的入口總壓為120 kPa，整個流域的出口壓力為 101.325 kPa，入口進(jìn)氣角和出口排氣角都為0°。研究圖2中翼型區(qū)域的流場特性。翼型區(qū)域和下游擾流柱區(qū)域的交界面取在x/Lc=1.42，下游擾流柱間距 0.5Lc。用邊界替代降階模型表征下游擾流柱區(qū)域?qū)σ硇偷挠绊憽ＤM采用雙精度理想氣體、Spallart-Allmaras模型、空氣和周期條件。收斂精度 1e -6。

圖2 計算域Fig.2 Flow domain

3.2 界面壓力與流量關(guān)系的降階模型

以下游擾流柱區(qū)域的進(jìn)口壓力和流量的關(guān)系建立界面降階模型。選定下游擾流柱流域進(jìn)口各點(diǎn)的基準(zhǔn)流量通量為190 kg/(m2·s)。對該基準(zhǔn)流量下的流動特性進(jìn)行模擬，流場的速度如圖3(a)所示。提取此狀態(tài)下的進(jìn)口壓力，將該壓力作為基準(zhǔn)壓力(圖4中baseline)。

圖3 不同位置流量擾動后下游流域的速度分布(部分?jǐn)?shù)據(jù))Fig.3 Velocity of the downstream domain with disturbances at several locations (part of data)

采用逐點(diǎn)擾動進(jìn)口流量的方法獲得不同流量擾動下的進(jìn)口壓力分布。本文流量的擾動量為 -30 kg/(m2·s)，即在擾動點(diǎn)處流量設(shè)置為190-30=160 kg/(m2·s)。圖3(b～d)給出了流量擾動下的流場速度。圖4給出了流量擾動下的進(jìn)口壓力分布。該算例的下游流域仿真模型進(jìn)口共有100個網(wǎng)格，因此要想完整獲得下游流域進(jìn)口流量和壓力的關(guān)系需要100組數(shù)據(jù)?？紤]到下游流域具有對稱性，因此僅需要25組數(shù)據(jù)就足夠降階模型辨識。另外，本文針對降階模型辨識的需要專門進(jìn)行了仿真，在實(shí)際問題中這些數(shù)據(jù)完全可以從已有的設(shè)計過程和試驗(yàn)結(jié)果中獲得。

圖4 不同位置流量擾動后下游流域界面的壓力分布(部分?jǐn)?shù)據(jù))

3.3 翼型仿真模型與下游界面降階模型的耦合

將翼型流域的仿真模型與下游界面降階模型耦合，形成混合模型。由該混合模型計算界面壓力和速度，結(jié)果如圖5和圖6的mixed所示。為了驗(yàn)證本文方法，還給出了翼型和下游流域都采用CFD模型的結(jié)果(圖5和圖6中CFD)，這里簡稱CFD模型結(jié)果?？梢钥闯?，混合模型結(jié)果和CFD模型結(jié)果一致，這說明本文的方法是可行的。

圖5 混合模型和CFD模型的界面壓力和速度Fig.5 Interface pressures and velocities of the mixed and CFD models

圖6 混合模型和CFD模型的翼型區(qū)域壓力(單位：Pa)Fig.6 Airfoil pressures of the mixed and CFD models (unit:Pa)

4 基于混合模型的翼型優(yōu)化

用混合模型優(yōu)化翼型，目標(biāo)是增大翼型升力。采用Hicks-Henne方法進(jìn)行翼型參數(shù)化[24,25]。本文只對上表面的6個形狀參數(shù)a1至a6進(jìn)行了優(yōu)化，其取值范圍為[-0.01,0.01]。采用序列二次規(guī)劃法SQP作為優(yōu)化算法。優(yōu)化共迭代230次，得到的翼型升阻比為0.07186/0.02424。

為了驗(yàn)證本文方法，還給出了采用CFD模型的優(yōu)化結(jié)果。優(yōu)化時采用完全相同的參數(shù)化方法、取值范圍、優(yōu)化目標(biāo)和優(yōu)化算法。采用CFD模型也迭代230次，得到翼型升阻比為0.07176/0.02417。

圖7和圖8比較了采用混合模型和CFD模型優(yōu)化得到的最優(yōu)翼型、壓力系數(shù)和壓力分布。結(jié)果表明，采用混合模型和CFD模型的優(yōu)化結(jié)果一致。因此，基于邊界替代降階模型法是可行的。

在優(yōu)化過程中每個混合模型需要13分鐘(2970次迭代，0.2626 s/次)收斂，每個CFD模型需要29 min(6262次迭代，0.2779 s/次)收斂，這表明降階模型提高了計算效率?？紤]到在算例中只替代了一個很小的流域，因此單次迭代計算效率提高不明顯。如果替代的區(qū)域很大，如發(fā)動機(jī)多級葉片的流域，則效率提升會非常顯著。從每個模型收斂所用的迭代次數(shù)來看，混合模型比CFD模型收斂速度快。分析可知，混合模型的網(wǎng)格比CFD模型少，因此，消除所有網(wǎng)格計算誤差的迭代次數(shù)少。

圖7 混合模型和CFD模型的優(yōu)化翼型和壓力系數(shù)Fig.7 Optimal airfoil and its pressures of the mixed and CFD models

圖8 優(yōu)化翼型的混合模型和CFD的壓力云圖(單位：Pa)Fig.8 Pressures contour of the optimal airfoil by the mixed and CFD model (unit:Pa)

5 總 結(jié)

為了降低航空發(fā)動機(jī)流道某個區(qū)域流動特性分析的計算量，提出了一種邊界替代降階模型方法。該方法用界面將研究區(qū)域與上下游區(qū)域分開，從上下游區(qū)域的模擬結(jié)果中辨識出界面的壓力和流量關(guān)系的邊界替代降階模型，將流域仿真模型與邊界替代降階模型耦合，舍棄了流道內(nèi)大部分區(qū)域的仿真模型，同時保留其對研究區(qū)域的影響。

用一排翼型和一排下游擾流柱的例子驗(yàn)證了本文的方法，結(jié)果表明,(1) 采用本文方法與完全采用CFD得到的翼型流域的壓力和流速完全一致； (2) 采用本文方法將減小計算量； (3) 采用本文方法更容易收斂?？偟膩碚f，本文方法可以極大地提高發(fā)動機(jī)流道氣動優(yōu)化和可靠性分析的效率。