基于DMD方法的翼型大迎角失速流動(dòng)穩(wěn)定性研究

葉 坤, 葉正寅, 武 潔, 屈 展

(西北工業(yè)大學(xué) 翼型葉柵空氣動(dòng)力學(xué)國(guó)家級(jí)重點(diǎn)實(shí)驗(yàn)室, 陜西 西安 710072)

0 引 言

失速是航空航天領(lǐng)域中一種普遍存在的流動(dòng)現(xiàn)象，失速研究具有重要的學(xué)術(shù)和應(yīng)用價(jià)值。傳統(tǒng)認(rèn)為靜態(tài)失速是當(dāng)迎角增大到某一定值時(shí)，在翼型上表面會(huì)出現(xiàn)大規(guī)模分離流動(dòng)，從而導(dǎo)致升力突然下降和阻力突然增加及非定常壓力脈動(dòng)放大等負(fù)面效應(yīng)，這對(duì)飛行器的整體氣動(dòng)特性和操縱性能產(chǎn)生嚴(yán)重影響[1]。隨著航空航天新技術(shù)的不斷發(fā)展，對(duì)現(xiàn)代飛行器的性能要求越來(lái)越高，如大迎角機(jī)動(dòng)性，蘇27過(guò)失速“眼鏡蛇”機(jī)動(dòng)時(shí)迎角超過(guò)了120°，用于改善飛機(jī)起降性能的增升襟翼也經(jīng)常在大迎角條件下工作，特別是利用迎角40°左右存在的二次升力峰值和阻力是著陸時(shí)增升增阻的一種理想選擇[2]。因此，只有深入理解大迎角復(fù)雜分離流動(dòng)演化規(guī)律，才能更好地控制失速現(xiàn)象，使飛行器可用迎角逾越靜態(tài)失速這一傳統(tǒng)禁區(qū)。

深入理解大迎角復(fù)雜分離流動(dòng)演化規(guī)律，首先需要獲取較準(zhǔn)確的流場(chǎng)。近年來(lái)，LES方法和DES方法被廣泛地應(yīng)用于分離流動(dòng)的模擬中[3-7]。在獲取較準(zhǔn)確的流場(chǎng)信息后，通常采用非線(xiàn)性理論、模態(tài)分解等方法對(duì)流動(dòng)結(jié)構(gòu)進(jìn)行分析。從非線(xiàn)性動(dòng)力學(xué)角度講，靜態(tài)失速中存在滯后、突跳及分岔現(xiàn)象。Mittal等[8]通過(guò)數(shù)值模擬對(duì)翼型靜態(tài)失速迎角附近的遲滯現(xiàn)象進(jìn)行了研究，并發(fā)現(xiàn)遲滯現(xiàn)象的發(fā)生是由于流動(dòng)具有對(duì)歷史行為的記憶能力。張家忠等[9]對(duì)翼型失速的非線(xiàn)性動(dòng)力學(xué)特性進(jìn)行研究，發(fā)現(xiàn)靜態(tài)失速是一種“鞍—結(jié)”分岔現(xiàn)象。另一方面，在獲得流動(dòng)結(jié)構(gòu)的眾多方法中，本征正交分解方法(Proper Orthogonal Decomposition， POD)被廣泛使用。自L(fǎng)umley[10]首先將POD方法引入湍流研究以來(lái)， Sirovich[11]引入了snapshot 方法來(lái)研究波動(dòng)流的動(dòng)力學(xué)問(wèn)題?？祩サ萚12]基于POD方法，采用非線(xiàn)性動(dòng)力學(xué)理論對(duì)翼型繞流的多模態(tài)耦合機(jī)制進(jìn)行分析，研究了模態(tài)耦合作用與流動(dòng)穩(wěn)定性的關(guān)系。

翼型大迎角分離渦的演化包含了渦量生成、積聚和耗散，其中的流動(dòng)機(jī)理非常復(fù)雜。本文采用DMD方法對(duì)分離渦的演化過(guò)程進(jìn)行穩(wěn)定性分析。首先基于SST-DES方法對(duì)NACA0012進(jìn)行數(shù)值模擬得到非定常流場(chǎng)數(shù)據(jù)。實(shí)際翼型大迎角分離流動(dòng)會(huì)存在明顯的三維效應(yīng)，且通常基于DES方法進(jìn)行計(jì)算時(shí)，對(duì)物面展向長(zhǎng)度和網(wǎng)格數(shù)量有一定的要求，文獻(xiàn)[22]中對(duì)此進(jìn)行了較詳細(xì)的分析。但受計(jì)算資源的限制，本文數(shù)值模擬中網(wǎng)格展向僅拉伸一層。事實(shí)上，對(duì)于一些經(jīng)典的三維分離流動(dòng)，例如圓柱、凹腔等流動(dòng)，大部分學(xué)者及文獻(xiàn)也均是從二維的角度進(jìn)行計(jì)算和分析。

