基于噴流擬序結(jié)構(gòu)預(yù)測的SGS 模型比較研究1)

劉琪麟 賴煥新

(華東理工大學(xué)機械與動力工程學(xué)院,上海 200237)

引言

擬序結(jié)構(gòu)普遍存在于湍流中[1-3],并影響湍流的發(fā)生、輸運與耗散等過程,對湍流擬序結(jié)構(gòu)的研究、對于理解湍流機理、分析流動相互作用有重要的意義.從大渦模擬來看,正確預(yù)測湍流擬序結(jié)構(gòu)的關(guān)鍵就在于算法和亞網(wǎng)格尺度(sub-grid-scale,SGS)模型的準(zhǔn)確性.

在算法方面,傳統(tǒng)的大渦模擬方法求解網(wǎng)格濾波的控制方程,并采用SGS 模型模化亞網(wǎng)格尺度流動與網(wǎng)格尺度流動之間的相互作用.雖然近年來部分學(xué)者利用數(shù)值離散的特點,也發(fā)展出不需要添加模型的隱式大渦模擬方法[4-5],但為了針對性地討論應(yīng)用中亞網(wǎng)格模型的選擇問題,本文仍采用高精度有限差分結(jié)合亞網(wǎng)格模型的大渦模擬經(jīng)典算法,并討論SGS 模型的性能.對大渦模擬算法更全面的介紹可以參考Sagaut 等[6]的綜述.

目前文獻中廣泛出現(xiàn)的亞網(wǎng)格尺度模型是經(jīng)典的常系數(shù)Smagorinsky 模型(Smagorinsky model,SM),已有眾多文獻報道該模型過度耗散了脈動動能[7-9].另一類是基于Germano 恒等式和Lilly 最小二乘法來動態(tài)調(diào)整模型系數(shù)的動態(tài)Smagorinsky 模型(dynamic Smagorinsky model,DSM),這類模型的局限性在于理論上要求在流場均勻的方向上計算平均的模型系數(shù),而這一平均過程缺少物理上合理的解釋,且造成分布式并行計算中不同處理器之間數(shù)據(jù)的頻繁調(diào)用與交換,降低計算效率.從計算的角度,較為理想的模型是局部模型,計算效率高且能較好地克服SM 對湍流黏性估計過高的缺點.Piomelli和Liu[9]提出了局部動態(tài)模型(localized dynamic Smagorinsky model,LDSM),去除了DSM 對流動中需存在各向同性方向的要求,而是利用流動的歷史信息來避免平均過程,計算效率得到提高.Lenormand等[10]提出的選擇多尺度模型(selective mixed-scale model,SMSM),在邊界層流動[10]、翼型繞流[11]和空腔流動[12]的預(yù)測中具有良好應(yīng)用.這一模型借助渦矢量方向的變化來判斷流動的狀態(tài)是否發(fā)生了轉(zhuǎn)捩,同時反映湍流的間歇性.Kobayashi 等[13-14]則利用流場中黏性耗散與擬序結(jié)構(gòu)的關(guān)聯(lián)性,提出了擬序結(jié)構(gòu)模型(coherent-structure Smagorinsky model,CSM),并應(yīng)用于各向同性湍流、湍流通道流以及橫向射流的計算中,其模型簡單且具有類似于動態(tài)模型的性能.Koyabashi[13]同時改進動能模型,提出擬序結(jié)構(gòu)動能模型(coherent-structure kinetic-energy model,CKM).這些局部模型都有望克服傳統(tǒng)模型的不足,但對于預(yù)測可壓縮噴流中的擬序結(jié)構(gòu),目前尚未看到有關(guān)它們預(yù)測性能的分析與評估的文獻.

典型的噴流擬序結(jié)構(gòu)有流向渦、高低速流體間的條帶[15]、剪切層邊沿的卷吸流動結(jié)構(gòu)[16]等,表現(xiàn)為長時間保持一定形狀的大尺度質(zhì)量團,且在時間和空間上具有強相關(guān)性[17].Lumley[18]最早提出采用本征正交分解(proper orthogonal decomposition,POD)方法提取擬序結(jié)構(gòu),且在自由剪切流、分離流、邊界層轉(zhuǎn)捩以及湍流通道流中得到廣泛的應(yīng)用.POD 方法的目標(biāo)是針對物理場尋找一組空間正交的基函數(shù),并按照基函數(shù)對物理量方差的貢獻排序,例如對于脈動速度場的研究就是基函數(shù)對湍動能的貢獻[19-20].本征正交分解的特點使得只需要少量的POD 模態(tài)就能夠提取出流場中能量占主導(dǎo)地位的流動模式.此外,POD 分解的基函數(shù)直接由湍流場確定,不同于傅里葉變換或者小波變換需要人為指定基函數(shù),因此表征的擬序結(jié)構(gòu)更為客觀[17].在POD 算法方面,雖然快照POD 算法能夠顯著減少計算量,但是有部分高階模態(tài)沒有參與計算.基于矩陣奇異值分解(singular value decomposition,SVD)的經(jīng)典POD 算法包含完整的模態(tài)信息,比快照POD算法對舍入誤差的魯棒性高[19-20].此外數(shù)學(xué)軟件Matlab 中提供了SVD 的庫函數(shù),也便于編制后處理程序.

