大迎角跨音速三角翼漩渦特性及激波與渦干擾的數(shù)值模擬研究

李自啟，梁斌，羅松，劉敏

（中航工業(yè)洪都，江西南昌330024）

大迎角跨音速三角翼漩渦特性及激波與渦干擾的數(shù)值模擬研究

李自啟，梁斌，羅松，劉敏

（中航工業(yè)洪都，江西南昌330024）

采用對流迎風(fēng)分裂格式改進(jìn)形式(AUSM+)和SST兩方程湍流模型結(jié)合求解三維雷諾平均Navier Stokes(RANS)方程。通過對跨音速65°后掠尖前緣三角翼的數(shù)值模擬，研究了三角翼上翼面正激波的移動情況以及渦核破裂點變化情況。計算結(jié)果與實驗對比表明：采用AUSM+格式和SST兩方程湍流模型結(jié)合能夠準(zhǔn)確模擬出激波結(jié)構(gòu)和大迎角下空間渦結(jié)構(gòu)以及激波與渦相互干擾情況。

迎風(fēng)分裂格式；SST湍流模型；三角翼；激波；空間渦

0 引言

隨著計算流體技術(shù)的發(fā)展，國內(nèi)外工程師對亞音速階段漩渦流動的計算研究已趨于成熟。然而，當(dāng)馬赫數(shù)增加到跨音速階段，激波出現(xiàn)以后，激波與附面層的相互作用以及與空間渦的相互干擾使流場變得極其復(fù)雜，呈現(xiàn)出與亞音速階段不同的流動現(xiàn)象。

由于跨音速、大攻角條件下三角翼上翼面漩渦的流動特別復(fù)雜，而當(dāng)前對漩渦的形態(tài)及破裂的變化過程的模擬還不夠成熟。關(guān)于翼面漩渦的初始形態(tài)以及初始分離位置,Sarpkaya[1]認(rèn)為,主要為氣泡型和螺旋型兩種形態(tài)。關(guān)于旋渦破裂的發(fā)展演變過程,前人在實驗和數(shù)值模擬兩方面做了大量工作。風(fēng)洞實驗顯示,旋渦破裂的型態(tài)是由泡型到螺旋型,最后又回到泡型。在跨音速流場中，隨著迎角的增加，三角翼上出現(xiàn)復(fù)雜的激波[7-10]，Elsenaar[11]和Hoeijmakers[12]將其歸納為2個激波系統(tǒng)：第一個是由于前緣渦旋轉(zhuǎn)而產(chǎn)生的，位于渦的底部，第二個垂直于對稱面，位于機翼的后緣，但是在靠近前緣渦的區(qū)域，受前緣渦發(fā)生的影響發(fā)生了彎曲。

在三角翼旋渦破裂的數(shù)值模擬方面,許多人采用通量差分分裂迎風(fēng)格式，迎風(fēng)耗散格式[2]主要分為矢通量分裂(FVS)和基于Riemann解的通量差分分裂(FDS)。FVS格式有較大的數(shù)值耗散、穩(wěn)定性好，但是激波分辨率不高,求解粘性流動的精度較低。FDS格式的代表格式是Roe-FDS格式,該格式具有較高的間斷分辨率和數(shù)值模擬精度,但是存在非物理的數(shù)值解,需要加入人工熵修正。在20世紀(jì)90年代,Liou和Stefen提出了AUSM[3](Advection Upstream Splitting Method)格式。AUSM格式的主要思想是將無粘通量項分裂為對流通量項和壓力項,基于馬赫數(shù)對其分別進(jìn)行特征分裂。AUSM格式在CFD研究領(lǐng)域得到重視,并且迅速發(fā)展和完善,相繼出現(xiàn)了AUSM+[4]。AUSM+一類格式具有較強的分辨激波、滑流等物理間斷的能力,而且具有良好的保正定性,計算高速流時不會產(chǎn)生“Carbuncle”現(xiàn)象。

本文采用AUSM+格式和三步Runge-Kutta時間離散方法以及SST二方程湍流模型[5,6]數(shù)值求解定常雷諾平均的RANS方程。對尖前緣三角翼跨音速擾流進(jìn)行了模擬，分析了跨音速條件下三角翼上表面渦結(jié)構(gòu)、激波分布以及前緣渦和激波的相互干擾。

