雷諾應(yīng)力與渦黏性模型的分離流預(yù)測對比分析

趙雅甜，邵志遠(yuǎn)，閻超，向星皓

1.中南大學(xué) 交通運輸工程學(xué)院，長沙 410083

2.北京航空航天大學(xué) 航空科學(xué)與工程學(xué)院，北京 100191

3.中國空氣動力研究與發(fā)展中心，綿陽 621000

流動分離與翼型和機翼氣動特性息息相關(guān)，飛行包線也時刻受到流動分離的影響?，F(xiàn)代飛行器愈發(fā)苛刻和復(fù)雜的飛行條件使得設(shè)計過程中對分離流動的準(zhǔn)確預(yù)測愈發(fā)關(guān)鍵。然而，對于較大逆壓梯度下的分離流動，精確的CFD（Com‐putational Fluid Dynamics）模擬仍十分困難［1-2］。包含展向流動的三維邊界層分離和高速下激波誘導(dǎo)邊界層分離等典型復(fù)雜流動現(xiàn)象對預(yù)測模型的精度提出了更高要求［3-4］。

湍流模型作為RANS（Reynalds Averaged Navier-Stokes Equations）方法的核心之一，對流場計算結(jié)果影響顯著。第2 次阻力預(yù)測會議（The Second DPW）的分析結(jié)果顯示，湍流模型對CFD 計算結(jié)果的影響最大，占到約15%［5］。湍流場中各參數(shù)脈動幾乎均與幾何邊界相關(guān)，因此，當(dāng)前尚不存在某種湍流模型可適用于全部流場。被廣泛使用的渦黏性模型由于Boussinesq假設(shè)的“先天缺陷”，導(dǎo)致對分離流動預(yù)測失準(zhǔn)［6］。不同于渦黏性模型，雷諾應(yīng)力模型從雷諾應(yīng)力輸運方程出發(fā)，對雷諾應(yīng)力直接求解，考慮了雷諾應(yīng)力的對流和擴散和流動歷史效應(yīng)［7］。因此該模型被認(rèn)為是最自然、最符合邏輯的模型［8］。NASA（National Aeronautics and Space Admin‐istration）在CFD 2030 遠(yuǎn)景中也將雷諾應(yīng)力模型作為主要發(fā)展的RANS 方法之一［9］?；谠擃惸Ｐ停琑umsey 等［10］通過對亞聲速駝峰算例和ONERA M6 機翼算例進行計算后，發(fā)現(xiàn)雷諾應(yīng)力模型存在再附點附近流線后彎，但預(yù)測精度優(yōu)于渦黏性模型。對于小攻角來流下的三維機翼激波位置預(yù)測較為準(zhǔn)確。王圣業(yè)等［11］利用基于雷諾應(yīng)力模型的分離渦模擬方法對 24.6°迎角鈍前緣三角翼進行數(shù)值仿真，發(fā)現(xiàn)模擬結(jié)果明顯好于傳統(tǒng)基于線性渦黏模型的分離渦模擬方法。舒博文等［12］發(fā)現(xiàn)雷諾應(yīng)力模型可在翼身接合區(qū)角區(qū)流動和三維強激波誘導(dǎo)分離等問題中得到正確的流動特征。

然而，雷諾應(yīng)力模型當(dāng)前也存在諸多問題。熊莉芳等［13］研究發(fā)現(xiàn)，雷諾應(yīng)力模型對邊界條件十分敏感，對計算網(wǎng)格要求苛刻，數(shù)值穩(wěn)定性較差。為解決數(shù)值剛性和邊界條件奇性問題，Men‐ter 和Egorov［14］提出Stress-BSL 模型，兼顧了計算效率，在常規(guī)流動預(yù)測中表現(xiàn)良好，但缺乏在強逆壓梯度下分離流動的具體研究及其與渦黏性模型的對比。本文基于Stress-BSL 模型，分析雷諾應(yīng)力模型在二維駝峰（Hump）逆壓梯度分離、跨聲速凸塊（Bump）激波-邊界層干擾分離、M6 跨聲速機翼三維邊界層分離中的預(yù)測表現(xiàn)。基于模型所反映的物理現(xiàn)象和數(shù)學(xué)構(gòu)造，對預(yù)測分離流動時出現(xiàn)誤差的原因及機制進行分析探究。通過與渦黏性模型進行對比，預(yù)測其性能優(yōu)勢，為未來湍流模型的改進和優(yōu)化提供參考。

