雷玉昌，張登成，張艷華

(1.空軍工程大學(xué) 研究生院，西安 710051；2.中國人民解放軍95034部隊(duì)，百色 533600；3.空軍工程大學(xué) 航空工程學(xué)院，西安 710038）

在航空技術(shù)應(yīng)用中，使用環(huán)量控制技術(shù)作為主動流動控制的解決方案越來越引起人們的注意［1-2］。環(huán)量控制技術(shù)是指在翼型后緣表面開縫以形成沿著物面切向的射流，用以增加沿著翼型表面的環(huán)量，進(jìn)而增加升力［3-4］。風(fēng)洞試驗(yàn)和數(shù)值仿真結(jié)果表明，環(huán)量控制技術(shù)能夠大幅度提高翼型升力，并在高升力條件下改善升阻比［5-6］。環(huán)量控制技術(shù)最初就是因?yàn)槠渥吭降母呱π阅芏痍P(guān)注，經(jīng)過八十多年的發(fā)展，該技術(shù)逐漸在無尾飛行器的氣動力控制、風(fēng)力渦輪機(jī)的載荷控制、降低噪聲、降低飛行器結(jié)構(gòu)質(zhì)量及提高隱身性能等方面表現(xiàn)出優(yōu)越的性能［7-10］。特別是在無尾飛行器的控制方面，逐漸表現(xiàn)出了替代傳統(tǒng)控制舵面的發(fā)展趨勢。英國BAE系統(tǒng)公司先后試飛了“惡魔(Demon）”、“Magma”等技術(shù)驗(yàn)證無人機(jī)，利用環(huán)量控制技術(shù)進(jìn)行了滾轉(zhuǎn)操作，驗(yàn)證了環(huán)量控制射流替代傳統(tǒng)舵面的可行性［11-12］。雖然已經(jīng)證明了環(huán)量控制技術(shù)的相關(guān)優(yōu)越性，但要實(shí)現(xiàn)可用的環(huán)量控制射流，還有很多關(guān)鍵理論和技術(shù)問題需要研究。其中考慮迎角改變時(shí)環(huán)量控制射流對于非定常氣動力的瞬態(tài)反應(yīng)是必不可少的。因?yàn)橐环矫姝h(huán)量控制射流并不能像傳統(tǒng)舵面一樣，在迎角改變時(shí)與翼型保持同步偏轉(zhuǎn)，而是存在時(shí)間遲滯效應(yīng)；另一方面高動量系數(shù)情況下會導(dǎo)致機(jī)翼提前失速，導(dǎo)致翼型在中小迎角下就會出現(xiàn)失速分離情況。因此，傳統(tǒng)飛行器的非定常氣動力模型將無法適用于環(huán)量控制飛行器。而飛機(jī)過失速機(jī)動時(shí)會產(chǎn)生流動分離和漩渦破碎，使得氣動力和氣動力矩呈現(xiàn)高度的非線性和非定常特性［13-14］。建立環(huán)量控制翼型在機(jī)動過程中的非定常氣動力模型，對于準(zhǔn)確理解射流作用機(jī)理及進(jìn)行飛行控制系統(tǒng)設(shè)計(jì)具有重要意義。

目前，關(guān)于環(huán)量控制技術(shù)的研究中并未建立起高精度的非定常氣動力模型。Loth和Englar等［15-16］給出了基于定常實(shí)驗(yàn)數(shù)據(jù)的曲線擬合公式，其中，文獻(xiàn)［15］給出的經(jīng)驗(yàn)公式在較低迎角、單射流情況下與可用實(shí)驗(yàn)數(shù)據(jù)有較好的一致性。Hoholis［17］針對二維環(huán)量控制翼型的氣動參數(shù)進(jìn)行了曲線擬合，并添加了穩(wěn)定性導(dǎo)數(shù)以凸顯其非定常效應(yīng)，結(jié)果表明，該方法對升力和力矩系數(shù)預(yù)測結(jié)果較好，對阻力系數(shù)預(yù)測較差，而且并不能準(zhǔn)確預(yù)測非定常效應(yīng)的氣動力遲滯效應(yīng)。Krukow和Dinkler［18］在進(jìn)行環(huán)量控制機(jī)翼的氣動彈性研究中利用高階多項(xiàng)式擬合得到了定常氣動力模型。Semaan等［19］使用LASSO方法辨識了周期性射流和脈沖型射流疊加后對機(jī)翼升力系數(shù)的增量影響。綜上所述，目前對于環(huán)量控制飛行器氣動力模型的建立多集中在定常模式下，對非定常氣動力模型的研究較為匱乏。而傳統(tǒng)飛行器非定常氣動力建模的發(fā)展已經(jīng)歷經(jīng)了幾十年，逐漸發(fā)展形成了氣動導(dǎo)數(shù)模型［20］、狀態(tài)空間模型［21］和微分方程模型［22］等數(shù)學(xué)模型和模糊邏輯模型［23］、神經(jīng)網(wǎng)絡(luò)模型［24］等人工智能模型。這些模型都能在一定程度上描述傳統(tǒng)飛行器在大迎角下的非定常氣動力變化趨勢。