本文對亞聲速可壓縮平行噴流進行大渦模擬.在驗證計算結(jié)果的基礎(chǔ)上,研究SM,CKM,SMSM,LDSM 和CSM 模型的特性,分析各個模型預(yù)測的瞬時渦結(jié)構(gòu).同時采用基于SVD 的POD 分解方法提取噴流擬序結(jié)構(gòu)并分析各模型預(yù)測結(jié)果的差異,為可壓縮噴流的大渦模擬提供參考.

1 數(shù)值計算與方法

1.1 亞網(wǎng)格尺度模型

控制方程組為Favre 質(zhì)量加權(quán)濾波的三維非穩(wěn)態(tài)可壓縮N-S 方程組[21].方程中的亞網(wǎng)格尺度未知量由亞網(wǎng)格模型給出,包括亞網(wǎng)格尺度應(yīng)力τij

以及亞網(wǎng)格尺度熱通量

其中ρ,T,νsgs以及ui(i=1,2,3) 分別為密度、溫度、模化的亞網(wǎng)格黏性系數(shù)以及速度矢量的3 個笛卡爾分量.γ 是量熱完全氣體的比熱比,空氣取γ=1.4.M為馬赫數(shù),湍流普朗特數(shù)常數(shù)Prt=1[22],而,頂標(biāo)“-”表示網(wǎng)格濾波,“～”表示Favre 質(zhì)量加權(quán)的網(wǎng)格濾波,所有計算變量均為網(wǎng)格濾波后的變量,因此除非特殊說明以下省去頂部符號.

Smagorinsky 模型(SM)中的亞網(wǎng)格尺度應(yīng)力為

δij=0.對于噴流,取模型系數(shù)C1=0.17[23].τkk的模化根據(jù)文獻[7]的建議得出.其中=Li j?δi jLkk/3,頂標(biāo)“^”表示測試濾波運算,Ct為上一時刻的模型系數(shù),初始值可以根據(jù)SM 模型取Ct=C1,當(dāng)前時刻的模型系數(shù)C2是未知量,且式(5)對C2是超定的.C2的估計值由式(5)兩邊取αij的縮并得出[9],如式(6)所示.LDSM 模型的亞網(wǎng)格尺度應(yīng)力與亞網(wǎng)格黏性系數(shù)通過將式(4)中的替換為得出

選擇多尺度模型(slective mixed-scale model,SMSM)是根據(jù)局部渦矢量方向的脈動來調(diào)整亞網(wǎng)格黏性系數(shù)

其中C3=0.06,α 取0.5[24].為了反映湍流的多尺度特性,引入小尺度脈動的動能=(1/2),其中脈動速度=?.進一步,對網(wǎng)格濾波尺度的瞬時渦量進行測試濾波,得出.和兩個向量之間的夾角 θ 反映了渦運動的脈動強弱.對于均勻各向同性湍流,在夾角 θ 的概率密度函數(shù)的峰值位置處有θ0=20°[24],據(jù)此,可判斷局部流動的狀態(tài).模型引入了選擇函數(shù)fθ0,當(dāng) θ> θ0時,取fθ0=1,表示渦方向在兩個濾波尺度上變化較大,是無規(guī)則的完全湍流狀態(tài);否則,當(dāng) θ<θ0時,取fθ0=tan(θ/2)/tan(θ0/2),意味著湍流逐步減弱,亞網(wǎng)格黏性也應(yīng)相應(yīng)減弱,直到在層流區(qū)僅保留分子黏性.因此,SMSM 模型的亞網(wǎng)格黏性系數(shù)為

擬序結(jié)構(gòu)模型是利用擬序結(jié)構(gòu)與湍流中耗散的相關(guān)性來調(diào)整亞網(wǎng)格黏性系數(shù)vsgs=CΔ2|S(u)|,模型系數(shù)C=C4|FCS|3/2中常數(shù)C4=0.05[13],FCS=Q/E為擬序結(jié)構(gòu)函數(shù)，其中