1 計(jì)算方法及驗(yàn)證算例

三維非定常N-S方程在直角坐標(biāo)系中的積分守恒形式為:

(1)

其中Ω為控制體，?Ω為控制體單元邊界，Q為守恒變量，Q=[ρ,ρu·ρv,ρw,ρE]T，F(xiàn)(Q)為無(wú)黏通量，F(xiàn)V(Q)為黏性通量。其中ρ、P、T、E分別為密度、壓強(qiáng)、溫度和單位質(zhì)量流體的總能，u、v、w分別為x、y、z方向的速度分量。

采用基于剪切應(yīng)力輸運(yùn)SSTk-w兩方程湍流模型的SST-DDES混合方法對(duì)流場(chǎng)進(jìn)行模擬，采用有限體積方法進(jìn)行空間離散，空間格式采用Roe，時(shí)間推進(jìn)采用LU-SGS。

目前公開(kāi)的翼型大迎角實(shí)驗(yàn)結(jié)果較少，本文采用文獻(xiàn)[23]中NACA0012的實(shí)驗(yàn)數(shù)據(jù)。實(shí)驗(yàn)條件為Ma=0.15、Re=6×106。通常采用DES方法進(jìn)行計(jì)算時(shí)，對(duì)物面展向長(zhǎng)度和網(wǎng)格數(shù)量有一定的要求，文獻(xiàn)[22]中對(duì)此進(jìn)行了較詳細(xì)的分析。但受計(jì)算資源的限制，數(shù)值模擬中網(wǎng)格展向僅拉伸1層。圖1為生成的結(jié)構(gòu)網(wǎng)格，翼型弦長(zhǎng)為1 m。為了避免使用拉伸比很大的“C”型網(wǎng)格，采用420×180的“O”型結(jié)構(gòu)網(wǎng)格。其中翼型上表面300個(gè)點(diǎn)，下表面120個(gè)點(diǎn)，物面網(wǎng)格第一層高度為5×10-6m，滿(mǎn)足y+<1，附面層徑向采用40層網(wǎng)格，增長(zhǎng)率為1.25。為準(zhǔn)確模擬分離，對(duì)上翼面附近區(qū)域的網(wǎng)格進(jìn)行加密，網(wǎng)格層數(shù)為100層，增長(zhǎng)率為1，之后采用40層過(guò)渡到遠(yuǎn)場(chǎng)，遠(yuǎn)場(chǎng)距離物面約為100倍弦長(zhǎng)。非定常計(jì)算中的時(shí)間步長(zhǎng)為5×10-4s。

圖2為升力系數(shù)(時(shí)間平均下的升力系數(shù))隨迎角的變化曲線(xiàn)。實(shí)驗(yàn)數(shù)據(jù)中，迎角為17.13°時(shí)，升力系數(shù)達(dá)到最大，當(dāng)迎角為18.21°時(shí)，升力系數(shù)突然大幅下降。本文計(jì)算結(jié)果中，當(dāng)迎角為17°時(shí)，升力系數(shù)達(dá)到最大。為進(jìn)一步細(xì)致了解17°至18°之間的流動(dòng)，在迎角17.5°至18°之間，每隔0.1°進(jìn)行一次非定常計(jì)算。當(dāng)迎角由17.7°進(jìn)一步增加到17.8°時(shí)，升力系數(shù)突然大幅下降，且開(kāi)始呈現(xiàn)周期性變化。圖3為17.7°和17.8°的升力系數(shù)隨時(shí)間的變化，圖4為流場(chǎng)流線(xiàn)和渦量云圖。上述結(jié)果表明本文的計(jì)算結(jié)果與實(shí)驗(yàn)數(shù)據(jù)吻合較好，說(shuō)明本文計(jì)算結(jié)果是可信的。

2 DMD方法

DMD方法只基于流動(dòng)快照而不受模型限制，通過(guò)提取流動(dòng)中的模態(tài)，準(zhǔn)確地描述流動(dòng)結(jié)構(gòu)。對(duì)于線(xiàn)化流動(dòng)而言，被提取的模態(tài)可以用來(lái)表征全局穩(wěn)定模態(tài)。對(duì)于非線(xiàn)性流動(dòng)而言，DMD的模態(tài)描述了在數(shù)據(jù)序列中起支配作用的流動(dòng)結(jié)構(gòu)。以下為其數(shù)學(xué)推導(dǎo)過(guò)程。

(2)

第一步，假設(shè)線(xiàn)性映射A是快照vi到vi+1的映射，即:

vi+1=Avi

(3)