因此，本文借助高精度的數(shù)值仿真技術(shù)，分析環(huán)量控制翼型非定常流動狀態(tài)。基于微分方程模型，添加環(huán)量控制射流參數(shù)源項(xiàng)。通過兩步線性回歸參數(shù)辨識方法建立適用于環(huán)量控制翼型的非定常氣動力數(shù)學(xué)模型，并對高動量系數(shù)下的非線性影響進(jìn)行修正。通過對任意條件下的非定常流動充分建模，用于環(huán)量控制技術(shù)的非定常狀態(tài)研究。

1 數(shù)值計(jì)算方法與驗(yàn)證

1.1 計(jì)算模型和網(wǎng)格

超臨界翼型有著鈍前緣、大厚弦比的特點(diǎn)，內(nèi)部能夠容納環(huán)量控制技術(shù)所需要的供氣機(jī)構(gòu)或管道等，因而成為環(huán)量控制技術(shù)良好的研究翼型。本文采用修形過的超臨界翼型進(jìn)行相關(guān)研究。翼型相對厚度為17%，弦長c＝240 mm。采用文獻(xiàn)［25］中的修形方法，將距后緣25%弦長的部分修形處理，按照后緣半徑r／c＝2%繪制柯恩達(dá)后緣曲面，射流口高度h＝0.001c。圖1為修形過的超臨界翼型示意圖。計(jì)算區(qū)域選取翼型弦長的30倍，網(wǎng)格剖分時(shí)采用結(jié)構(gòu)網(wǎng)格，生成O型網(wǎng)格拓?fù)浣Y(jié)構(gòu)，邊界層內(nèi)第1層的網(wǎng)格高度約為1×10－5m，射流口附近第1層網(wǎng)格高度為2×10－6m，保證第1層高度的y＋均小于1，以滿足黏性底層的計(jì)算要求。網(wǎng)格總數(shù)約為58萬，進(jìn)行非定常計(jì)算時(shí)采用FLUENT網(wǎng)格動態(tài)重構(gòu)技術(shù)，圖2為翼型網(wǎng)格劃分及非定常計(jì)算過程迎角α從－5°變化到5°翼型后緣網(wǎng)格變形情況?？梢钥闯觯W(wǎng)格適應(yīng)情況較好。

圖1 修形后的超臨界翼型Fig.1 Supercritical airfoil after modification

圖2 翼型計(jì)算網(wǎng)格Fig.2 Airfoil’s computational grid

1.2 計(jì)算方法和驗(yàn)證

本文的數(shù)值模擬方法采用二維雷諾平均Navier-Stokes方程，其積分形式為

式中：Q為流動變量；Ω為控制體；FC和FV分別為無黏性通量和黏性通量；n為控制體面的外法線向量。采用k-εSST湍流模型，該模型對于有較大逆壓梯度的邊界層流動、分離預(yù)測性能較好。采用有限體積法離散控制方程，黏性通量采用二階迎風(fēng)格式離散。遠(yuǎn)場邊界為壓力遠(yuǎn)場，壁面邊界采用無滑移壁面條件。

本文采用的數(shù)值仿真軟件為商用軟件FLUENT。采用的求解器為基于Simple方法的Pressure-Based Solver。對求解不可壓縮流動及微可壓流動效果較好。該方法對于非定常問題的流動計(jì)算，主要借助隱式時(shí)間積分方案，在每個(gè)計(jì)算時(shí)間步內(nèi)進(jìn)行迭代，在計(jì)算步之內(nèi)的內(nèi)迭代計(jì)算與普通穩(wěn)態(tài)問題的內(nèi)迭代計(jì)算是相類似的，在得到本時(shí)間步內(nèi)的收斂解之后，轉(zhuǎn)入下一個(gè)時(shí)間步進(jìn)行繼續(xù)計(jì)算。

為了驗(yàn)證本文所采用數(shù)值計(jì)算方法的準(zhǔn)確性，分別對標(biāo)模NACA0012翼型強(qiáng)迫俯仰振動下的非定常氣動力數(shù)據(jù)及環(huán)量控制翼型在射流吹氣下的定常氣動力數(shù)據(jù)進(jìn)行對比驗(yàn)證。