1 湍流模型與算法介紹

雷諾應(yīng)力模型直接求解雷諾應(yīng)力輸運方程，本文選取 Stress-BSL 模型作為雷諾應(yīng)力模型的代表，將長度尺度ε與長度尺度ω進行混合，該模型的構(gòu)造為

雷諾應(yīng)力輸運方程為

式中：ρ為密度；為平均速度；Reij為雷諾應(yīng)力；μ為分子黏度性；μt為渦黏性系數(shù)。雷諾應(yīng)力再分配項作為輸運方程的重中之重，該項的構(gòu)造為

式中：k為湍動能；δij為克羅內(nèi)克符號；C1、α、β、γ為模式常數(shù)；Sij為平均應(yīng)變率；Wij為角速度；aij為各向異性張量。Sij、Wij、aij3個參數(shù)的計算公式為

耗散項采用各向同性假設(shè)進行?；?，可表示為

ω輸運方程為

本文空間離散方法針對有限體積法求解采用二階中心差分格式，針對黏性通量插值算法采用二階迎風(fēng) MUSCL（Monotone Upstreamcentred Schemes for Conservation Laws）重構(gòu)格式，無黏通量采用 HLLC（Harten-Lax-van Leer-Contact）格式。時間離散格式采用隱式LU-SGS（Lower-Upper Symmetric Gauss-Seidel）方法［15］。

2 數(shù)值計算與結(jié)果分析

2.1 二維駝峰（Hump）逆壓梯度分離

二維駝峰（Hump）模型［16-17］的前半部分上凸，流道收縮，流體在此處加速，壓強減小，后半部分突然下凹，流道擴張，出現(xiàn)較強的逆壓梯度，邊界層發(fā)生分離。由于駝峰下游的壓力梯度接近于0，因此邊界層可以再附于下壁面，分離區(qū)位于駝峰后半部分下凹位置。

圖1 給出了Hump 模型的計算網(wǎng)格，對邊界層和駝峰后半部分進行了網(wǎng)格加密，網(wǎng)格分辨率為1 633×433（流向×法向）。來流馬赫數(shù)為0.1，來流總溫與參考溫度的比值為1.007，來流總壓與參考壓強的比值為1.002；出口靜壓與參考壓強的比值為0.999 62。

圖1 Hump 模型計算網(wǎng)格Fig.1 Calculation grids of hump

圖2 給出了2 種模型計算的駝峰后緣流線和無量綱流向速度（U/U∞）分布，其中，U為當(dāng)前位置流速；U∞為自由來流速度；X、Y方向位置采用弦長（c）無量綱化?？梢园l(fā)現(xiàn)分離區(qū)上緣流動存在明顯的剪切，Stress-BSL 模型預(yù)測結(jié)果誤差更小，具體的流動分離和再附位置如表1。但是在再附點附近可以明顯地發(fā)現(xiàn)，Stress-BSL 模型預(yù)測的流線出現(xiàn)了不符合物理規(guī)律的向前彎曲，其他雷諾應(yīng)力模型也被發(fā)現(xiàn)存在類似的問題［4］。這可能是由于雷諾應(yīng)力模型在此處預(yù)測的湍流長度尺度增長率過大，導(dǎo)致剪切區(qū)域本應(yīng)向下游發(fā)展的流線受到分離區(qū)剪切影響突然增大，出現(xiàn)流線彎曲。

圖2 Hump 模型后緣無量綱流向速度與流線圖Fig.2 Nondimensional velocity and streamline behind hump

2.2 跨聲速凸塊（Bump）激波誘導(dǎo)邊界層分離