且在全部區(qū)域和整個(gè)時(shí)間段內(nèi)的采樣都滿(mǎn)足上述映射關(guān)系。對(duì)于緩慢變化的系統(tǒng)，乘法準(zhǔn)則是上述假設(shè)的基礎(chǔ)。

={Av1,Av2,Av3,…,AvN-1}

(4)

A≈UFU*

(5)

(6)

其中，Σ是一個(gè)r×r的對(duì)角矩陣，有非0的對(duì)角元素{σ1,…,σr}，并且，U和V為酉矩陣。

U∈CM×r，并且U*U=I

V∈Cr×N，并且V*V=I

將式(5)和式(6)代入式(4)中，可得：

(7)

F即為A的最優(yōu)低維估計(jì)矩陣，由此通過(guò)求解F的特征值和特征向量即可得到DMD分析結(jié)果。

假設(shè)F是一個(gè)r維子空間內(nèi)的精確映射矩陣而不是對(duì)A的估計(jì)，那么

xk+1=Fxk

(8)

則A的POD模態(tài)U可以看作從xk到vk的近似映射，則

vk≈Ux

(9)

如果F的特征向量個(gè)數(shù)與秩相等且都線(xiàn)性獨(dú)立，那么F可以寫(xiě)做如下形式

(10)

(11)

(12)

(13)

由此可得：

[v1,v2,…,vN-1]≈[φ1,φ2,…，φr]·

(14)

動(dòng)態(tài)模態(tài)分解獲得的最重要的結(jié)果是模態(tài)φi和模態(tài)對(duì)應(yīng)的特征值μi,lg(μi/Δt)的實(shí)部為對(duì)應(yīng)模態(tài)的

放大率，其虛部表示對(duì)應(yīng)模態(tài)的頻率。在進(jìn)行穩(wěn)定性判斷時(shí)，放大率為正，則對(duì)應(yīng)的模態(tài)發(fā)散；放大率為負(fù)，則對(duì)應(yīng)的模態(tài)收斂；若放大率為0，則對(duì)應(yīng)的模態(tài)為穩(wěn)定極限環(huán)模態(tài)。同樣也可以將μi置于復(fù)平面的單位圓上進(jìn)行判斷，若在單位圓內(nèi)，則收斂；單位圓外，則發(fā)散；單位圓上，則為穩(wěn)定周期性模態(tài)。

3 計(jì)算結(jié)果分析

本文在驗(yàn)證算例中實(shí)驗(yàn)條件的馬赫數(shù)和雷諾數(shù)下對(duì)NACA0012大迎角失速翼型流場(chǎng)的渦穩(wěn)定性演化規(guī)律進(jìn)行研究。首先，對(duì)迎角在17.8°至40°范圍內(nèi)的狀態(tài)進(jìn)行非定常計(jì)算，其中22°至24°的計(jì)算中，迎角間隔為1°。圖5為Strouhal數(shù)(St)隨迎角的變化曲線(xiàn)，

St=fC/U∞

(15)

式中f為頻率，C為弦長(zhǎng)，U∞為來(lái)流速度。

可以看出，St隨迎角的增加并不是呈現(xiàn)單調(diào)減小的規(guī)律。當(dāng)迎角在17.8°至23°范圍內(nèi)，隨迎角的增加，St基本呈線(xiàn)性減小；當(dāng)迎角由23°增加到24°時(shí)，St大幅減小；當(dāng)迎角在24°至29°范圍內(nèi)，St隨迎角增加緩慢減小，并且迎角為29°時(shí)，St達(dá)到最小，為0.166168；迎角由29°增加到30°時(shí)，St有微幅的增加；當(dāng)迎角進(jìn)一步由30°增加到31°時(shí)，St大幅增加；迎角從31°增加到40°時(shí)，St再次基本呈線(xiàn)性減小。因此，采用DMD方法分別對(duì)以下3種迎角狀態(tài)下的分離流場(chǎng)進(jìn)行穩(wěn)定性分析：(1) 迎角為17.8°時(shí)，流動(dòng)剛剛進(jìn)入失速狀態(tài)，稱(chēng)為失速初始狀態(tài)；(2) 迎角增加到24°，流動(dòng)進(jìn)一步分離，St突然大幅減小，稱(chēng)為淺失速狀態(tài)；(3) 迎角進(jìn)一步增加到29°時(shí)，St達(dá)到最小，稱(chēng)為深失速狀態(tài)。

由于在上述三種迎角下，流場(chǎng)最后均表現(xiàn)為穩(wěn)定極限環(huán)狀態(tài)，升力系數(shù)均呈現(xiàn)周期性變化。因此，采用DMD方法進(jìn)行穩(wěn)定性分析時(shí)，DMD采樣中僅提取一個(gè)周期的瞬態(tài)流場(chǎng)數(shù)據(jù)，三種迎角下的樣本數(shù)量分別為65、182、234。