為速度導(dǎo)數(shù)Dij=?ui/?xj的第二不變量,它主要由Di j中反映當(dāng)?shù)匦D(zhuǎn)的反對稱部分Aij和反映當(dāng)?shù)刈冃蔚膶ΨQ部分Bi j組成,而為Di j的范數(shù).其中

基于直接數(shù)值模擬(direct numerical simulation,DNS),流場的相關(guān)性分析表明,Q的絕對值與擬序小尺度渦附近能量耗散的強弱存在關(guān)聯(lián)[13],因此|FCS|反映了流場中黏性與耗散的局部變化,而|FCS|的冪次“3/2”是根據(jù)不可壓縮壁面流動中C∝y3,Q∝y2以及E∝ 常數(shù)(y是壁面的法向)而確定的.因此,CSM 模型的亞網(wǎng)格黏性系數(shù)為

擬序結(jié)構(gòu)動能模型是在擬序結(jié)構(gòu)模型的基礎(chǔ)上,在亞網(wǎng)格黏性系數(shù)中顯式地包含亞網(wǎng)格尺度動能ksgs,并做量綱調(diào)整后得到

至此,可以通過式(2)、式(4)來計算亞網(wǎng)格尺度熱通量和亞網(wǎng)格尺度應(yīng)力.雖然文獻[12]指出亞網(wǎng)格尺度應(yīng)力中各向同性的部分 τkk相對于熱力學(xué)壓力可以忽略,但根據(jù)Vreman[7]的建議,在基于SM 的LDSM 以及CSM 模型中仍然模化 τkk.

1.2 計算對象與數(shù)值離散

本文采用Hu 等[22]的亞聲速平行噴流為考核算例.噴流出口雷諾數(shù)ReJ=DJρ∞UJ/μ∞=2000,馬赫數(shù)MJ=UJ/c∞=0.9,其中DJ=1、UJ=1 分別為歸一化的噴口寬度和噴口速度,參考值ρ∞,μ∞和c∞分別為來流密度、動力黏度以及環(huán)境聲速.計算域與網(wǎng)格劃分如圖1 所示,流向(x)、橫向(y)與展向(z)的尺寸分別為Lx×Ly×Lz=24×15×3,采用笛卡爾網(wǎng)格進行離散,對應(yīng)方向的結(jié)點數(shù)分別為 181×181×16.為了提高近噴口區(qū)域和噴流唇線附近的網(wǎng)格分辨率,采用雙曲正切函數(shù)對x和y方向的網(wǎng)格進行拉伸,而噴流在z方向采用周期邊界條件,劃分均勻的網(wǎng)格.基于Celik 等[25]提出的大渦模擬分辨率指標(biāo),由于課題組前期工作已證明該網(wǎng)格能夠分辨湍動噴流中75%以上的湍動能,考慮到噴流的雷諾數(shù)較低,該計算網(wǎng)格可以滿足大渦模擬的網(wǎng)格質(zhì)量要求[26].

圖1 計算域與網(wǎng)格(3 個方向都每隔3 條網(wǎng)格線畫一條線)Fig.1 Computional domain and the Grid(One in every three lines is shown in each direction)

采用有限差分方法開展數(shù)值模擬,本文的計算源程序以Sandham 課題組[27-28]的(shock-boundary layer interaction,SBLI)代碼為基礎(chǔ),計算精度和可靠性經(jīng)過了充分的考核驗證[22,29].在此基礎(chǔ)上,本文進一步采用OpenMP 并行化,使該代碼可用于內(nèi)存共享式的小型計算工作站.

空間離散的1 階與2 階導(dǎo)數(shù)均采用5 點4 階中心差分格式計算,邊界點采用Carpenter 提出的單側(cè)4 階穩(wěn)定差分格式[30],從而保持空間離散的整體高精度.為了抑制差分計算引起的非物理振蕩,算法上采用了熵分裂方法,通過將對流項分裂為守恒形式與非守恒形式的加權(quán)求和來實現(xiàn)[31].同時,在x和y方向的計算邊界上,采用無反射邊界條件來抑制非物理振蕩波在計算域邊界上的反射,算法上通過保留走出計算域的特征波而將走入計算域的特征波幅值置零來實現(xiàn)[32].在z方向上采用周期邊界條件來模擬平行噴流在展向的無限延伸.

時間離散采用低存儲顯式3 階龍格庫塔算法,這一算法完成3 個子時間步的推進僅需要兩個存儲流場的數(shù)組,對內(nèi)存資源的消耗低[21].為便于模型特性的研究比較,時間步長取定值Δt=0.003 075,且滿足數(shù)值穩(wěn)定性條件CFL<1.在初始條件方面,給定噴口平面上流向速度u、溫度T以及壓強p的形線,并用噴口平面上的物理量來初始化全場.噴口的速度分布采用雙曲正切函數(shù)給定,并施加伴流速度U0=0.1,可得