跨聲速凸塊（Bump）流動模型［18］由一個環(huán)形圓弧凸起構(gòu)成，氣體流經(jīng)凸起速度增加產(chǎn)生激波，在凸起的后半部分發(fā)生激波誘導(dǎo)分離，隨后邊界層再附并在凸起后緣產(chǎn)生分離區(qū)。圖3 給出了軸對稱跨聲速凸塊的計算網(wǎng)格，對邊界層、激波和環(huán)形凸塊后緣部分進行了網(wǎng)格加密，網(wǎng)格分辨率為721×321（流向×法向）。為模擬軸對稱條件，需將二維平面網(wǎng)格繞X軸旋轉(zhuǎn)1°，并將橫向邊界面設(shè)置為周期邊界。來流馬赫數(shù)為0.875，來流單位雷諾數(shù)為2.76×106m?1。

圖3 Bump 模型計算網(wǎng)格Fig.3 Calculation grids of bump

圖4 Bump 模型后緣無量綱流向速度與流線圖Fig.4 Nondimensional velocity and streamline behind bump

2 種模型預(yù)測所得流場流線和無量綱流向速度（U/U∞）分布見圖 4，可以明顯地發(fā)現(xiàn) SST 模型預(yù)測的流動再附點位置過于靠后，分離區(qū)過長，具體的流動分離點和再附點位置見表2。2 種模型的計算結(jié)果均出現(xiàn)了流動分離點提前，但Stress-BSL 模型的誤差更小。對于再附點，SST 模型預(yù)測的結(jié)果靠后，而Stress-BSL 模型的結(jié)果則較實驗值靠前。由于凸塊模型高度較小，再附點附近流線曲率很小，因此前述的原因?qū)е铝骶€向前彎曲，易使流動提前再附。整體來看 Stress-BSL 模型預(yù)測效果更好。

表2 Bump 模型分離、再附點位置Table 2 Separation and reattachment points of bump

表3 網(wǎng)格無關(guān)性驗證Table 3 Grid independence verification

2.3 跨聲速 M6 機翼三維邊界層分離

ONERA M6 機翼是跨聲速三維激波誘導(dǎo)邊界層分離的經(jīng)典算例［19］。在跨聲速來流條件下，機翼的上表面會出現(xiàn)激波，并且由于三維邊界層存在展向流動，在靠近翼尖的區(qū)域激波增強并出現(xiàn)激波誘導(dǎo)分離。研究表明，2 種湍流模型對于小攻角下的流動預(yù)測均與實驗吻合良好，因此本文采用大攻角6.06°進行強激波誘導(dǎo)分離的預(yù)測能力對比，來流馬赫數(shù)為0.84，來流單位雷諾數(shù)為1.46×107m?1。

圖5 給出了跨聲速三維 ONERA M6 機翼模型的計算網(wǎng)格，采用 O-H 型拓?fù)?，由于沒有側(cè)滑角，故采用半模進行計算，并對邊界層、機翼前緣部分和機翼翼尖段進行了網(wǎng)格加密，總共設(shè)計了粗、中、細(xì)、極細(xì)4 種網(wǎng)格，半模網(wǎng)格量分別約為133 萬、265 萬、530 萬、795 萬，網(wǎng)格細(xì)節(jié)見表 3。

圖5 M6 計算網(wǎng)格Fig.5 Calculation grids of M6

圖6 為采用 Stress-BSL 湍流模型計算的4 套 ONERA M6 機翼網(wǎng)格在6.06°攻角下20%展長位置處的壓力系數(shù)（Cp）分布計算結(jié)果與實驗值的對比結(jié)果。發(fā)現(xiàn)細(xì)網(wǎng)格和極細(xì)網(wǎng)格對激波位置的捕捉基本一致，且與實驗值誤差較小，為節(jié)省計算資源，后續(xù)均采用細(xì)網(wǎng)格進行數(shù)值模擬。

圖6 網(wǎng)格無關(guān)性驗證Fig.6 Grid independence verification