3.1 DMD穩(wěn)定性分析

圖6為升力系數(shù)隨時(shí)間的變化。三種狀態(tài)下，升力系數(shù)在一個(gè)周期內(nèi)的變化相差較大。失速初始狀態(tài)下，升力系數(shù)一個(gè)周期內(nèi)呈現(xiàn)類(lèi)似正弦曲線(xiàn)的變化，周期性變化規(guī)律簡(jiǎn)單，一個(gè)周期中僅有一個(gè)極大值和一個(gè)極小值，最大值與最小值相對(duì)均值的振蕩幅值分別為28.6%和32.2%。淺失速狀態(tài)中，升力系數(shù)一個(gè)周期內(nèi)的變化較復(fù)雜，一個(gè)周期中存在四個(gè)極大值和四個(gè)極小值，在5/8T時(shí)刻附近達(dá)到的極大值略小于1/8T時(shí)刻的最大值，最大值與最小值相對(duì)均值的振蕩幅值分別為24.4%和18.6%。深失速狀態(tài)中，升力系數(shù)在一個(gè)周期中存在兩個(gè)極大值和兩個(gè)極小值，且在5/8T時(shí)刻附近達(dá)到的極大值小于1/8T時(shí)刻的最大值。最大值與最小值相對(duì)均值的振蕩幅值分別為38.3%和17.6%。

圖7為升力系數(shù)FFT分析結(jié)果?？梢?jiàn)FFT可提取其主頻及高階倍頻，圖中僅標(biāo)出前五階頻率，按照St對(duì)應(yīng)的峰值對(duì)St進(jìn)行排序，三種狀態(tài)的順序各不相同。失速初始狀態(tài)，St的順序?yàn)椋篠t1，St2，St3，St4，St5；淺失速狀態(tài)，St的順序?yàn)椋篠t2，St4，St1，St3，St5；深失速狀態(tài)，St的順序?yàn)椋篠t2，St1，St4，St3，St5。

表1為DMD方法提取頻率與FFT分析結(jié)果的對(duì)比，兩者之間頻率最大差異為0.15037%，說(shuō)明DMD方法可以非常準(zhǔn)確地提取流動(dòng)中的頻率。

(a)α=17.8° (b)α=24.0° (c)α=29.0°

圖6升力系數(shù)隨時(shí)間的變化
Fig.6Liftcoefficientvaryingwithtimeforthethreestates

(a)α=17.8° (b)α=24.0° (c)α=29.0°

圖7 升力系數(shù)FFT分析結(jié)果Fig.7 FFT results of lift coefficient

基于DMD方法分別對(duì)三種狀態(tài)下的流場(chǎng)進(jìn)行穩(wěn)定性分析。圖8為翼型三種迎角狀態(tài)下DMD特征值的實(shí)部和虛部在單位圓上的分布，圖中黑色點(diǎn)表示在單位圓上或單位圓內(nèi)，紅色點(diǎn)表示在單位圓外。圖9為特征值取對(duì)數(shù)之后的實(shí)部和虛部的分布，圖中黑色點(diǎn)表示實(shí)部等于0或小于0，紅色的點(diǎn)表示實(shí)部大于0。圖10為模態(tài)幅值與放大率之間的關(guān)系，圖中紅色的表示放大率正，黑色表示放大率為負(fù)或0。由前面DMD方法介紹中，這三個(gè)圖中黑色點(diǎn)對(duì)應(yīng)的模態(tài)為收斂模態(tài)或穩(wěn)定極限環(huán)模態(tài)，而紅色點(diǎn)對(duì)應(yīng)的模態(tài)為發(fā)散模態(tài)。三種迎角狀態(tài)下，不穩(wěn)定模態(tài)占總模態(tài)數(shù)量的比例分別為29.23%、29.28%和62.52%。

從圖8看，所有模態(tài)的特征值基本都落在單位圓上，且從圖9和圖10中可以看出，所有發(fā)散模態(tài)對(duì)應(yīng)的放大率都非常小，最大放大率分別為3.53×10-5，3.56×10-3和6.95×10-4，這些不穩(wěn)定模態(tài)處于弱發(fā)散狀態(tài)，其他模態(tài)處于收斂或極限環(huán)狀態(tài)。

圖11為模態(tài)系數(shù)幅值隨St的變化，可以看出，模態(tài)系數(shù)幅值隨St的增加，基本呈現(xiàn)逐漸減小的趨勢(shì)，說(shuō)明St較低的模態(tài)在流場(chǎng)中占據(jù)主要影響。同時(shí)，DMD提取的St范圍為0～20，覆蓋了脫落渦主頻St的0～33倍，說(shuō)明DMD方法提取到流動(dòng)頻率范圍很寬。