在數(shù)值結(jié)果方面,雖然DNS 能夠準(zhǔn)確給出流場結(jié)構(gòu)的信息,但是DNS 對計算資源的消耗量高,而且流場后處理的數(shù)據(jù)量大.以文獻[22-33]中提供的DNS 數(shù)據(jù)為基礎(chǔ),課題組前期工作已開展了計算代碼的考核[34-35]、LES 的網(wǎng)格分辨率分析[26]以及局部模型的研究[21,36],算法的可靠性與模型都已得到驗證.本文在已有工作的基礎(chǔ)上,進一步分析局部模型揭示的噴流擬序結(jié)構(gòu).

1.3 POD 方法

POD 方法采用一組正交基函數(shù) ψ(x) 的線性組合來表示物理場 φ(x,t),組合系數(shù)為c(t),即

所得的基函數(shù)也稱為POD 模態(tài).POD 分解的目的是使這種表示的誤差在平方平均的意義下最小,并且在給定誤差下最小化表示所需要的模態(tài)個數(shù),從而得到最佳的基函數(shù)[20].這組最佳基函數(shù)能夠最優(yōu)地表征原物理量的方差,對于脈動速度場，這種最優(yōu)是模態(tài)反映的湍動能隨模態(tài)階數(shù)增大而迅速的減小[19],脈動速度的第一階POD 模態(tài)對湍動能的貢獻最大.詳細(xì)的理論介紹可參考文獻[20].

在計算中,可對流場的速度分量進行時間采樣,得到關(guān)于標(biāo)量 φ 的N個流場樣本,對于網(wǎng)格點數(shù)為M=Mx×My的二維樣本,可將每個樣本都寫成M×1的列向量,由此得到M×N的樣本矩陣

在離散求解時,尋找最佳基函數(shù)的過程在數(shù)學(xué)上等價于求樣本矩陣F的SVD 問題[20].采用Matlab 中低存儲的“svd()”庫函數(shù)計算F的奇異值分解

其中,左矩陣 Ψ 的各列是正交的,Ψ=(ψ1,ψ2,···,ψN),它的第i列 ψi即為物理場的第i階POD 模態(tài),右矩陣C的各列也是正交的,C=(c1,c2,···,cN),它的第i列ci即為模態(tài) ψi對應(yīng)的系數(shù),反映POD 模態(tài)與時間有關(guān)的特征[19],而S為奇異值矩陣,S=dig(s1,s2,···,sN),S的對角元為從大到小排列的奇異值si(i=1,2,···,N),其他矩陣元為零,λi=代表了與第i階POD 模態(tài)對應(yīng)的能量,且對各階模態(tài)可以定義能量貢獻率

以及累計能量貢獻率

另外,由于各階POD 模態(tài) ψi是相互正交的,所以λ衰減越快,當(dāng)前模態(tài)所含有的能量相比下一階模態(tài)的差就越大,當(dāng)前模態(tài)所表征的物理場的相關(guān)性就越強.

2 SGS 模型的特性

2.1 平均流場

在分析湍流統(tǒng)計特征時,為排除初場對統(tǒng)計結(jié)果的影響,取無量綱時間t=48～ 216 的流場數(shù)據(jù),對應(yīng)流體以噴口速度流過2～ 9 倍的計算域長度,以計算時間步長Δt為采樣間隔,得到54 640 個流場樣本.圖2 以LDSM 模型為例,顯示了在唇線上無量綱位置x=2,6,8,10,20 處的流向脈動速度的功率譜,譜線具有一段固定的斜率且能觀察到慣性區(qū),反映了湍流的能量級串現(xiàn)象,也說明流動已經(jīng)達到了完全湍動狀態(tài).

圖2 LDSM 模型預(yù)測的脈動速度無量綱功率譜Fig.2 Dimensionless power spectra of velocity signals predicted by the LDSM