圖7、表4 給出了2 種湍流模型的計算結(jié)果對比。由圖7（a）、圖7（b）可以發(fā)現(xiàn)，在機翼的背風(fēng)面存在明顯的“λ”型激波。對比2 種模型的預(yù)測結(jié)果，Stress-BSL 模型預(yù)測的激波結(jié)構(gòu)更為合理。對于激波誘導(dǎo)邊界層分離而言，波后高壓前傳導(dǎo)致波前邊界層增厚，流動受阻進一步導(dǎo)致了激波的增強，因此流動分離的位置與激波位置幾乎重合。隨著逐漸靠近翼尖，激波增強，Stress-BSL 模型預(yù)測的結(jié)果始終與實驗值符合較好，而 SST 模型預(yù)測的激波位置開始前移，分離區(qū)過大（見圖7（c）），與實驗值的誤差增大。過大的分離區(qū)導(dǎo)致流動再附受阻，在機翼尾緣附近壓力系數(shù)明顯偏大（見圖8）。需要注意的是，圖8（a）、圖8（b）中Stress-BSL 模型預(yù)測的激波位置反而較SST 模型誤差偏大。這是由于風(fēng)洞實驗中翼根平面并不是對稱面，而CFD 計算中將翼根平面設(shè)置為對稱面，因此產(chǎn)生了邊界條件干擾，使得激波預(yù)測位置較實驗值后移。但由于SST 模型預(yù)測的激波位置靠前，因此導(dǎo)致了圖8（a）、圖8（b）的SST 模型預(yù)測激波位置更準(zhǔn)的巧合性結(jié)果。這種邊界條件干擾會隨著逐漸遠(yuǎn)離翼根平面而消失。Rumsey［4］研究也提到了這一點。

表4 升力系數(shù)與阻力系數(shù)對比Table 4 Comparison of lift and drag coefficient

圖7 不同湍流模型計算結(jié)果對比Fig.7 Comparison of simulation results with different turbulence models

圖8 不同站位壓力系數(shù)Fig.8 Pressure coefficient at different stations

另外，由于Stress-BSL 模型預(yù)測的流動分離更加準(zhǔn)確，分離區(qū)更小，背風(fēng)面后緣的壓力系數(shù)較SST 模型的預(yù)測結(jié)果更小，使得機翼上下表面間的壓強差更大，且壓差阻力更大，因此相較于SST 模型的預(yù)測結(jié)果，升力系數(shù)和阻力系數(shù)均更大（見表4）。

3 湍流模型性能機制

基于計算結(jié)果，本節(jié)著手從湍流模型的構(gòu)造機制出發(fā)，探明預(yù)測分離點提前、再附點滯后、分離區(qū)過大的具體原因。雷諾應(yīng)力的物理含義為流體微元表面上脈動動量輸運的平均值，可以認(rèn)為對于湍流邊界層而言，更小的雷諾應(yīng)力意味著邊界層抵抗逆壓梯度的能力越差，越容易出現(xiàn)流動分離。因此，湍流模型對雷諾應(yīng)力的預(yù)測能力直接決定了對分離流動預(yù)測的準(zhǔn)確性。

3.1 Stress-BSL 模型

圖9 給出了Hump 模型、Bump 模型特點站位的Stress-BSL 模型雷諾應(yīng)力計算值與風(fēng)洞實驗值對比?？梢园l(fā)現(xiàn)，Stress-BSL 模型的預(yù)測值較實驗值偏低，導(dǎo)致了預(yù)測分離提前和再附滯后。雖然通過直接描述雷諾應(yīng)力的輸運過程避免了引入假設(shè)的誤差，但雷諾應(yīng)力輸運方程中難以實驗觀測和?；捻棧o建模帶來了未知性和誤差，使得 Stress-BSL 模型雖然具有優(yōu)勢，但仍存在較大誤差。分析雷諾應(yīng)力輸運方程中的各項，可以得到 Stress-BSL 模型預(yù)測不準(zhǔn)的原因。

圖9 Stress-BSL 模型雷諾應(yīng)力計算結(jié)果與實驗值對比Fig.9 Comparison of stress-BSL and experimental Reynolds stress

1）生成項

由于雷諾應(yīng)力模型直接求解雷諾應(yīng)力輸運方程，因此可以對生成項精確求解。但精確并不是絕對準(zhǔn)確，在求解雷諾應(yīng)力時由于模型建模誤差導(dǎo)致的計算誤差會發(fā)生累積。因而，所謂的精確是指在生成項求解過程中沒有引入新的變量和假設(shè)，不會產(chǎn)生新的誤差。

2）擴散項