(a)α=17.8° (b)α=24.0° (c)α=29.0°

圖8特征值的實(shí)部與虛部在單位圓上的分布
Fig.8Eigenvaluesonunitcircleforthethreestates

(a)α=17.8° (b)α=24.0° (c)α=29.0°

圖9特征值取對(duì)數(shù)后的實(shí)部與虛部的分布
Fig.9Logarithmiceigenvaluesonunitcircleforthethreestates

(a)α=17.8° (b)α=24.0°(c)α=29.0°

圖10放大率與幅值之間的關(guān)系
Fig.10Amplitudevaryingwithgrowthrateforthethreestates

(a)α=17.8° (b)α=24.0° (c)α=29.0°

圖11幅值隨St的變化
Fig.11AmplitudevaryingwithStforthethreestates

3.2 分離渦演化及DMD模態(tài)分析

3.2.1 失速初始狀態(tài)

圖12為失速初始狀態(tài)下一個(gè)周期中8個(gè)瞬態(tài)流場(chǎng)的渦量圖以及流線(xiàn)圖，其中1/8T和5/8T時(shí)刻分別對(duì)應(yīng)升力系數(shù)最大和最小狀態(tài)。在1/8T時(shí)刻，翼型上表面存在一個(gè)順時(shí)針的大分離渦(主渦)，同時(shí)后緣存在一個(gè)較小的逆時(shí)針尾渦，并且翼型中段物面附近處存在一個(gè)由于剪切層分離形成的分離泡。當(dāng)進(jìn)入2/8T時(shí)刻時(shí)，前緣處生成一個(gè)前緣渦，同時(shí)翼型中段處的一個(gè)分離泡演化為兩個(gè)分離泡。當(dāng)進(jìn)入3/8T時(shí)刻時(shí)，2/8T時(shí)刻中靠近后緣處的的分離泡消失，同時(shí)靠近前緣處的分離泡和前緣渦進(jìn)一步發(fā)展，尾渦逐漸變小。當(dāng)進(jìn)入4/8T時(shí)刻時(shí)，在后緣處形成一個(gè)新的小分離渦，同時(shí)翼型中段處的分離泡和前緣渦進(jìn)一步發(fā)展，尾渦進(jìn)一步變小。當(dāng)進(jìn)入5/8T時(shí)刻時(shí)，尾渦脫落消失，同時(shí)主渦，前緣渦、翼型中段處的分離泡和后緣處分離渦進(jìn)一步發(fā)展。當(dāng)進(jìn)入6/8T時(shí)刻時(shí)，在后緣處分離渦發(fā)展成為新的尾渦，同時(shí)前緣渦、主渦和翼型中段處的分離泡均變小，且主渦呈現(xiàn)脫落的趨勢(shì)。當(dāng)進(jìn)入7/8T時(shí)刻時(shí)，主渦脫落消失，尾渦變大，前緣渦進(jìn)一步向后發(fā)展，翼型中段處的分離泡消失。當(dāng)進(jìn)入T時(shí)刻時(shí)，前緣渦進(jìn)一步發(fā)展成為新的主渦，尾渦也進(jìn)一步變大，流動(dòng)開(kāi)始向下一個(gè)周期發(fā)展?？傮w來(lái)看，分離渦在一個(gè)周期的演化過(guò)程中，前緣渦、主渦、翼型中段處的分離泡以及尾渦都經(jīng)歷了一次形成、發(fā)展、消失、再形成。

圖13為前六階DMD模態(tài)的速度云圖，其中Mode1對(duì)應(yīng)的St為0，Mode2～Mode6分別對(duì)應(yīng)St1～St5(淺失速狀態(tài)和深失速狀態(tài)與此相同)?？梢钥闯?，Mode1主要表現(xiàn)為整個(gè)翼型上表面的大分離渦結(jié)構(gòu)，其模態(tài)系數(shù)幅值不隨時(shí)間變化，說(shuō)明此模態(tài)的分離渦結(jié)構(gòu)為翼型分離渦演化過(guò)程中最主要的靜止模態(tài)。Mode2主要表現(xiàn)為前緣渦、翼型中段處的分離泡,以及后緣處的尾渦結(jié)構(gòu)，且其對(duì)應(yīng)的頻率為升力系數(shù)的主頻，結(jié)合一個(gè)周期中流場(chǎng)渦結(jié)構(gòu)的演化過(guò)程可以發(fā)現(xiàn)，此模態(tài)對(duì)應(yīng)的渦結(jié)構(gòu)與流動(dòng)中以此頻率進(jìn)行演化的渦結(jié)構(gòu)比較一致。Mode3～Mode6主要表現(xiàn)為后緣附近小尺度的尾渦以及尾跡區(qū)的渦結(jié)構(gòu)，其對(duì)應(yīng)的頻率為高階倍頻。圖14 為模態(tài)系數(shù)隨時(shí)間的變化，可以看出不同模態(tài)系數(shù)之間不僅頻率不同，并且不同系數(shù)之間存在明顯的相位差，說(shuō)明分離渦流動(dòng)中不同頻率對(duì)應(yīng)的渦結(jié)構(gòu)運(yùn)動(dòng)不同步。按照模態(tài)系數(shù)最大幅值對(duì)Mode2～Mode6進(jìn)行排序，發(fā)現(xiàn)排序后的模態(tài)對(duì)應(yīng)的St為St1、St2、St3、St4、St5，這與圖7(a)中頻率按峰值排序一致。