1 控制方程及離散方法

1.1 控制方程

笛卡爾坐標(biāo)系積分形式的Navier-Stokes方程為：

其中：

式中：Q為守恒變量；E、F、G為無黏通量；Ev、Fv、Gv為黏性通量；S為源項。

1.2 離散方法

本文中無粘性的數(shù)值離散采用AUSM+格式[4]，AUSM格式的主要思想是將無粘通量項分裂為對流通量項和壓力項,基于馬赫數(shù)對其分別進(jìn)行特征分裂。在超聲速下,根據(jù)特征變量,格式具有正確的迎風(fēng)方向。在亞、跨聲速下,格式中需要包含所有的特征波。為了改進(jìn)AUSM格式的不足,如：激波后形成“紅寶石”現(xiàn)象,1994年Liou和Steffen[3]在AUSM格式的基礎(chǔ)上,對該格式進(jìn)行了改進(jìn),構(gòu)成了后來的AUSM+格式。通量離散形式如下：

AUSM+格式中采用的馬赫數(shù)分裂方法如下：

左右馬赫數(shù)定義如下：

式中，C1/2為界面聲速；cL為左側(cè)網(wǎng)格聲速；cR為右側(cè)網(wǎng)格聲速。

控制面處左右狀態(tài)決定了格式的精度,本文采用三階MUSCL方法插值得到計算時所用左右狀態(tài)，針對于控制方程（1）中的粘性項采用二階中心差分格式，控制方程（1）中的時間項采用三步Runge-Kutta時間離散方法。

2 湍流方程

湍流模型采用SST模型[5,6]，k-ω模型求解了兩個運輸方程，一個關(guān)于湍動能k，另一個關(guān)于頻率ω。k-ω模型的優(yōu)點是可以很好的處理近壁處低雷諾數(shù)的數(shù)值計算。此模型不涉及k-ε模型中復(fù)雜的非線性衰減函數(shù)，因此計算結(jié)果更加準(zhǔn)確，收斂性更好。低雷諾數(shù)k-ε要求近壁面處y+＜0.2。而低雷諾數(shù)的k-ω模型要求近壁面處y+＜2。

k方程為：

ω方程為：

式中：Pk為湍流生成速率；模型中的常數(shù)由以下給出，即：

基于SST的k-ω方程考慮了剪切應(yīng)力的傳輸，可以精確的預(yù)測流動的開始和負(fù)壓力梯度條件下流體的分離量。SST模型的最大優(yōu)點就在于考慮了湍流剪切力，從而不會對湍流粘度造成過度預(yù)測。其傳輸行為可由包含限制數(shù)的渦流粘度方程求得：

式中：F2是一個混合函數(shù)，其功能是用來約束壁面附面層流動的限制數(shù)。

3 幾何模型和計算網(wǎng)格[13]

幾何模型采用參考文獻(xiàn)[13]中第四章提到的標(biāo)模，機翼為65°切尖三角翼，實驗數(shù)模翼型為在NACA64A005基礎(chǔ)上將前緣用尖頭的雙圓弧修形。本文直接采用NACA64A005，平面參數(shù)無鴨翼，模型如圖1（a）所示，根弦長為c=1.2m，后掠角為65°，展長為0.452m，稍弦長為0.239m。

計算使用的網(wǎng)格采用C-O型網(wǎng)格，網(wǎng)格首層高度0.001mm，網(wǎng)格空間增長率1.2，網(wǎng)格總數(shù)421萬，網(wǎng)格拓?fù)浣Y(jié)構(gòu)如圖1（b）。

圖1 幾何模型與網(wǎng)格

4 計算結(jié)果及分析