Stress-BSL 模型中擴散項將雷諾應(yīng)力的擴散速率與其在流場中的梯度表示為正比。De‐muren 和Sarkar［20］通過對平面通道湍流研究后發(fā)現(xiàn)：擴散項在對數(shù)層對雷諾正應(yīng)力的計算結(jié)果影響不大，但可以對邊界層外層湍流轉(zhuǎn)化為各向同性產(chǎn)生影響。這一項并不是雷諾應(yīng)力模型預(yù)測帶有強烈剪切的分離流動準(zhǔn)確與否的決定性因素，絕大多數(shù)文獻并未對其影響進行專門研究。

3）耗散項

湍流耗散的情況十分復(fù)雜，包括由大渦和小渦拉伸、壓強作用、湍流輸運、分子擴散造成的耗散作用，以及湍動能耗散自身的耗散作用。由于耗散作用基本發(fā)生在湍流的最小尺度中，此時的耗散作用除近壁區(qū)等特殊區(qū)域外，可以認(rèn)為是各向同性的。Stress-BSL 模型中擴散項采用了各向同性的耗散模型［21］。由圖9 也可以發(fā)現(xiàn)，在近壁區(qū)雷諾應(yīng)力的計算值和實驗誤差很小，因此該假設(shè)基本合理。但需要注意的是，耗散項中的模式封閉參數(shù)事實上需要根據(jù)雷諾數(shù)不同進行經(jīng)驗修正，在不同的流動區(qū)域內(nèi)取值不同，且耗散項與湍流長度尺度息息相關(guān)，對雷諾應(yīng)力模型的預(yù)測能力影響較大。

4）再分配項

雷諾應(yīng)力再分配項是脈動壓強和應(yīng)變率張量相關(guān)的平均值［22-24］。由于湍動能輸運方程中并不存在再分配項，說明再分配項對于湍動能的產(chǎn)生沒有貢獻，而只是在湍流脈動速度的各個分量之間進行調(diào)節(jié)。不同于擴散作用依賴于明確的梯度，再分配項是一種能量轉(zhuǎn)移，通過將量值較大的速度分量轉(zhuǎn)移給量值較小的速度分量，以達到各方向分量的平衡。為明確其作用，給出二維湍流場的雷諾應(yīng)力輸運方程

可以發(fā)現(xiàn)，Y方向的雷諾正應(yīng)力生成項為0。假設(shè)在均勻剪切場中初始時只有正應(yīng)力，沒有切應(yīng)力。由于流場的剪切作用，雷諾切應(yīng)力開始增長，導(dǎo)致X方向正應(yīng)力的生成項大于0，X方向正應(yīng)力開始增長。但Y方向的雷諾正應(yīng)力生成項為0，這將使得Y方向正應(yīng)力被耗散進而導(dǎo)致切應(yīng)力生成項消失。再分配項則可以將X方向分量的一部分能量轉(zhuǎn)移給Y方向，使其保持一定強度，從而保持了切應(yīng)力的生成和強度，形成了一個維持湍流場的閉環(huán)機制。由于再分配項的重要作用，同時與生成項有相同的量級，且缺乏模化依據(jù)，因此是預(yù)測誤差的主要來源。Stress-BSL 模型將其表示為快速和慢速項之和。快速項中包含平均速度梯度，這說明雷諾應(yīng)力再分配不是局部過程，不能僅取決于當(dāng)?shù)匾稽c的參數(shù)。但模型構(gòu)造卻基于同一物理位置的湍流波動量的相關(guān)性完成模型閉合，這會導(dǎo)致一部分湍流信息的丟失，使得快速項的建模失準(zhǔn)，在速度梯度越大的區(qū)域，誤差也會越大。圖9 的結(jié)果證明，在邊界層內(nèi)剪切最強烈（速度梯度大）的區(qū)域，雷諾應(yīng)力計算值與實驗值誤差最大。其中，圖中縱坐標(biāo)Y?Y0、Z?Z0分別為Hump 模型和Bump 模型的當(dāng)?shù)乇诿娣ㄏ蚓嚯x，并使用弦長進行無量綱化。由于邊界層外側(cè)的速度梯度小于內(nèi)側(cè)，因此邊界層內(nèi)側(cè)的雷諾應(yīng)力預(yù)測誤差遠(yuǎn)大于外側(cè)。另外，當(dāng)前快速項的?；腔诳焖倩兘?，非線性模型可認(rèn)為是利用張量函數(shù)性質(zhì)在原線性模型中添加了各項異性張量的二次項［25］，而不是在物理現(xiàn)象層面的對流動信息的補充。