對于標(biāo)模NACA0012翼型而言，計(jì)算過程中強(qiáng)迫翼型繞重心做俯仰振蕩，其迎角變化規(guī)律為

式中：減縮頻率k＝ωc／(2V）∞，c為弦長；V∞為來流速度為無量綱時(shí)間；α0為基準(zhǔn)迎角；αm為振幅。減縮頻率用來描述翼型的振蕩頻率。計(jì)算條件：k＝0.081 4，α0＝0.016°，αm＝2.51°。本文計(jì)算時(shí)采用的物理時(shí)間步長t＝0.01T，內(nèi)迭代步數(shù)為50，其中T表示翼型非定常運(yùn)動周期。為加快收斂，初始解從中心迎角處的定常收斂解開始計(jì)算。將計(jì)算結(jié)果與文獻(xiàn)［26］中Landon的實(shí)驗(yàn)結(jié)果和Uygun的數(shù)值仿真結(jié)果進(jìn)行對比，結(jié)果如圖3所示。圖中：CL為升力系數(shù)，α為迎角。在部分迎角下，仿真結(jié)果與實(shí)驗(yàn)數(shù)據(jù)存在差異，但是在總體上與實(shí)驗(yàn)數(shù)據(jù)是一致的，能夠反映非定常流動時(shí)的氣動力變化趨勢。

圖3 NACA0012翼型非定常升力系數(shù)仿真與實(shí)驗(yàn)對比Fig.3 Comparison of simulation and experimental results of unsteady lift coefficient of NACA0012 airfoil

對于環(huán)量控制翼型而言，動量系數(shù)是影響環(huán)量控制中后緣射流控制效果的重要參數(shù)，升阻特性受動量系數(shù)影響較大，其定義為

式中：m˙為射流出口處的質(zhì)量流量；Vjet為射流速度；ρ∞為來流密度；V∞為來流速度；S為二維翼型單位展長的面積，在數(shù)值上等于弦長c。假定氣流等熵膨脹至射流出口，遠(yuǎn)前方來流靜壓為射流出口靜壓。通過調(diào)節(jié)射流出口速度，調(diào)節(jié)射流出口處的質(zhì)量流量和動量系數(shù)。本文將使用上述參數(shù)作為后緣射流的邊界條件，調(diào)整參數(shù)以進(jìn)行相關(guān)研究。對0°迎角下不同吹氣動量系數(shù)下的環(huán)量控制翼型進(jìn)行數(shù)值模擬。計(jì)算條件為：來流速度為30 m／s，來流溫度為293.15 K，雷諾數(shù)為5×105，動量系數(shù)在0～0.06之間變化。將仿真結(jié)果與文獻(xiàn)［25］中的實(shí)驗(yàn)數(shù)據(jù)進(jìn)行對比，仿真采用穩(wěn)態(tài)方法，結(jié)果如圖4所示，仿真結(jié)果與實(shí)驗(yàn)數(shù)據(jù)基本吻合。綜上所述，本文采用的數(shù)值計(jì)算方法應(yīng)當(dāng)能夠較好地模擬環(huán)量控制翼型的非定常氣動力數(shù)據(jù)。

圖4 環(huán)量控制翼型升力系數(shù)仿真與實(shí)驗(yàn)對比Fig.4 Comparison of simulation and experimental results of circulation control airfoil’s lift coefficient

2 定常氣動力插值

借助數(shù)值仿真或風(fēng)洞試驗(yàn)手段得到的翼型以迎角作為自變量的氣動參數(shù)往往具有極大的離散性和不完整性，因此，在進(jìn)行氣動力和力矩建模時(shí)需要對離散的氣動力數(shù)據(jù)進(jìn)行擬合或插值，使之成為光滑的連續(xù)函數(shù)。本節(jié)基于Kriging模型對定常狀態(tài)下的環(huán)量控制翼型氣動力進(jìn)行插值處理，并與傳統(tǒng)氣動導(dǎo)數(shù)模型的預(yù)測精度進(jìn)行對比。

2.1 K riging模型

Kriging模型最初是由南非采礦工程師Krige在1951年提出的一種空間估計(jì)技術(shù)，最初應(yīng)用于礦床儲量計(jì)算和誤差估計(jì)，經(jīng)過幾十年的發(fā)展，逐漸在地質(zhì)、氣象、航空航天、汽車等領(lǐng)域得到發(fā)展和應(yīng)用［27］。Kriging模型由于其卓越的非線性函數(shù)插值預(yù)測能力和誤差估計(jì)功能，正在受到越來越多的關(guān)注。

設(shè)m1，m2，…，mn為數(shù)量集上的一系列變量值。z(m1），z(m2），…，z(mn）為相應(yīng)的參數(shù)值。則Kriging模型定義變量 m0處的插值結(jié)果z*(m0）為已知樣本值的線性加權(quán)，即

因此，通過求解加權(quán)系數(shù)λi就可以得到數(shù)量集中任意位置的參數(shù)估計(jì)值。關(guān)于Kriging模型的詳細(xì)介紹和推導(dǎo)參見文獻(xiàn)［27］。

影響環(huán)量控制翼型定常氣動力的主要因素包括迎角α和動量系數(shù)Cμ。仿真計(jì)算了α為－5°～15°，Cμ為0～0.06之間的氣動力參數(shù)。

環(huán)量控制翼型的氣動力參數(shù)呈現(xiàn)高度的非線性特征，迎角α和動量系數(shù)Cμ對氣動力的影響高度耦合，翼型失速與動量系數(shù)大小密切相關(guān)，動量系數(shù)越大，翼型失速迎角越小。以其氣動力參數(shù)可作為兩者的空間分布變量。選取樣本數(shù)量集α為－5°～15°，Cμ為0～0.06之間的離散氣動力參數(shù)。對樣本數(shù)量集進(jìn)行歸一化后利用Kriging模型分不同方向進(jìn)行函數(shù)計(jì)算，其預(yù)測結(jié)果和均方根誤差結(jié)果如圖5和圖6所示。

圖5 Kriging模型預(yù)測結(jié)果Fig.5 Prediction result of Kriging model

圖6 Kriging模型預(yù)測誤差Fig.6 Prediction error of Kriging model

從整體上來看，Kriging模型對于環(huán)量控制翼型的定常氣動力插值結(jié)果精度較高，特別是對于采樣點(diǎn)范圍內(nèi)的內(nèi)插結(jié)果誤差值低于1×104。對于邊界值誤差較高，這反映了模型對氣動力參數(shù)的外插能力需要進(jìn)一步提高。Kriging模型能夠?qū)崿F(xiàn)環(huán)量控制翼型氣動力參數(shù)設(shè)計(jì)的光滑連續(xù)插值，且具有較高的預(yù)測精度。

2.2 傳統(tǒng)氣動導(dǎo)數(shù)模型

傳統(tǒng)的氣動導(dǎo)數(shù)模型是基于線性疊加原理，假設(shè)氣動力可由飛行器迎角、舵面偏角、角速度等狀態(tài)變量線性疊加，其本質(zhì)是氣動力的泰勒展開式分解后只保留其一階表達(dá)式，保留的各狀態(tài)變量系數(shù)分別被稱為氣動靜、動導(dǎo)數(shù)。對于定常氣動力而言，僅涉及氣動靜導(dǎo)數(shù)。為了提高氣動導(dǎo)數(shù)模型在大迎角非線性情況下的預(yù)測精度，一般采用高階展開式進(jìn)行擬合，以升力系數(shù)為例可以表示為

將射流的影響分為兩部分：一部分設(shè)置為初始狀態(tài)量，包含在系數(shù)f0中，另一部分用于修正迎角帶來的氣動力影響，包含在系數(shù)f1和f2中。將系數(shù)fi設(shè)置為射流動量系數(shù)的三次多項(xiàng)式函數(shù)：

分別在不同的吹氣動量系數(shù)下驗(yàn)證Kriging模型與傳統(tǒng)氣動導(dǎo)數(shù)模型的準(zhǔn)確性，將2種模型的預(yù)測結(jié)果與CFD數(shù)值計(jì)算值進(jìn)行比較，圖7給出了不同動量系數(shù)下2種模型的預(yù)測情況。

圖7 Kriging模型與傳統(tǒng)氣動導(dǎo)數(shù)模型對比結(jié)果Fig.7 Comparison between Kriging model and aerodynamic derivative model

從圖7中可以看出，在采樣范圍的動量系數(shù)(Cμ為0～0.06）作用下，Kriging模型較傳統(tǒng)的氣動導(dǎo)數(shù)模型預(yù)測精度有明顯提升，預(yù)測結(jié)果更加接近CFD計(jì)算值，并能準(zhǔn)確預(yù)測失速迎角范圍；對于采樣范圍之外的參數(shù)，Kriging模型仍然能保持一定的預(yù)測精度，預(yù)測效果較傳統(tǒng)氣動導(dǎo)數(shù)模型有較大提升。

3 非定常氣動力建模

環(huán)量控制射流相較于傳統(tǒng)舵面而言可以大幅度提高翼型升力、改善升阻比，這使得環(huán)量控制飛行器往往在小迎角下就能夠獲得傳統(tǒng)飛行器在大迎角下的氣動力特性。因此，本文主要進(jìn)行中低迎角(－5°～10°）下的環(huán)量控制翼型非定常氣動力建模。

3.1 環(huán)量控制翼型非定常氣動力模型

在機(jī)翼小振幅強(qiáng)迫振蕩風(fēng)洞試驗(yàn)中，若還未出現(xiàn)深動態(tài)失速現(xiàn)象，則根據(jù)小擾動假設(shè)，氣動力參數(shù)可以由高階線性狀態(tài)方程表示，包括定常和非定常兩部分［28］。其中定常氣動力反映了機(jī)翼振蕩過程中的氣動力快速響應(yīng)部分，非定常氣動力部分則反映了機(jī)翼運(yùn)動導(dǎo)致的流場內(nèi)渦系移動產(chǎn)生的氣動力遲滯效應(yīng)。已有的研究表明［29］，環(huán)量控制射流通過改變翼型表面壓力分布、控制分離點(diǎn)移動而改變氣動力系數(shù)，這與迎角對氣動力系數(shù)的作用機(jī)理一致。對于環(huán)量控制翼型而言，動量系數(shù)與迎角共同影響翼型的氣動力系數(shù)，并且呈現(xiàn)高度耦合。以縱向運(yùn)動為例，氣動力參數(shù)可以表示為

式中：Ci為氣動力力矩系數(shù)；Catt(α，Cμ）為氣流在附著流動的狀態(tài)下，即無分離假設(shè)下的氣動力系數(shù)；Cα˙att(α，Cμ）為機(jī)翼定常俯仰下的氣動力導(dǎo)數(shù)，為迎角變化率；Cdyn為非定常氣動力增量。

上述模型適用于尚未出現(xiàn)深動態(tài)失速時(shí)的氣動力系數(shù)變化情況，其中，Cdyn采用以下線性微分方程式表示：

式中：τ為非定常運(yùn)動中的特征時(shí)間常數(shù)；Cst(α，Cμ）為實(shí)際定常狀態(tài)下氣動力系數(shù)，即Kriging模型預(yù)測值；ΔC(α，Cμ）為由于氣流分離導(dǎo)致的定常氣動力增量。

對于環(huán)量控制翼型而言，假設(shè)在迎角俯仰過程中射流動量系數(shù)不變，這種假設(shè)是合理的，目前對于傳統(tǒng)舵面飛行器在大迎角下非定常氣動力模型的建立也是基于舵面偏角不變的情況。

3.2 非定常氣動力模型參數(shù)辨識

借助經(jīng)典的強(qiáng)迫振動法計(jì)算諧振運(yùn)動過程中的非定常氣動力，強(qiáng)迫翼型繞重心做小幅度的俯仰振蕩，其迎角變化規(guī)律為

式中：Δα為迎角改變幅度。

非定常氣動力系數(shù)可以看做時(shí)間ˉt的連續(xù)函數(shù)，因此，可以將翼型小幅度俯仰振蕩過程中的氣動力系數(shù)變化按照減縮頻率進(jìn)行傅里葉分解，以升力系數(shù)為例，保留一階形式可以表示為

式(10）表示了氣動力系數(shù)可類似于“小擾動”假設(shè)，由基準(zhǔn)系數(shù)加上擾動系數(shù)構(gòu)成。其中，C0(α0，Cμ0）表示中心迎角對應(yīng)下的基準(zhǔn)氣動力系數(shù)；Cα(α0，Cμ0，k）稱為同相導(dǎo)數(shù)；Cα˙(α0，Cμ0，k）稱為異相導(dǎo)數(shù)。

基于小振幅情況下的線性假設(shè)：

該假設(shè)表示在小振幅情況下實(shí)際靜態(tài)氣動力系數(shù)與氣流不分離假設(shè)下靜態(tài)氣動力系數(shù)之量與迎角變化近似呈現(xiàn)線性關(guān)系。ΔCα(α0，Cμ0）表示該差值與迎角變化之間的導(dǎo)數(shù)。

借助歐拉方程求解微分方程式(7），即求解下述方程組：

使用微分算子法容易解得方程組(12）的解，將式(8）、式(9）的迎角變化規(guī)律代入式(6）中，與式(10）的相關(guān)系數(shù)作對比后整理得到以下關(guān)系：

由上述推導(dǎo)可知，同相導(dǎo)數(shù)并不等同于靜導(dǎo)數(shù)，異相導(dǎo)數(shù)也并不等同于動導(dǎo)數(shù)。

但是當(dāng)τ非常小時(shí)，此時(shí)同相導(dǎo)數(shù)和異相導(dǎo)數(shù)分別退化為傳統(tǒng)線性氣動導(dǎo)數(shù)模型的靜導(dǎo)數(shù)Cαst(α0，Cμ0）和動導(dǎo)數(shù)Cα˙att(α0，Cμ0），即如下關(guān)系式：

式中：Cα(α0，Cμ0）＝Cαst(α0，Cμ0）。

此時(shí)，Cdyn(0）＝0，C0(α0，Cμ0）＝Cst(α0，Cμ0）。

τ定量地表征了翼型非定常氣動力效應(yīng)的顯著程度。解式(13）可得

由于在同一中心迎角、不同減縮頻率對應(yīng)下的同相／異相導(dǎo)數(shù)數(shù)據(jù)點(diǎn)落在同一條直線上，而直線的斜率則表示特征時(shí)間常數(shù)的負(fù)值：－τ。

對于某一特定中心迎角下的小振幅強(qiáng)迫振動數(shù)據(jù)而言，采用最小二乘法可得到對應(yīng)中心迎角下的同相／異相導(dǎo)數(shù)數(shù)據(jù)點(diǎn)，部分?jǐn)?shù)據(jù)點(diǎn)如圖8所示。

圖8 不同迎角、減縮頻率下的同相／異相導(dǎo)數(shù)Fig.8 In-phase and out-phase derivatives at different angles of attack and reduced frequencies

由上述分析可知，對于特定的小振幅流動狀態(tài)(特定的迎角和動量系數(shù)）下，τ是固定的，計(jì)算得到的不同中心迎角、動量系數(shù)下的特征時(shí)間常數(shù)如圖9所示。

圖9 時(shí)間常數(shù)隨動量系數(shù)的變化曲線Fig.9 Variation curve of time constantwith momentum coefficient

特征時(shí)間常數(shù)值受迎角和動量系數(shù)影響較大。總體來看，τ隨Cμ增大呈現(xiàn)先降低后增加的趨勢，這說明低動量系數(shù)下，翼型的非定常氣動效應(yīng)減弱，但是當(dāng)動量系數(shù)增加到一定值后，翼型的非定常氣動力效應(yīng)增強(qiáng)，并且隨著迎角的增大，該臨界動量系數(shù)逐漸降低。當(dāng)α0＝0°時(shí)，在全動量系數(shù)范圍內(nèi)時(shí)間常數(shù)趨于0，氣動力模型可近似用傳統(tǒng)線性氣動導(dǎo)數(shù)模型代替。當(dāng)α0＝10°時(shí)，在全動量系數(shù)范圍內(nèi)時(shí)間常數(shù)都較大，翼型處于臨近動態(tài)失速或失速后狀態(tài)，此時(shí)翼型的非定常氣動力特性主要由后緣渦系的遲滯運(yùn)動導(dǎo)致，傳統(tǒng)線性氣動導(dǎo)數(shù)模型不再適用。

對于τ較小的情況，基于傳統(tǒng)線性氣動導(dǎo)數(shù)模型，使用不同迎角和動量系數(shù)下非定常運(yùn)動的若干離散點(diǎn)，借助最小二乘法辨識靜導(dǎo)數(shù)Cαst(α0，Cμ0）和動導(dǎo)數(shù)Cα˙att(α0，Cμ0）。

對于τ較大的情況，由于無法使用傳統(tǒng)的線性氣動導(dǎo)數(shù)模型，在獲得τ(α，Cμ）后，仍需要辨識同相／異相導(dǎo)數(shù)中的Cαatt(α0，Cμ0）、ΔCα(α0，Cμ0）和Cα˙att(α0，Cμ0）。

為了便于描述，將式(13）用以下方程描述：

式中：y1i、y2i、b0、a0、c0分別為式(13）中的相關(guān)參數(shù)。為了辨識相關(guān)參數(shù)，引入以下誤差函數(shù)：

式中：n為同一中心迎角和動量系數(shù)、不同頻率下的仿真實(shí)驗(yàn)次數(shù)；y1i和y2i分別為特定小振幅流動狀態(tài)下的同相／異相導(dǎo)數(shù)；b0、a0和c0為參數(shù)辨識后的相應(yīng)值。該方法的核心思想是兩步線性回歸法，即先根據(jù)非定常運(yùn)動的離散點(diǎn)辨識得到y(tǒng)1i和y2i，再根據(jù)得到的y1i和y2i辨識b0、a0和c0。辨識過程中，因?yàn)榇嬖谥行挠?、動量系?shù)、減縮頻率3個(gè)自變量，直接進(jìn)行參數(shù)辨識需要耗費(fèi)極大的仿真計(jì)算量，所以本節(jié)針對射流動量系數(shù)與迎角的相互作用特點(diǎn)，對b0、a0和c0做合理優(yōu)化，提高計(jì)算效率。

對于b0而言，Jones［25］曾指出動量系數(shù)對氣動力的影響呈現(xiàn)2個(gè)比較明顯的線性階段，其流動狀態(tài)滿足無分離假設(shè)，因此在低動量系數(shù)下b0可以用以下關(guān)系式代替：

式中：CαCμst(α0）為Cαst(α0）對Cμ的導(dǎo)數(shù)，可直接由Kriging模型插值得到，待辨識參量變?yōu)橛?。采用這樣的簡化方法在低迎角線性階段可直接退化為Cαst(α0，Cμ0），與無分離情況下的氣動導(dǎo)數(shù)一致。

對于c0而言，當(dāng)自由來流不發(fā)生變化時(shí)，迎角俯仰導(dǎo)致翼型表面產(chǎn)生繞重心的旋轉(zhuǎn)速度，與外流相對流動產(chǎn)生下洗效應(yīng)，而后緣射流并不會直接隨著翼型俯仰，而是由于后緣壓力場的變化時(shí)刻處于動態(tài)平衡之中。由定常旋轉(zhuǎn)產(chǎn)生的準(zhǔn)定常氣動力增量主要在小迎角時(shí)起作用，與迎角變化率近似呈線性關(guān)系。在較大迎角范圍內(nèi)，其引起的準(zhǔn)定常氣動力增量遠(yuǎn)小于非定常氣動力增量。將動量系數(shù)作為修正參數(shù)，可以將c0表示為

對于a0而言，a0＝Cst(α，Cμ）－Catt(α，Cμ），因此只需辨識b0即可。在辨識得到τ和ΔC(α，Cμ）后，根據(jù)非定常氣動力增值周期性變化的特性，即Cdyn(0）＝Cdyn(2π／ω），即可求得初值Cdyn(0）。

以誤差函數(shù)E最小為優(yōu)化目標(biāo)函數(shù)，以升力系數(shù)為例，得到的Cαatt(α0，Cμ0）曲線如圖10所示。

圖10 Cαatt(α0，Cμ0）的辨識結(jié)果Fig.10 Identification results of Cαatt(α0，Cμ0）

3.3 非定常氣動力模型與仿真數(shù)據(jù)對比驗(yàn)證

為了驗(yàn)證上述非定常氣動力模型的準(zhǔn)確性，分別對環(huán)量控制翼型在不同迎角、動量系數(shù)、振幅下的預(yù)測結(jié)果進(jìn)行對比分析。預(yù)測結(jié)果取升力系數(shù)CL和俯仰力矩系數(shù)Cm進(jìn)行對比分析。計(jì)算條件為：來流速度為30 m／s，來流溫度為293.15 K。中心迎角、動量系數(shù)、振幅、減縮頻率驗(yàn)證如圖11、圖12所示。

如圖11所示，為了得到非定常氣動力模型對大振幅情況下的預(yù)測結(jié)果，將大振幅的流動情況分段處理，每一段的氣動力變化情況都可以用小振幅情況下迎角對應(yīng)的同相／異相導(dǎo)數(shù)表示，分段處理后可以將大振幅情況下的預(yù)測結(jié)果近似表示為

圖11 小振幅流動狀態(tài)下結(jié)果對比Fig.11 Comparison of results under small-amplitude flow

將分段處理的αm取為1°，得到非定常氣動力模型對于大振幅情況下的預(yù)測結(jié)果，如圖12所示。

圖12 大振幅流動狀態(tài)下結(jié)果對比Fig.12 Comparison of results under large-amp litude flow

可以看出，對于升力系數(shù)而言，線性非定常氣動力模型對小振幅情況下流動狀態(tài)有較好的預(yù)測精度，特別是在0°迎角下，此時(shí)升力系數(shù)線性特征明顯，在較大振幅情況下(見圖11(a）、圖11(b））也能夠取得較好的預(yù)測精度。但是在較大迎角、較高動量系數(shù)的流動狀態(tài)下，如圖12(c）所示，此時(shí)升力系數(shù)非線性特征明顯，線性非定常氣動力模型難以取得較高的預(yù)測精度。對于俯仰力矩系數(shù)而言，線性非定常氣動力模型在各個(gè)動量系數(shù)下都能較好地預(yù)測力矩系數(shù)變化趨勢。這是因?yàn)楦┭隽叵禂?shù)主要受翼型上表面分離渦影響，當(dāng)翼型上表面分離渦移動時(shí)，對應(yīng)的負(fù)壓區(qū)逐漸擴(kuò)大或者縮小，導(dǎo)致俯仰力矩系數(shù)變化。高動量系數(shù)下翼型上表面不再分離，導(dǎo)致遲滯效應(yīng)逐漸減弱。值得一提的是，在出現(xiàn)動態(tài)失速分離前，高動量系數(shù)具有減弱氣動力遲滯效應(yīng)的優(yōu)勢，如圖11(f）、圖12(f）所示，此時(shí)俯仰力矩的變化幾乎不存在遲滯效應(yīng)，與定常狀態(tài)下的變化幾乎一致。在出現(xiàn)動態(tài)失速后，同樣能在線性小迎角階段減弱遲滯效應(yīng)，如圖12(b）、圖12(c）所示。

3.4 非定常氣動力模型非線性修正

當(dāng)振蕩幅度增大并出現(xiàn)深動態(tài)失速的情況時(shí)，上述線性非定常氣動力模型將會在失速段失去預(yù)測精度。這是因?yàn)樵诖笳穹闆r下進(jìn)行分段處理完全忽略了洗流時(shí)差的累計(jì)效應(yīng)，在非線性較為明顯的階段將完全喪失其預(yù)測精度。為了改善非定常氣動力模型在動態(tài)失速情況下的預(yù)測精度，采用文獻(xiàn)［30］的方法將非定常氣動力增量式(7）改寫為以下非線性形式：

式中：ki(α，Cμ）＝1／τi(α，Cμ），該方法采用多次逼近的思想，避免了大振幅過程中僅用單一時(shí)間常數(shù)τ線性表示增量的誤差，同時(shí)充分考慮了累計(jì)效應(yīng)。利用上述線性非定常氣動力模型辨識得到的Catt(α，Cμ）和Cα˙att(α，Cμ），可以獲得線性假設(shè)下CL的變化，與實(shí)際CL的差值即非定常氣動力增量，與線性假設(shè)下得到的Cdyn相疊加即實(shí)際情況下的Cdyn，保持ΔC(α，Cμ）和k1(α，Cμ）不變，利用最小二乘法即可辨識得到ki(α，Cμ），采用這樣的假設(shè)能夠有效降低參數(shù)辨識的難度。

上述線性非定常氣動力模型已經(jīng)能夠滿足無射流情況下中低迎角的大振幅預(yù)測精度，因此，非線性項(xiàng)主要用來提高高動量系數(shù)情況下的預(yù)測精度。為了提高辨識效率，本文僅對包含Kriging模型中靜態(tài)失速點(diǎn)附近的流動狀態(tài)進(jìn)行參數(shù)辨識。采用三次非線性形式描述非定常氣動力增量，辨識得到的k2(α，Cμ）和k3(α，Cμ）如圖13所示。

圖13 ki(α，Cμ）辨識結(jié)果Fig.13 Identification result of ki(α，Cμ）

非線性項(xiàng)修正后的非定常氣動力模型進(jìn)行大振幅流動情況下的升力系數(shù)預(yù)測，結(jié)果如圖14所示，其結(jié)果較線性模型有明顯提升。

圖14 非線性非定常氣動力模型結(jié)果對比(Cμ ＝0.05）Fig.14 Comparison of results of nonlinear unsteady aerodynamic models(Cμ ＝0.05）

4 結(jié)論

本文主要進(jìn)行了環(huán)量控制翼型非定常氣動力的建模，通過上述分析，可以得到以下結(jié)論：

1）基于Kriging模型完成了環(huán)量控制翼型的定常氣動力插值，能夠精確的反映其定常氣動力變化，插值精度要優(yōu)于傳統(tǒng)氣動導(dǎo)數(shù)模型。

2）基于微分方程模型，通過添加射流參數(shù)源項(xiàng)，完成了適用于環(huán)量控制翼型的線性微分方程建模。對Cαatt(α0，Cμ0）和Cα˙att(α0，Cμ0）等參數(shù)簡化處理，提高辨識效率。

3）對高動量系數(shù)情況下的非定常氣動力模型進(jìn)行非線性修正，能夠較為準(zhǔn)確的反映不同動量系數(shù)下的非定常流動特性。

4）對于升力系數(shù)而言，線性非定常氣動力模型能夠滿足小振幅和低動量系數(shù)下大振幅流動狀態(tài)；對于俯仰力矩系數(shù)而言，線性非定常氣動力模型能夠滿足大部分流動狀態(tài)。

環(huán)量控制翼型的非定常氣動力耦合特性復(fù)雜，仍需要在以下方面進(jìn)一步研究：①高動量系數(shù)下缺乏足夠精度的數(shù)學(xué)模型來描述后緣渦的分離變化，需要對流場特性進(jìn)行更深入的研究，進(jìn)而優(yōu)化模型結(jié)構(gòu)。提高在高動量系數(shù)下非線性非定常氣動力模型預(yù)測效果。②關(guān)于環(huán)量控制翼型的非定常氣動力建模研究仍處于起步階段，如何有效的結(jié)合已經(jīng)較為成熟的大迎角非定常氣動力數(shù)學(xué)模型是未來研究的途徑之一。