噴流主流流體與環(huán)境流體存在速度差,形成剪切層.圖3 所示為噴流剪切層的發(fā)展.考慮噴流中心線速度Uc與伴流速度U0,圖3(a)采用速度半值寬δ0.5,即流向速度u= 0.5(Uc+U0) 時的橫向位置y,顯示出剪切層的發(fā)展呈現(xiàn)兩個線性階段,即初期的緩慢增長段與下游的快速增長段,這對應(yīng)噴流勢流核末端x≈6 處的流動轉(zhuǎn)捩現(xiàn)象.對比5 種模型的結(jié)果,SM 模型預(yù)測的剪切層橫向發(fā)展程度較其他模型大,表明該模型過度耗散的特點.圖3(b)以LDSM 模型為例,在相似坐標(biāo)中繪出了噴流不同流向位置的速度形線,相似坐標(biāo)的橫向位置采用速度半值寬歸一化,平均噴流速度差〈u〉?Uc采用主、伴流的速度差ΔUc=Uc?U0歸一化,〈·〉 表示時間平均.圖3(b)中,速度形線在入口呈平頂狀分布,而在下游則快速發(fā)展為自相似的速度形線,可以看到x=4 時,自相似形線還沒有建立,但當(dāng)x≥8 時,速度形線與實驗噴流[37-39]的自相似數(shù)據(jù)點吻合良好,也驗證了計算結(jié)果.

圖3 噴流剪切層的發(fā)展Fig.3 Developments of the jet shear layer

噴流流體與環(huán)境流體的強烈摻混伴隨著能量的耗散.對比5 種SGS 模型的亞網(wǎng)格尺度平均黏性耗散率,定義為.圖4 為在展向中心平面(z=1.5)上的分布.由于噴流入口施加的平均來流條件,流體噴出后在近噴口區(qū)基本為層流狀態(tài),隨后在勢流核末端(x≈6)附近發(fā)生流動轉(zhuǎn)捩,并在下游發(fā)展為充分湍流.圖2 中,近噴口的脈動水平比勢流核及其下游的區(qū)域低1～ 2 個數(shù)量級,進一步表明來流的層流特征,此處與湍流脈動相關(guān)的亞網(wǎng)格尺度黏性耗散應(yīng)趨近于零.但在圖4(a)中,SM 模型預(yù)測的在近噴口區(qū)內(nèi)呈兩條平行的條帶,表明SM 模型對層流擾動敏感,不能正確反映湍流黏性耗散率的分布區(qū)域.圖4(b)～ 圖4(e)分別為CKM,SMSM,LDSM 以及CSM 的預(yù)測結(jié)果,這4 種局部模型雖然能夠正確地反映近噴口區(qū)的分布特征,在轉(zhuǎn)捩區(qū)(x=6～9)內(nèi)呈現(xiàn)一對峰值,但是以CSM 模型為參考,CKM,LDSM,SMSM 以及SM 模型預(yù)測的峰值耗散不同,分別為CSM 模型預(yù)測值的3.0,1.8,1.3 和0.9 倍.其中CKM 預(yù)測的峰值接近其他模型的1～ 2 倍,稍后分析的分布對轉(zhuǎn)捩區(qū)內(nèi)渦結(jié)構(gòu)的影響.

圖4 對比在中心平面(z=1.5)上的分布Fig.4 Comparisons of the distributi on of ?in the center plane of the jet (z=1.5)

2.2 瞬時流場的渦結(jié)構(gòu)

分析無量綱時間t=207 的瞬時流場,此時噴流以噴口速度流過8 倍以上的計算域長度,已經(jīng)充分排除初場的影響.參考Hu 等[22]的部分DNS 結(jié)果,以噴流中心線速度〈u〉 衰減為0.99 時的流向位置得出平均的勢流核長度為x=6.43.圖5 為Q=0.1 的等值面,圖5 中紅色切片位于Hu 等DNS 結(jié)果[22]的勢流核末端.在預(yù)測流動轉(zhuǎn)捩方面,圖5(c)～ 圖5(e)中,SMSM,LDSM 以及CSM 模型都能預(yù)測大尺度展向渦在平均勢流核末端附近的破碎,以及下游區(qū)域內(nèi)湍流多尺度渦結(jié)構(gòu)的生成.然而對于SM 模型,圖5(a)在切片位置下游仍然能夠觀察到流向的大尺度渦結(jié)構(gòu)(如圖5 中箭頭所示),結(jié)合圖4(a)可知,湍流中小尺度渦的發(fā)展可能受到近噴口過度的亞網(wǎng)格尺度耗散的抑制,這與文獻[7-9]對SM 模型過度耗散特點的報道一致.圖5(b)中,CKM 模型的預(yù)測結(jié)果與SM 類似,同樣可能由于該模型中局部亞網(wǎng)格尺度耗散較強.圖6 為展向渦 ωz的等值面,同樣采用紅色切片標(biāo)記DNS 結(jié)果[22]的勢流核末端位置.可以看到,ωz的等值面在勢流段內(nèi)呈平行的渦片,隨著流動向下游的發(fā)展,ωz等值面在勢流核末端附近快速卷起,隨后渦結(jié)構(gòu)被拉伸并在下游發(fā)展出復(fù)雜的三維結(jié)構(gòu).圖6(a)中,SM 模型預(yù)測的渦卷起位置位于切片之后,渦片偏長,對應(yīng)圖4(a)中近噴口區(qū)的不合理耗散,SM 模型抑制了流場中展向渦的發(fā)展.圖6(b),圖6(d),圖6(e) 中,CKM,LDSM 以及CSM 模型在勢流核末端位置附近都能夠觀察到展向渦的卷起(如圖6 中箭頭標(biāo)記)以及下游的渦破碎,但LDSM 與CSM 模型在切片后預(yù)測的小尺度渦結(jié)構(gòu)比CKM 模型豐富,再次表明CKM 模型的局部強耗散可對小尺度渦結(jié)構(gòu)產(chǎn)生抑制.圖6(c) 中,SMSM 預(yù)測的渦卷起不顯著,但切片下游的預(yù)測結(jié)果與LDSM,CSM類似.相比SM 與CKM 模型,SMSM,LDSM 以及CSM 模型更能反映勢流核下游展向渦的多尺度特性.