3.2 SST 模型

相較于雷諾應(yīng)力模型，渦黏性模型對雷諾應(yīng)力的求解是基于Boussinesq 假設(shè)，通過渦黏性系數(shù)（μt）完成的，因此SST 與Stress-BSL 模型對分離流動預(yù)測失準(zhǔn)的原因大相徑庭。SST 湍流模型的具體構(gòu)造細(xì)節(jié)可參見文獻［26］。

1）流動分離

起初的兩方程模型（k-ε模型、k-ω模型）對于逆壓梯度的預(yù)測能力很大程度上局限于邊界層對數(shù)區(qū)域［27-29］，預(yù)測效果并不好。Menter 和Egorov［14］發(fā)現(xiàn)邊界層尾跡區(qū)域的渦黏性水平最終決定了渦黏性模型對逆壓梯度流動的預(yù)測能力，因此在構(gòu)造SST 模型時引入了 Bradshaw 假設(shè)，在渦黏性系數(shù)表達式中強制保證湍流剪切應(yīng)力與邊界層尾跡區(qū)的湍動能成正比，從而獲得了很好的附著流和弱逆壓梯度流動的預(yù)測結(jié)果。其中渦黏性系數(shù)表達式為

對于存在很大速度梯度剪切很強的邊界層內(nèi)，渦黏性系數(shù)的表達式中分母則由平均應(yīng)變率（S）代替，使得湍動能生成項保持與耗散項相等，這也對應(yīng)了Bradshaw 假設(shè)中包含的湍流能量平衡。但對比雷諾應(yīng)力計算結(jié)果與實驗值（見圖10），顯然模型對于雷諾應(yīng)力的限制過度了，導(dǎo)致雷諾應(yīng)力預(yù)測值較實際情況偏小。這是由于Bradshaw 假設(shè)是通過Klebanoff 零壓力梯度平板邊界層試驗數(shù)據(jù)總結(jié)得出的［30］，其基礎(chǔ)為無壓力梯度的附著流動，當(dāng)逆壓梯度非常大時，流動的平均應(yīng)變率很高，使得湍動能生成項明顯大于耗散項，但 SST 模型中的控制器限制其生成項仍等于耗散項，使得湍動能和渦黏性系數(shù)小于真實值，造成了預(yù)測的雷諾應(yīng)力偏小。

圖10 SST 模型雷諾應(yīng)力計算結(jié)果與實驗值對比Fig.10 Comparison of SST and experimental Reyn‐olds stress

還需要注意的是，基于二維平板邊界層數(shù)據(jù)提出的Bradshaw 假設(shè)在推廣到三維邊界層時引入了“主雷諾應(yīng)力”的概念，將三維邊界層內(nèi)雷諾應(yīng)力張量其他各向分量的變化描述為對主雷諾應(yīng)力的變化的貢獻。但在強逆壓梯度下，雷諾切應(yīng)力各方向分量（UV、VW、UW Stress）間的差值很小，雷諾正應(yīng)力（UU、VV、WW Stress）與切應(yīng)力達到同一量級（見圖11）。顯然主雷諾應(yīng)力的的概念不再成立。

圖11 雷諾應(yīng)力對比Fig.11 Reynolds stress comparison

考慮到BSL 模型與SST 模型的構(gòu)造基本相同，唯一區(qū)別在于BSL 模型未對湍動能生成項進行限制，因此可以通過BSL 模型，對基于Brad‐shaw 假設(shè)的控制器所帶來的影響進行對比分析。通過BSL 模型與SST 模型的計算結(jié)果對比（見圖8、圖12），發(fā)現(xiàn) BSL 模型的計算結(jié)果不同程度上抑制了流動分離的提前，證實了上述分析的SST 模型預(yù)測失準(zhǔn)的原因。