(a) Frst order mode (b) Second order mode (c) Third order mode

(d) Fourth order mode (e) Fifth order mode (f) Sixth order mode

圖13DMD模態(tài)速度云圖
Fig.13ContoursofvelocityforDMDmodes

3.2.2 淺失速狀態(tài)

圖15為淺失速狀態(tài)下一個(gè)周期中8個(gè)瞬態(tài)流場(chǎng)的渦量圖以及流線(xiàn)圖?？梢钥闯觯?/8T與5/8T、2/8T與6/8T、3/8T與7/8T、4/8T與T，每一對(duì)時(shí)

刻的流場(chǎng)均非常相似。從圖6(b)中可以看出，1/8T和5/8T時(shí)刻對(duì)應(yīng)升力系數(shù)的兩個(gè)較大的極值點(diǎn)；2/8T、4/8T、6/8T以及T時(shí)刻對(duì)應(yīng)升力系數(shù)的值較接近，但其中2/8T與6/8T時(shí)刻升力系數(shù)處于下降階段，4/8T與T時(shí)刻升力系數(shù)處于上升階段，3/8T和7/8T時(shí)刻的升力系數(shù)接近且均處于升力系數(shù)下降階段。

流場(chǎng)中渦的演化過(guò)程為：在1/8T時(shí)刻，翼型上表面存在一個(gè)順時(shí)針的大分離渦(主渦)，其后存在一個(gè)逆時(shí)針的二次渦，同時(shí)后緣處存在一個(gè)較小的逆時(shí)針?lè)蛛x渦；當(dāng)進(jìn)入2/8T時(shí)刻時(shí)，主渦和二次渦脫落消失，后緣處的小分離渦進(jìn)一步發(fā)展；3/8T時(shí)刻，主渦再次形成，后緣處的小分離渦發(fā)展為較大的二次渦；當(dāng)進(jìn)入4/8T時(shí)刻時(shí)，主渦和二次渦向下游發(fā)展；5/8T時(shí)刻時(shí)，二次渦進(jìn)一步向下游發(fā)展，后緣處再次形成一個(gè)小分離；之后分離渦演化過(guò)程與2/8T至4/8T的過(guò)程類(lèi)似。

圖16為前六階DMD模態(tài)的速度云圖。Mode1主要表現(xiàn)為整個(gè)翼型上表面的大分離渦結(jié)構(gòu)，這與失速初始狀態(tài)類(lèi)似，該分離渦結(jié)構(gòu)為分離渦演化過(guò)程中最主要的靜止模態(tài)。與失速初始狀態(tài)不同的是，Mode2、Mode3和Mode4主要表現(xiàn)為翼型后緣分離渦和主渦上方區(qū)域的一些渦結(jié)構(gòu)，Mode5和Mode6主要表現(xiàn)為尾跡區(qū)的渦結(jié)構(gòu)。實(shí)際上，隨著迎角的增加，翼型上表面逐漸成為一個(gè)背風(fēng)區(qū)，在空間上，流場(chǎng)中渦不再僅僅從靠近后緣區(qū)域脫落，整個(gè)翼型上表面上方區(qū)域會(huì)有渦脫落。例如圓柱和豎直平板的流動(dòng)，可看成迎角為90°的翼型。文獻(xiàn)[4]和文獻(xiàn)[9]中分別采用LES和DES方法，在對(duì)翼型分離流動(dòng)的模擬中，均捕捉到了這種分離渦結(jié)構(gòu)。本文DMD方法提取的Mode2、Mode3和Mode4對(duì)應(yīng)的渦結(jié)構(gòu)與這一流動(dòng)規(guī)律吻合。

圖17為模態(tài)系數(shù)隨時(shí)間的變化，可以看出，不同系數(shù)之間也存在明顯的相位差。按照模態(tài)系數(shù)最大幅值對(duì)Mode2～Mode6進(jìn)行排序，發(fā)現(xiàn)排序后的模態(tài)對(duì)應(yīng)的St為St2、St4、St1、St3、St5，這與圖7(b)中的頻率按峰值排序一致。