圖5 對比 Q=0.1 的等值面,采用流向速度著色,紅色切片位置標(biāo)記了文獻[22]中DNS 結(jié)果的勢流核末端:x=6.43Fig.5 Comparisons of the iso-surface of Q=0.1,colored by streamwize velocity.The red slice marks the end of the potential core x=6.43 of the DNS results in Ref.[22]

圖6 對比 ωz 的等值面,紅色切片位置標(biāo)記了文獻[22]中DNS 結(jié)果的勢流核末端:x=6.43Fig.6 Comparisons of the iso-surface of ωz,the red slice marks the end of the potential core x=6.43 of the DNS results in Ref.[22]

2.3 速度場的POD 模態(tài)分析

對5 種SGS 模型的預(yù)測結(jié)果進行POD 分解.平行噴流主要表現(xiàn)為x-y平面內(nèi)的二維流動特征,因此在湍流統(tǒng)計的流場數(shù)據(jù)中,截取三維流場的中心平面(z=1.5)為樣本,以10Δt(=0.030 75) 為采樣間隔,得到4685 個樣本.針對5 種SGS 模型,提取脈動速度分量u′(=u?〈u〉),v′(=v?〈v〉) 和w′(=w?〈w〉)的POD 模態(tài).脈動速度u′,v′和w′的第一階模態(tài)分別對原始脈動速度場的脈動強度:〈u′u′〉,〈v′v′〉 和〈w′w′〉的貢獻最大,并在這一意義下脈動速度的第一階模態(tài)與原始脈動速度場最為接近.圖7 為流向脈動速度u′模態(tài)的累計能量貢獻率 η,可以看到,累計能量貢獻率曲線在前200 階內(nèi)就達到了90%以上,這200 階模態(tài)僅為總模態(tài)數(shù)的4%,因此曲線 η 反映了POD 的快速收斂特性.但5 種SGS 模型的 η 具有不同的增長速度,SM 模型最快,而SMSM 模型最慢.脈動速度v′與w′的 η 曲線具有與圖7 相似的快速收斂特征,因此不再繪出.圖8 為脈動速度分量在前20 階模態(tài)的能量貢獻率衰減曲線 ξ.可以看到,5 種SGS 模型得出的 ξ 曲線都在前5 階內(nèi)迅速衰減,說明低階模態(tài)占主導(dǎo)且模態(tài)表征的流場具有強的相關(guān)性.對比圖8 中第1 階模態(tài)的能量貢獻率,圖8(a)中u′模態(tài)的衰減曲線顯示SM 模型的第1 階模態(tài)貢獻率最高,SMSM 模型最低;圖8(b)中橫向速度脈動v′模態(tài)的衰減曲線顯示CKM 模型的貢獻率最高,而SMSM 模型依然排最低;圖8(c)展向速度脈動w′模態(tài)的衰減曲線中SMSM 模型的能量貢獻率排序有所提高.圖8 表明,第1 階模態(tài)的能量貢獻率已經(jīng)體現(xiàn)出了不同模型之間的差異,接下來針對第1 階模態(tài)展開分析.

圖7 u′ 模態(tài)的累計能量貢獻率的收斂曲線Fig.7 Convergency curves of the cumulative energy contribution for the modes of u′