圖12 BSL 與SST 模型計算結(jié)果對比Fig.12 Comparison of BSL and SST models

究其原因，是因為渦黏性模型所基于的Boussinesq 假設(shè)包含的線性關(guān)系、各向同性等“先天缺陷”，導(dǎo)致了渦黏性模型在預(yù)測存在突然的平均應(yīng)變率變化、額外應(yīng)變率和明顯的各向異性等特征的流動時出現(xiàn)失準(zhǔn)。

2）流動再附

如上所述，SST 模型對邊界層雷諾應(yīng)力的低估導(dǎo)致了流動分離的提前。相應(yīng)的，湍流模型對分離區(qū)上緣與主流之間的雷諾應(yīng)力的預(yù)測則影響了流動的再附位置的計算結(jié)果［31］。由于雷諾應(yīng)力表征脈動動量的輸運，因此在分離區(qū)上緣預(yù)測的雷諾應(yīng)力越大，動量輸運越多，湍流能量越高，流動再附位置越靠前。相反地，預(yù)測的雷諾應(yīng)力越小，則流動再附位置的計算結(jié)果越靠后?？紤]到SST 模型中雷諾切應(yīng)力主要由渦黏性系數(shù)（μt）和平均應(yīng)變率（S）決定。因此希望通過對比 Stress-BSL 模型與 SST 模型的計算結(jié)果找出主要影響因素（見圖13、圖14）。

圖13 渦黏性系數(shù)對比Fig.13 Comparison of eddy viscosity coefficient

圖14 平均應(yīng)變率對比Fig.14 Comparison of mean strain rate

發(fā)現(xiàn)相較于Stress-BSL 模型，SST 模型預(yù)測的μt偏小，S則偏大，但圖9、圖10 的結(jié)果表明SST 模型預(yù)測的雷諾應(yīng)力更小。據(jù)此推測，相較于S，μt的計算結(jié)果對于雷諾應(yīng)力的預(yù)測影響更大?？疾煊嬎悝蘴的SST 模型輸運方程，由于雷諾應(yīng)力被低估，導(dǎo)致湍動能輸運方程的生成項被低估，擴散項同樣偏小，上述共同作用使得湍動能k偏小，但表達式為隱式，影響較為復(fù)雜。對于比耗散率ω，利用量綱分析，可以得到湍流特征長度尺度δl的表達式

發(fā)現(xiàn)δl與ω之間為倒數(shù)關(guān)系。由于雷諾應(yīng)力表征湍流內(nèi)的動量輸運，因此雷諾應(yīng)力的預(yù)測與湍流能量緊密相關(guān)，而湍流能量又與湍流長度尺度直接相關(guān)。因此雷諾應(yīng)力被低估，說明湍流能量和長度尺度被低估，則ω被高估。考慮ω輸運方程，則耗散項同樣會被高估，這就要求生成項偏高更多，由于已經(jīng)分析過湍動能生成項被低估，因此可以認(rèn)為幾乎完全由位于分母位置的μt導(dǎo)致了ω生成項被高估。這說明，雖然k、ω輸運方程相互耦合，且均與μt相關(guān)，但顯然μt、ω的相互影響更為直接。據(jù)此推測，與ω相關(guān)的模型長度尺度的建模失準(zhǔn)、忽略耗散過程的各向異性特征是造成μt預(yù)測結(jié)果偏小的主要原因。

3）關(guān)鍵參數(shù)再標(biāo)定

基于Bradshaw 假設(shè)的限制器強制湍動能生成項與耗散項相等，使得SST 模型低估了雷諾應(yīng)力，但該假設(shè)豐富了模型的物理信息，因此并不能直接摒棄該假設(shè)。Zhao 等［32］研究后發(fā)現(xiàn)，系數(shù)a1對于分離激波具有很高的敏感度，并且a1是激波誘導(dǎo)分離的分離區(qū)內(nèi)壓力系數(shù)不確定性的主要來源，因此可以通過增大a1放寬對雷諾應(yīng)力的限制。經(jīng)過試算，發(fā)現(xiàn)a1=0.35 可以取得相對滿意的結(jié)果。