為了更形象說明跨音速條件下三角翼空間渦的分布情況，圖2和圖3給出了三角翼上翼面典型空間渦的結(jié)構(gòu)和三角翼上翼面表面壓力分布情況。由圖2可以看出，本文方法清晰撲捉出主渦、次渦以及后激波的位置。文獻(xiàn)中采用三階迎風(fēng)格式，由圖3可以看出，該方法能更清晰的撲捉正激波的位置。圖4給出了截面處表面壓力系數(shù)分布，由圖3可知，在70%弦長以前，能夠準(zhǔn)確捕捉出橫向激波的位置，70%弦長以后計算結(jié)果在機翼后緣出現(xiàn)了激波，導(dǎo)致主次渦的破裂。為了研究跨音速三角翼空間渦的變化以及空間渦與激波的相互干擾情況，選取Ma= 0.85,迎角α=19°、20°、21°、23°作為計算狀態(tài)進(jìn)行數(shù)值模擬。

圖2 三角翼空間流動情況

圖3 計算結(jié)果與迎風(fēng)格式結(jié)果對比

4.1 跨音速條件下空間渦隨迎角的變化

圖4給出了三角翼不同迎角下空間渦變化情況，主渦和次渦起始位置靠近三角翼尖點處，主渦和次渦在20°以前比較穩(wěn)定，當(dāng)迎角達(dá)到21°時，在x軸向80%左右處主渦能量變?nèi)酰螠u逐漸破裂，并且隨著迎角的增加，主次渦破裂位置向前移動。

圖4不同迎角空間渦變化情況

圖5 給出了主渦渦核處空間流線，在21°迎角以前，渦核處的空間流線比較緊湊，此時空間渦能量強，當(dāng)迎角達(dá)到21°以后，渦核處的空間流線在根弦線80%以后出現(xiàn)發(fā)散，此時空間渦能量變?nèi)?，并且隨著迎角的增加，渦能量變?nèi)醯奈恢弥饾u向前移動，渦軸的位置逐漸向三角翼內(nèi)側(cè)移動，渦的起點位置逐漸向三角翼尖點移動。由三角翼表面壓力系數(shù)分布可以看出，隨著迎角的增加，后激波位置和渦核發(fā)散點逐漸前移，由此可知，后激波將導(dǎo)致渦能量降低甚至破裂。

圖6給出了不同截面處表面壓力系數(shù)分布，在根弦線50%以前，截面壓力系數(shù)分布顯示，在主渦處存在較強的吸力峰值，主渦到前緣之間存在較弱的二次吸力峰值。隨迎角增加，兩個吸力隨著渦的移動向三角翼內(nèi)側(cè)移動。

4.2 跨音速條件下空間激波隨迎角的變化

圖7給出了X=0.5m截面處空間馬赫分布，1處激波由于主渦與二次分離區(qū)之間的橫向流動，形成了一個類似收縮-擴張管的流動，局部加速到超音速，并且流動不能一直保持超音速，最終以激波減速。2處激波由于次渦順時針旋轉(zhuǎn)產(chǎn)生，迎角大于23°時此處出現(xiàn)分離。由于氣流是向上加速產(chǎn)生的，故3處激波出現(xiàn)在剪切層上方，靠近三角翼前緣位置。

圖5 不同迎角三角翼上空間渦線

圖6三角翼不同截面表面壓力系數(shù)分布

圖8 給出了不同迎角下三角翼上表面流線，由圖可以看出，氣流在三角翼前緣形成穩(wěn)定的漩渦，渦跡線如1處所示，漩渦的轉(zhuǎn)動在2處形成橫向激波，使氣流在此形成二次分離流，分離氣流又形成漩渦，漩渦在貼近上翼面處形成激波，使氣流在3處產(chǎn)生三次分離流。由于面積率和翼型的綜合作用,在機翼后緣4處形成正激波。又由于激波的作用，在機翼后緣形成氣流分流區(qū)，如5處。圖9給出了正激波隨馬赫數(shù)變化情況，隨著迎角的增加，激波位置逐漸前移。

圖7 x=0.5m處截面展向激波分布（1、2、3為激波）

圖8 不同迎角三角翼上表面流線

圖9 正激波分布

5 結(jié)論

1）用AUSM+格式結(jié)合SST二方程湍流模型求解RANS方程，模擬鈍前緣三角翼擾流跨音速特性。該方法準(zhǔn)確地模擬出三角翼上表面的主渦、次渦吸力峰值以及流場中復(fù)雜流動情況，并準(zhǔn)確地模擬出展向流動激波和正激波的位置。