噴流主要表現(xiàn)為流向的動量傳遞.對比5 種SGS 模型預(yù)測的流向速度脈動u′的第一階POD 模態(tài),圖9 顯示了模態(tài)幅值在中心平面(z=1.5)上的分布,采用接近零的模態(tài)幅值等值線u′=±0.003 來顯示模態(tài)的輪廓.圖9(a)中,SM 得到兩條沿流向的帶狀模態(tài),且在噴流下游x=12～24 的范圍內(nèi)一直保持高的模態(tài)幅值,對應(yīng)圖5(a)中大尺度流向渦結(jié)構(gòu)對流場的攪動.圖9(b)中,CKM 模型預(yù)測的模態(tài)也在下游保持沿流向的高幅值帶狀模態(tài),在轉(zhuǎn)捩區(qū)(x=6～9)內(nèi)的模態(tài)強度呈現(xiàn)一對低谷區(qū),位置恰好匹配圖4(b)中的峰值耗散區(qū)位置.這表明流向速度脈動的POD 模態(tài)對亞網(wǎng)格尺度的黏性耗散率敏感,局部的強耗散導(dǎo)致脈動強度〈u′u′〉 的降低,從而改變了POD 模態(tài)的形狀.圖9(d)、圖9(e)分別為LDSM與CSM 模型預(yù)測的模態(tài),兩者都沿流向兩側(cè)發(fā)散開來,并且在噴流下游呈現(xiàn)分岔或破碎的形態(tài),模態(tài)幅值也逐漸減小,反映下游的多尺度結(jié)構(gòu)以及脈動能量的逐漸消散.圖9(c)為SMSM 模型預(yù)測的模態(tài),模態(tài)輪廓在流向上呈斷續(xù)的條狀,局部的模態(tài)幅值低,這可能與圖8(a)中SMSM 模型第一階POD 模態(tài)的能量貢獻率低有關(guān).而在下游,SMSM 的模態(tài)呈現(xiàn)與CSM,LDSM 模型預(yù)測結(jié)果相似的分岔形態(tài).因此,與SM 和CKM 模型相比,CSM,LDSM 以及SMSM 模型都合理地反映了u′模態(tài)的形態(tài).

圖8 脈動速度分量的POD 模態(tài)的能量貢獻率Fig.8 Energy contribution rate of POD modes for the fluctuate velocity components

圖9 流向脈動速度 u′ 的第一階POD 模態(tài)云圖,實、虛等值線分別表示 u′=+0.003和u′=?0.003Fig.9 Contours of the first-order POD modes of u′,with solid contour lines for u′=+0.003 and dashed contour lines foru′=?0.003

圖10 為橫向脈動速度v′的第一階POD 模態(tài),等值線v′=±0.003 顯示出模態(tài)的輪廓.同時,圖10中采用面內(nèi)矢量 (u′,v′) 的POD 模態(tài)來顯示噴流中的流動模式.觀察圖10,v′的POD 模態(tài)呈正負(fù)相間的肋狀,并沿流向排列.在肋狀模態(tài)上下端部的附近,矢量箭頭顯示出環(huán)形的流動模式,在噴流上部的紅圓標(biāo)出了逆時針環(huán)流,下部的藍(lán)圓標(biāo)出了順時針環(huán)流.環(huán)流模式穿過了主流與環(huán)境流體,在圓環(huán)的下游側(cè),高速的主流流體被送入兩側(cè)的低速環(huán)境流體中,而在環(huán)流的上游側(cè),低速的環(huán)境流體被帶入主流流體中.同時,v′模態(tài)輪廓的橫向尺寸也沿著流向逐漸增大.可見,(u′,v′) 的POD 模態(tài)反映了噴流的流動卷吸現(xiàn)象[16].對比圖10(a)～ 圖10(e),SM,SMSM,LDSM 以及CSM模型預(yù)測的模態(tài)都反映了噴流的卷吸現(xiàn)象,其中SMSM 模型預(yù)測的旋渦尺度小,這與圖8(a) 和圖8(b) 中的低模態(tài)貢獻率對應(yīng).而圖10(b)中,CKM模型預(yù)測的模態(tài)基本沒有反映該區(qū)域內(nèi)的流動卷吸,結(jié)合圖4(b),脈動速度分量u′和v′可能受到CKM模型在x=6 附近局部的強耗散抑制,導(dǎo)致未能預(yù)測出明顯的流動卷吸.

圖10 橫向脈動速度v'的第一階POD 模態(tài)云圖與面內(nèi)矢量 (u′,v′) 的第一階POD 模態(tài)矢量圖(實、虛等值線分別表示v′=+0.003和 v′=?0.003,紅圈標(biāo)示逆時針環(huán)流、藍(lán)圈標(biāo)示順時針環(huán)流)Fig.10 Contours of the first-order POD modes of v' and the in plane vector (u′,v′)(Solid contour lines for v′=+0.003 and dashed contour lines for v′=?0.003.Red circle marks counterclockwise circular flow and blue circle indicates clockwise circular flow)