由圖15 可以發(fā)現(xiàn)，對a1進行調(diào)整可以在一定程度上取得更加接近實驗數(shù)據(jù)的結(jié)果，再標(biāo)定后的模型對于M6 機翼的背風(fēng)面激波預(yù)測更加合理，兩道激波合并為一道的位置向翼尖移動，且合成后的強激波的位置向下游移動。由壁面流線可以發(fā)現(xiàn)再標(biāo)定后模型預(yù)測的分離區(qū)更小，沒有出現(xiàn)分離起始位置過于靠近翼根位置的現(xiàn)象，整體對激波位置的捕捉和激波后區(qū)域的壓力系數(shù)的預(yù)測均更為準(zhǔn)確。

圖15 a1標(biāo)定前后壓力對比Fig.15 Comparison of pressure with different a1

4 結(jié)論

經(jīng)過試算，對比篩選后選用了SST 模型、Stress-BSL 模型作為渦黏性模型和雷諾應(yīng)力模型的代表對湍流分離流動進行了計算。選取了NACA4412 翼型、二維駝峰、二維跨聲速凸塊、跨聲速三維ONERA M6 機翼，從數(shù)學(xué)構(gòu)造和物理現(xiàn)象等方面入手對這2 個湍流模型的分離流動的預(yù)測性能進行了對比，給出了較為深入的機制分析，并基于分析所得原因?qū)δＰ瓦M行了常數(shù)再標(biāo)定，提出了改進方向，得到如下主要結(jié)論：

1）對于強逆壓梯度下的分離流動，Stress-BSL 模型和SST 模型均出現(xiàn)預(yù)測失準(zhǔn)，但Stress-BSL 模型的預(yù)測結(jié)果與實驗值誤差更小、性能更好。在激波誘導(dǎo)邊界層分離的復(fù)雜流動中，Stress-BSL 模型預(yù)測的激波特征和壓力系數(shù)均與實驗值吻合較好，而SST 模型則明顯預(yù)測分離區(qū)過大。在三維邊界層分離中，物理基礎(chǔ)更加可靠的Stress-BSL 模型的計算結(jié)果與實驗值誤差很小，而SST 模型計算的壓力系數(shù)與雷諾應(yīng)力則與實驗值偏差較大。

2）Stress-BSL 模型對分離流動預(yù)測的結(jié)果存在誤差，主要是由于雷諾應(yīng)力再分配項的模化不準(zhǔn)確，其基于快速畸變近似和張量函數(shù)性質(zhì)，缺乏與實際物理過程的聯(lián)系。雷諾應(yīng)力輸運方程中耗散項同樣復(fù)雜，影響重要但程度不如再分配項。其他如擴散項等對雷諾應(yīng)力計算的影響不大。另外，計算中發(fā)現(xiàn)，數(shù)值穩(wěn)定性仍是雷諾應(yīng)力模型亟待解決的問題，對條件要求較為苛刻，需根據(jù)其特點開發(fā)相應(yīng)的數(shù)值方法。

3）SST 模型對分離流動預(yù)測的結(jié)果取決于對雷諾應(yīng)力的?；?。分離提前是因為基于Brad‐shaw 假設(shè)的控制器，使得本應(yīng)遠(yuǎn)大于耗散項的湍動能生成項被限制，導(dǎo)致雷諾應(yīng)力預(yù)測偏低。而耗散尺度建模失準(zhǔn)導(dǎo)致的主導(dǎo)雷諾應(yīng)力預(yù)測結(jié)果的渦黏性系數(shù)被低估則是再附滯后的主要原因。

4）對于流動分離提前的問題，本文通過對湍動能生成項限制放寬后，預(yù)測結(jié)果得到改善。另外，渦黏性模型對分離流動計算不準(zhǔn)的問題很大程度上可以通過對湍流尺度的建模進行改善。對于湍流各向異性特性的補充，則應(yīng)聚焦于對雷諾應(yīng)力本構(gòu)關(guān)系添加非線性項進行修改。

5）模型的長度尺度對預(yù)測準(zhǔn)確性起到關(guān)鍵作用，可以在消除邊界奇性的基礎(chǔ)上，補充對長度尺度的控制，以獲得更好的預(yù)測結(jié)果。