(a) First order mode (b) Second order mode (c) Third order mode

(d) Fourth order mode (e) Fifth order mode (f) Sixth order mode

圖16DMD模態(tài)速度云圖
Fig.16CountersofvelocityforDMDmodes

3.2.3 深失速狀態(tài)

圖18為深失速狀態(tài)下一個(gè)周期中8個(gè)瞬態(tài)流場(chǎng)的渦量圖以及流線(xiàn)圖，流場(chǎng)中分離渦的演化規(guī)律與淺失速時(shí)類(lèi)似，與之不同的是，翼型上表面的分離區(qū)域進(jìn)一步增大。圖19為前六階DMD模態(tài)的速度云圖。

(a) First order mode (b) Second order mode (c) Third order mode

(d) Fourth order mode (e) Fifth order mode (f) Sixth order mode

圖19DMD模態(tài)速度云圖
Fig.19ContoursofvelocityforDMDmodes

其中Mode1依然反應(yīng)的是整個(gè)翼型上表面的大分離渦結(jié)構(gòu)。與淺失速狀態(tài)相比，Mode2、Mode3和Mode4進(jìn)一步表現(xiàn)為了隨迎角增大而形成的后緣分離渦和主渦上方區(qū)域的一些渦結(jié)構(gòu)，同時(shí)Mode5和Mode6主要表現(xiàn)為后緣分離渦和尾跡區(qū)的渦結(jié)構(gòu)。

圖20為模態(tài)系數(shù)隨時(shí)間的變化，可以看出，不同系數(shù)之間依然存在明顯的相位差。按照模態(tài)系數(shù)最大幅值對(duì)Mode2～Mode6進(jìn)行排序，發(fā)現(xiàn)排序后的模態(tài)對(duì)應(yīng)的St為St2、St1、St4、St3、St5，這與圖7(c)中的頻率按峰值排序一致。這說(shuō)明DMD提取的模態(tài)頻率不僅與FFT分析結(jié)果相差極小，并且兩者提取的頻率在流動(dòng)中的主導(dǎo)作用順序也一致。

4 結(jié) 論

本文基于SST-DDES混合方法對(duì)NACA0012大迎角分離流動(dòng)進(jìn)行非定常數(shù)值模擬，采用DMD方法對(duì)三種典型大迎角失速狀態(tài)下流場(chǎng)的渦穩(wěn)定性演化規(guī)律進(jìn)行研究。對(duì)于本文的計(jì)算模型及條件，得到以下結(jié)論：

1) DMD方法可以準(zhǔn)確地提取翼型大迎角流動(dòng)中的主頻和高階倍頻及對(duì)應(yīng)的流場(chǎng)模態(tài)結(jié)構(gòu)。與FFT分析結(jié)果相比，頻率最大相差小于0.16%，且可通過(guò)特征值對(duì)相應(yīng)的模態(tài)進(jìn)行穩(wěn)定性分析。由于三種迎角下所有模態(tài)的放大率均接近0，因此所有模態(tài)處于弱發(fā)散、弱收斂或極限環(huán)狀態(tài)。

2) 三種流動(dòng)狀態(tài)下，DMD提取的一階模態(tài)主要表現(xiàn)為分離渦演化過(guò)程中最主要的靜止分離渦結(jié)構(gòu)。前三階低頻對(duì)應(yīng)的模態(tài)渦結(jié)構(gòu)與流動(dòng)中以此頻率進(jìn)行演化的渦結(jié)構(gòu)比較一致。更高階的倍頻主要表現(xiàn)為尾渦和尾跡區(qū)的渦結(jié)構(gòu)。

3) DMD提取的模態(tài)頻率不僅與FFT分析結(jié)果相差極小，兩者提取的頻率在流動(dòng)中的主導(dǎo)作用順序也一致。不同DMD模態(tài)系數(shù)之間存在相位差， 說(shuō)明分離渦流動(dòng)中不同頻率對(duì)應(yīng)的渦結(jié)構(gòu)運(yùn)動(dòng)不同步。

參 考 文 獻(xiàn)：

[1]Alferez N, Mary I, Lamballais E. Study of stall development around an airfoil by means of high fidelity large eddy simulation[J]. Flow Turbulence Combust, 2013, 91: 623-641.

[2]Chen Pei. Simulation and analysis of airfoil deep stall flow at high angle of attack[C]//8th Cross-straits Conference on Aeronautics and Astronautics. BeiJing, 2012: 218-225. (in Chinese)陳培. 翼型大迎角深失速流動(dòng)的模擬和分析[C]//第八屆海峽兩岸航空航天學(xué)術(shù)研討會(huì). 北京, 2012: 218-225.