圖11 為展向脈動速度w′的POD 模態(tài)云圖,等值線w′=±0.003 顯示出模態(tài)的輪廓.圖11(b)～圖11(e)分別為CKM,SMSM,LDSM 以及CSM 模型預(yù)測的結(jié)果,w′的POD 模態(tài)在轉(zhuǎn)捩區(qū)(x=3～9)內(nèi)呈現(xiàn)規(guī)則的脊?fàn)钆帕?紅色表示脈動速度w′指向平面外,藍(lán)色表示脈動速度w′指向平面內(nèi),強弱相間的分布反映了w′模態(tài)對渦結(jié)構(gòu)的展向拉伸,如圖11(f)所示.圖11 中的脊?fàn)钅B(tài)顯示渦拉伸主要發(fā)生在轉(zhuǎn)捩區(qū),而在下游的模態(tài)輪廓呈碎塊狀,強度也逐漸減弱,對應(yīng)渦結(jié)構(gòu)逐漸發(fā)展為三維結(jié)構(gòu),展向的拉伸現(xiàn)象不再顯著.圖11(a)中,SM 模型預(yù)測的模態(tài)沒有反映轉(zhuǎn)捩區(qū)內(nèi)的脊?fàn)钆帕?在下游呈大團塊,位置與圖4(a) 中SM 模型噴口附近的過度耗散區(qū)域相對應(yīng).圖11(b)中,CKM 模型預(yù)測的模態(tài)雖然呈現(xiàn)了脊?fàn)钆帕?但是模態(tài)輪廓在下游x=21 處仍然呈大尺度塊狀,與SM 模型的預(yù)測結(jié)果相似.相比而言,SMSM,LDSM 以及CSM 模型的預(yù)測結(jié)果則能夠合理地反映展向拉伸的流動模式.

圖11 展向脈動速度 w′ 的第一階POD 模態(tài)云圖Fig.11 Contours of the first-order POD modes of w′,with solid contour lines for w′=+0.003and dashed contour lines forw′=?0.003

圖11 展向脈動速度 w′ 的第一階POD 模態(tài)云圖 (續(xù))Fig.11 Contours of the first-order POD modes of w′,with solid contour lines for w′=+0.003and dashed contour lines for w′=?0.003 (continued)

在計算效率方面,表1 對比了5 種SGS 模型的CPU 時間,以SM 模型為基準(zhǔn),LDSM,SMSM,CKM,CSM 的計算時間分別為SM 模型的:1.64,1.28,1.22 和1.00 倍,其中LDSM,SMSM 以及CKM 模型需要進行濾波運算,因此需要更多的CPU 時間,CSM 模型的計算時間最接近SM 模型.

表1 每時間步內(nèi)每網(wǎng)格點所占用的CPU 時間Table 1 CPU time per grid point per time step

3 結(jié)論

本文開展了亞聲速可壓縮平行噴流的大渦模擬.在分析平均流動與耗散的基礎(chǔ)上,針對瞬時渦流以及脈動速度POD 模態(tài)表征的擬序結(jié)構(gòu),對SM,CKM,SMSM,LDSM 和CSM,5 種SGS 模型的性能進行研究,得出以下主要結(jié)論:

(1) 初步揭示了u′模態(tài)從流向的帶狀到下游分岔或破碎形態(tài)反映的多尺度特征、v′模態(tài)沿流向的肋狀排列與尺寸增長、矢量 (u′,v′) 模態(tài)的環(huán)流模式表征的流動卷吸、以及w′模態(tài)表征的流場結(jié)構(gòu)在展向上受拉伸的模式.

(2) SMSM,LDSM 以及CSM 均較好地反映出湍流的多尺度特性,清晰地分辨了小尺度結(jié)構(gòu)的發(fā)展過程.而SM 與CKM 模型則未能有效地預(yù)測湍流中的小尺度渦結(jié)構(gòu)及部分流場特征模態(tài).

(3) 脈動速度場的POD 模態(tài)對SGS 模型的亞網(wǎng)格尺度耗散敏感,表現(xiàn)為CKM 預(yù)測的峰值耗散區(qū)對應(yīng)u′模態(tài)的低谷區(qū)因而環(huán)流模式不顯著,以及SM 未能反映轉(zhuǎn)捩區(qū)內(nèi)w′模態(tài)的脊?fàn)罾炷Ｊ?SMSM,LDSM 和CSM 克服了以上模型的不足之處,其中CSM 模型同時兼有計算效率較高的優(yōu)勢.