2）迎角達(dá)到21°以后，三角翼后激波位置突然向前緣移動，并且與主渦發(fā)生干擾，使渦的破裂位置突然前移。

3）前緣分離渦的強度、后激波的位置以及渦軸處能量大小直接影響空間渦破裂。

4）渦能量減弱的原因是,流動空間中存在粘性，并且沿渦核移動的方向具有逆壓梯度，尤其是渦核處，逆壓梯度更大，會導(dǎo)致渦能量減低甚至破裂。

[1]Sarpkaya T.On Stationary and T raveling Vortex Breakdown[J].J.Fluid Mech.,1977,45(3)：545-559.

[2]Steger J L,Warming R L.Flux vector splitting of the inviscid gas dynamic equations with application tofinitedifferencemethods[J].Journalof Computational Physics,1981.40(2)：263-293.

[3]Liou M S.Ten years in the making AUSM family[R].AIAA Paper 2001225212CP,2001.

[4]Liou M S.A sequel to AUSM：AUSM+[J]. Journal of Computational Physics,1996,129(2)：364-382.

[5]Hellsten,A.New Two-Equation Turbulence Model for Aerodynamics Applications[D].Department of Mechanical Engineering,HelsinkiUniversity of Technology(Espoo,Finland),2004.

[6]Wallin,S.and Johansson,A.V.An Explicit andAlgebraicReynoldsStressModelfor incompressible and Compressible Turbulent Flows[J]. Journal of Fluid Mechanics,Vol.403,2000,pp.89-132.

[7]Drougge,G.The International Vortex Flow Experiment for Computer Code Validation[J].ICAS Proceedings,1988,Vol.1,pp.35-41.

[8]Elsenaar,A.Hjelmberg,L.Bütefisch,K.-A., and Bannink,W.J.,The International Vortex Flow Experiment[？].AGARDCP 437,Vol.1,1988,pp.9-1 to 9-23.

[9]Londenberg,W.K.,Transonic Navier-Stokes Calculations About a 65 Degree Delta Wing,NASACR-4635,Nov.1994.

[10]Chiba,K.and Obayashi,S.,CFD Visualization of Second Primary Vortex Structure on a 65-Degree Delta Wing,AIAA Paper 2004-1231,Jan.2004.

[11]Erickson G E,Rogers L W.Experimental study of the vortex flow behaviour on a generic fighter wing at subsonic and transonic speeds.AIAA-1987-1262.1987

[12]Donohoe S R,Bannink W J.Surface reflective visualisation of shock-wave/vortex interactions above a delta wing.AIAA Journal,1995,33(4)：680-687.

[13]方寶瑞.飛機氣動布局設(shè)計[M].北京：航空工業(yè)出版社,1997.

[14]朱自強.現(xiàn)代飛機空氣動力學(xué)[M].北京：北京航空航天大學(xué)出版社,1982.

>>>作者簡介

李自啟，男，1985年1月出生，2013年畢業(yè)于南京航空航天大學(xué)，現(xiàn)從事氣動布局設(shè)計工作。

Study on Vortex Features of Big-elevation Angle Transonic Delta Wing and Numerical Simulation of Interference of Shock Wave and Vortex

Li Ziqi,Liang Bin,Luo Song,Liu Min
(AVIC Hongdu Aviation Industry Group,Nanchang,Jiangxi,330024)

The combination of AUSM+(Advection Upstream Splitting Method)and SST two-equation turbulence model are used to solve three-dimension Navier Stokes(RANS)equation.With the numerical simulation of transonic sweptback leading-edge delta wing of 65°,study on movement of normal shock wave on the wing upper surface of delta wing and change of vortex fracture point has been done.The calculation result and test comparison show that the combination of AUSM+and SST two-equation turbulence model can accurately simulate structures of shock wave and space vortex at a big elevation angle,as well as interference of shock wave and vortex each other.

upstream splitting method;SST turbulence model;delta wing;shock wave;space vortex

2016-01-13）