[3]Mellen C P, Fr-Ograve J, Hlich et al. Lessons from LESFOIL project on large-eddy simulation of flow around an airfoil[J]. AIAA Journal, 2012, 41(4): 573-581.

[4]Almutairi J H, Alqadi I M. Large-eddy simulation of natural low-frequency oscillations of separating-reattaching flow near stall conditions[J]. AIAA Journal, 2013, 51(4): 981-991.

[5]Spalart P R. Detached-eddy simulation[J]. Annual Review of Fluid Mechanics, 2009, 41: 181-202.

[6]Liu Z, Yang Y J, Zhou W J, et al. Study of unsteady separation flow around airfoil at high angle of attack using hybrid RANS-LES method[J]. Acta Aeronautica et Astronautica Sinica, 2014, 35(2): 372-380. (in Chinese)劉周, 楊云軍, 周偉江, 等. 基于RANS-LES混合方法的翼型大迎角非定常分離流動(dòng)研究[J]. 航空學(xué)報(bào), 2014, 35(2): 372-380.

[7]Richez F, Pape A L, Costes M. Zonal detached-eddy simulation of separated flow around a finite-span wing[J]. AIAA Journal, 2015, 53(11): 3157-3166.

[8]Mittal S, Saxena P. Prediction of hysteresis associated with static stall of an airfoil[J]. AIAA Journal, 2000, 38(5): 933-935.

[9]Zhang J Z, Li K L, Chen L Y. Nonlinear dynamics of static stall of airfoil and its control[J]. Acta Aeronautica et Astronautica Sinica, 2011, 32 (12): 2163-2173. (in Chinese)張家忠, 李凱倫, 陳麗鶯. 翼型失速的非線(xiàn)性動(dòng)力學(xué)特性及其控制[J]. 航空學(xué)報(bào), 2011, 32(12): 2163-2173.

[10]Lumley J L.The structure of inhomogeneous turbulence flows[J]. Atmospheric Turbulence and Radio Wave Propagation, 1967: 166-176.

[11]Sirovich L. Turbulent and the dynamics of coherent structures. Parts I-III[J]. Quarterly of Applied Mathematics, 1987, XLV(3): 561-582.

[12]Kang W, Dai X Y, Liu N. Multi-modal interaction and flow instability of flow around an airfoil at low reynolds number[J]. Journal of Northwestern Polytechnical University, 2015, 33(3): 382-387. (in Chinese)康偉, 代向艷, 劉凝. 低速翼型繞流的多模態(tài)耦合與流動(dòng)穩(wěn)定性研究[J]. 西北工業(yè)大學(xué)學(xué)報(bào), 2015, 33(3): 382-387.

[13]Schmid P J. Dynamic mode decomposition of numerical and experimental data[J]. Journal of Fluid Mechanical, 2010, 656: 5-28.

[14]Rowley C W, Mezic I, Bagheri S, et al. Spectral analysis of nonlinear flows[J]. Journal of Fluid Mechanical, 2009, 641: 115-127.

[17]Seena A, Sung H J. Dynamic mode decomposition of turbulent cavity flows for self-sustained oscillations[J]. International Journal of Heat and Fluid Flow, 2011, 32(6): 1098-1110.

[18]Muld T W, Efraimsson G, Henningson D S. Flow structures around a high-speed train extracted using proper orthogonal decomposition and dynamic mode decomposition[J]. Computers & Fluids, 2012, 57: 87-97.

[19]Pan C, Chen H, Wang J J. Dynamical mode decomposition of complex flow field[C]//8th National Conference on Experimental Fluid Mechanics. Guangzhou, 2010: 77-82. (in Chinese)潘翀, 陳皇, 王晉軍. 復(fù)雜流場(chǎng)的動(dòng)力學(xué)模態(tài)分解[C]//第八屆全國(guó)實(shí)驗(yàn)流體力學(xué)學(xué)術(shù)會(huì)議論文集. 廣州, 2010: 77-82.

[20]Mariappan S, Gardner A D, Richter K, et al. Analysis of dynamic stall using dynamic mode decomposition technique[J]. AIAA Journal, 2014, 52 (11): 2427-2439.

[21]Mohan A T, Gaitonde D V, Visbal M R. Model reduction and analysis of deep dynamic stall on a plunging airfoil using dynamic mode decomposition[R]. AIAA 2015-1058.

[22]Spalart P R. Young person’s guide to detached-eddy simulation grids[R]. NASA CR-2001-211032.

[23]Mccroskey W J. A critical assessment of wind tunnel results for the NACA 0012 airfoil[R]. NASA TM 1987-100019.