基于數(shù)值模擬的流場(chǎng)附面層邊緣識(shí)別方法

王 磊 張堃元 司江濤 劉凱禮 蔡北京 楊心宇

(1.上海飛機(jī)設(shè)計(jì)研究院，上海 201210；2.南京航空航天大學(xué)，南京 210016)

0 引言

附面層是壁面附近流動(dòng)因黏性而受阻滯的一部分區(qū)域，其邊緣通常取在速度達(dá)到主流0.99倍的位置。對(duì)附面層的研究是進(jìn)氣道、噴管、機(jī)翼等內(nèi)外流流場(chǎng)設(shè)計(jì)中的重要內(nèi)容，例如對(duì)無黏設(shè)計(jì)結(jié)果進(jìn)行黏性修正、分析附面層對(duì)主流的影響、判斷流動(dòng)是否有分離的趨勢(shì)等。獲取附面層厚度以及附面層內(nèi)參數(shù)分布是進(jìn)行這些研究的基礎(chǔ)。

工程設(shè)計(jì)中有一些簡(jiǎn)單的經(jīng)驗(yàn)公式可用于估算附面層的發(fā)展。例如，文獻(xiàn)[1-4]計(jì)算噴管附面層位移厚度δ1采用發(fā)展長度x的線性函數(shù)：

δ1=xtanβ

(1)

其中，角度β與噴管出口馬赫數(shù)有關(guān)；文獻(xiàn)[5-6]根據(jù)公式(2)計(jì)算高超聲速進(jìn)氣道附面層位移厚度：

δ1=Ax+Bxe-x

(2)

其中，A、B為根據(jù)經(jīng)驗(yàn)確定的系數(shù)；文獻(xiàn)[7]中應(yīng)用了一種形式更復(fù)雜的經(jīng)驗(yàn)公式。這些公式計(jì)算量小、使用方便，但顯然其計(jì)算精度以及能夠獲取的信息是比較有限的。

求解附面層方程是較為傳統(tǒng)的附面層近似計(jì)算方法。其中多數(shù)是對(duì)附面層動(dòng)量積分方程補(bǔ)充經(jīng)驗(yàn)關(guān)系式后求解，例如文獻(xiàn)[8-11]在進(jìn)氣道、噴管、葉柵等流道設(shè)計(jì)的黏性修正中應(yīng)用的方法。也有些直接針對(duì)附面層偏微分方程組進(jìn)行數(shù)值計(jì)算，例如Cebeci[12]、杜國樑[13]所提出的方法，其中Cebeci提供的計(jì)算程序在許多進(jìn)氣道、噴管設(shè)計(jì)中得到應(yīng)用[14-18]。

隨著CFD技術(shù)的飛速發(fā)展，基于黏性N-S方程組的數(shù)值模擬已得到普遍應(yīng)用，與上述兩種方法相比，其精度更高，適應(yīng)性更廣。但N-S方程組求解中一般不特意區(qū)分附面層與主流，這使得如何在流場(chǎng)數(shù)據(jù)中準(zhǔn)確認(rèn)定附面層邊緣成了新的問題。按照通常的定義，附面層邊緣取在速度達(dá)到主流0.99倍的位置，而在較復(fù)雜的流場(chǎng)中主流本身并不均勻，因此需要根據(jù)經(jīng)驗(yàn)人為給定主流速度的數(shù)值[19-21]，有較大的主觀性并且難以實(shí)現(xiàn)自動(dòng)化識(shí)別；有文獻(xiàn)中根據(jù)附面層和主流中法向速度梯度大小的差異判斷附面層邊緣[22]，也存在類似的問題。

針對(duì)在流場(chǎng)數(shù)值計(jì)算結(jié)果中識(shí)別附面層的困難之處，本文根據(jù)附面層的概念，提出以排除了附面層干擾的“參考主流流場(chǎng)”中的速度代替實(shí)際流動(dòng)主流速度確定附面層邊緣的方法，以典型的超聲速壓縮型面為例對(duì)該識(shí)別方法進(jìn)行了驗(yàn)證，并根據(jù)識(shí)別結(jié)果考察了應(yīng)用附面層修正的設(shè)計(jì)效果。

1 附面層邊緣識(shí)別方法

針對(duì)在流場(chǎng)數(shù)值計(jì)算結(jié)果中識(shí)別附面層的困難之處，本文根據(jù)附面層的概念，提出以排除了附面層干擾的“參考主流流場(chǎng)”中的速度代替實(shí)際流動(dòng)主流速度確定附面層邊緣的方法，并以超聲速壓縮型面為例對(duì)該識(shí)別方法進(jìn)行了驗(yàn)證。

如引言所述，根據(jù)數(shù)值模擬結(jié)果識(shí)別附面層困難之處在于復(fù)雜流場(chǎng)中主流速度的確定，主流應(yīng)是附面層外、沒有受到壁面與流體之間剪切力影響的流動(dòng)。滿足這種要求的流場(chǎng)通過零剪切滑移壁面條件下的數(shù)值模擬可以獲得，這里稱之為“參考主流”，其中流動(dòng)沒有附面層，完全排除了壁面與流體之間剪切力的影響。

在流場(chǎng)中黏性相互作用較弱時(shí)(對(duì)其強(qiáng)弱的判斷參見文獻(xiàn)[23]的說明)，可以認(rèn)為實(shí)際流動(dòng)中主流流動(dòng)與這種參考主流相同。因此根據(jù)附面層的定義判斷其邊緣時(shí)，可以用這種參考主流中對(duì)應(yīng)位置的速度代替所需的實(shí)際流場(chǎng)主流速度，這種方法避免了人為指定主流速度的主觀性，因而更為科學(xué)合理，也有利于識(shí)別過程的自動(dòng)化。

圖1為該方法識(shí)別過程的示意圖。在壁面上的指定位置(例如圖1中x)建立沿壁面法向的單位向量n和沿壁面切向的單位向量s，沿向量n方向新建一系列網(wǎng)格點(diǎn)，依次驗(yàn)證各點(diǎn)參數(shù)是否滿足公式(3)：

圖1 附面層邊緣識(shí)別過程示意圖

(V-0.99Vref)·s≥0

(3)

式中，Vref為上述參考主流中該網(wǎng)格點(diǎn)的速度向量，V為實(shí)際流動(dòng)中的速度向量(可以是數(shù)值計(jì)算或?qū)嶒?yàn)測(cè)量的結(jié)果)，第一個(gè)滿足該式的點(diǎn)即判定為該截面附面層邊緣。對(duì)于壁面所有位置重復(fù)該過程，能夠獲取一系列邊緣點(diǎn)，依次連接就得到了整個(gè)壁面上的附面層邊緣曲線。

值得說明的是，對(duì)于主流湍流度較小并且不存在強(qiáng)剪切的流動(dòng)，也可以用無黏條件下的模擬結(jié)果近似代替滑移壁面條件下的計(jì)算作為參考主流，這有利于減小計(jì)算量。

本文在商業(yè)軟件ANSYS FLUENT的Scheme語言環(huán)境中編制了能夠自動(dòng)完成整個(gè)計(jì)算和識(shí)別過程的程序。

2 流場(chǎng)數(shù)值模擬方法及算例校驗(yàn)

本文研究中采用軟件ANSYS FLUENT進(jìn)行數(shù)值模擬，首先利用文獻(xiàn)[24-25]中超聲速壓縮拐角的實(shí)驗(yàn)結(jié)果考核了所采用的軟件和算法對(duì)超聲速流動(dòng)中附面層計(jì)算的可靠性。

圖2為該實(shí)驗(yàn)?zāi)Ｐ偷氖疽鈭D。風(fēng)洞實(shí)驗(yàn)段尺寸203 mm×203 mm，來流馬赫數(shù)2.87，總壓0.68 MPa，總溫282 K，實(shí)驗(yàn)?zāi)Ｐ蜑?°的斜楔，模型寬度152 mm，貼于風(fēng)洞底面安裝，斜楔拐角距離噴管出口1.2 m。實(shí)驗(yàn)采集了壓縮拐角前后若干截面上的參數(shù)分布，本文采用了拐角前25.4 mm(記為截面1)和拐角后66.04 mm(記為截面2)兩個(gè)垂直于壁面的截面上時(shí)均速度分布。實(shí)驗(yàn)的其他數(shù)據(jù)見文獻(xiàn)[25]。

圖2 實(shí)驗(yàn)?zāi)Ｐ褪疽鈭D[24]

針對(duì)該實(shí)驗(yàn)進(jìn)行了二維數(shù)值計(jì)算。計(jì)算模型參考了文獻(xiàn)[26-27]的做法，在截面1前補(bǔ)充長度為2.0 m的附面層發(fā)展平板以保證截面1附面層速度分布與實(shí)驗(yàn)相同，圖3為計(jì)算模型示意圖。計(jì)算的總網(wǎng)格量約400 000，壁面y+約25。進(jìn)口邊界條件為壓力遠(yuǎn)場(chǎng)，出口為壓力出口，壁面為絕熱無滑移固壁邊界。

圖3 計(jì)算模型示意圖[24]

采用基于密度的耦合隱式算法求解RANS方程，對(duì)流通量采用基于二階迎風(fēng)格式的AUSM方法計(jì)算。湍流模型選用標(biāo)準(zhǔn)k-ε模型，采用二階迎風(fēng)格式離散，近壁采用標(biāo)準(zhǔn)壁面函數(shù)處理。比熱比取定值1.4，黏性系數(shù)采用Sutherland公式計(jì)算。以來流條件初始化開始計(jì)算，當(dāng)各項(xiàng)殘差不再下降并且進(jìn)出口流量相對(duì)誤差低于10-3時(shí)完成計(jì)算。

圖4顯示了截面1和截面2速度分布的實(shí)驗(yàn)與計(jì)算結(jié)果，圖中η為沿壁面法向的坐標(biāo)，uξ為平行于壁面方向的速度大小。其中截面1計(jì)算與實(shí)驗(yàn)最大誤差為5.5%，均方根誤差0.35%，表明拐角前來流滿足實(shí)驗(yàn)條件；截面2計(jì)算結(jié)果與實(shí)驗(yàn)相吻合，速度最大誤差為4.4%，均方根誤差0.27%，表明所采用的計(jì)算方法是可靠的，后續(xù)研究中均采用類似的算法進(jìn)行計(jì)算。

圖4 兩個(gè)截面上速度分布的實(shí)驗(yàn)和計(jì)算結(jié)果

3 附面層識(shí)別方法的驗(yàn)證

取來流條件為：M=6，p=2 511 Pa，T=221.6 K，在無黏條件下設(shè)計(jì)了以下兩種超聲速壓縮型面：

1)三級(jí)斜楔壓縮型面：三個(gè)壓縮角分別為5.07°、5.71°和6.47°，無黏情況下三道激波匯聚點(diǎn)距離型面前緣高度為100 mm，壁面長度為433 mm；

2)彎曲壓縮型面：根據(jù)壓力分布反設(shè)計(jì)得到，設(shè)計(jì)方法參見文獻(xiàn)[28]，無黏情況下激波高度100 mm，壁面長度L=394.5 mm，壁面壓力分布為分段函數(shù)：

(4)

其中,下標(biāo)i為分段節(jié)點(diǎn)序號(hào)，p1/p=1.67，式中其他參數(shù)取值見表1。

表1 壁面壓力分布函數(shù)中的參數(shù)取值

對(duì)這兩個(gè)壓縮型面進(jìn)行數(shù)值模擬得到流場(chǎng)參數(shù)，計(jì)算中來流湍流度為1%，來流湍流黏性與分子黏性之比為1。之后對(duì)兩個(gè)壓縮型面進(jìn)行無黏計(jì)算作為參考主流，應(yīng)用第1節(jié)所述識(shí)別算法提取了附面層邊緣。

圖5顯示了兩個(gè)壓縮流場(chǎng)的馬赫數(shù)等值線以及提取出的附面層邊緣曲線。從圖中可見，除有激波穿過的局部范圍外，所提取的邊緣曲線光滑連續(xù)，均位于速度急劇變化區(qū)域的外緣，與附面層的物理意義相符。

(a)三級(jí)斜楔壓縮

圖6為兩個(gè)算例出口截面近壁區(qū)域的速度和馬赫數(shù)分布曲線。其中三級(jí)斜楔壓縮出口截面主流參數(shù)均勻，按照0.99倍的定義，附面層厚度應(yīng)為6.76 mm；而本算法識(shí)別的結(jié)果為7.04 mm(參見圖6(a))，兩者相對(duì)誤差僅4.1%，表明識(shí)別的結(jié)果是比較準(zhǔn)確的。

從圖6(b)可見，彎曲壓縮型面出口截面主流不均勻，馬赫數(shù)和速度分布在遠(yuǎn)離壁面時(shí)均有明顯增加。按照傳統(tǒng)的做法，需要人為選定主流速度的數(shù)值，不僅存在較大的主觀性，也難以實(shí)現(xiàn)自動(dòng)化的識(shí)別。而本算法所識(shí)別的附面層邊緣位于速度梯度有明顯變化的位置，符合附面層的概念。

(a)三級(jí)斜楔壓縮

4 應(yīng)用附面層修正的型面設(shè)計(jì)

在工程應(yīng)用中設(shè)計(jì)進(jìn)氣道、噴管等流道型面時(shí)，常采用先進(jìn)行無黏設(shè)計(jì)、再修正附面層影響的方法。在進(jìn)行修正設(shè)計(jì)時(shí)，通過本文所述的附面層邊緣識(shí)別方法，能夠更準(zhǔn)確、方便地獲得各截面附面層邊緣位置及附面層內(nèi)參數(shù)分布，之后容易根據(jù)定義計(jì)算其位移厚度δ1：

(5)

其中，δ為附面層厚度，ρ、uξ分別為密度和沿壁面切向的速度分量，下標(biāo)δ表示附面層邊緣位置的參數(shù)。對(duì)于無黏設(shè)計(jì)的壁面坐標(biāo)向外移動(dòng)該距離可修正附面層對(duì)流量的影響。

應(yīng)用此方法對(duì)上述彎曲壓縮型面進(jìn)行了修正。圖7為修正前后壁面形狀、激波形狀的對(duì)比，圖8為修正前后壁面壓力分布及與設(shè)計(jì)值的對(duì)比?？梢娦拚蠹げǜ叨扰c特征線法設(shè)計(jì)結(jié)果的相對(duì)誤差大約從2.0%降低到0.3%，壁面壓力與給定壓力分布之間的誤差也明顯降低，末端壓力的相對(duì)誤差大小|ε|從6.6%降低至2.3%。

圖7 修正前后的壁面形狀和激波形狀

圖8 修正前后壁面壓力分布及與設(shè)計(jì)值間的誤差

5 結(jié)論

1)研究了依據(jù)數(shù)值模擬結(jié)果識(shí)別附面層的方法，與依靠經(jīng)驗(yàn)公式、附面層方程的近似計(jì)算相比有更高的精度和更廣的適應(yīng)性；

2)提出了用“參考主流”代替實(shí)際主流識(shí)別流場(chǎng)附面層邊緣的方法：以零剪切力滑移壁面邊界條件下數(shù)值模擬得到的不受附面層干擾的流場(chǎng)作為參考主流，在根據(jù)附面層定義確定其邊緣時(shí)，以該參考主流流場(chǎng)中的速度代替實(shí)際的主流速度進(jìn)行判斷。這種方法避免了傳統(tǒng)做法中人為指定主流速度的主觀性，因而更為科學(xué)合理，也有利于識(shí)別過程的自動(dòng)化；

3)對(duì)兩個(gè)超聲速壓縮流場(chǎng)附面層識(shí)別的測(cè)試表明，該識(shí)別方法獲得的附面層邊緣比較準(zhǔn)確，其中斜楔壓縮出口截面上附面層厚度與采用實(shí)際主流速度判斷得到的厚度相對(duì)誤差僅4.1%；

4)根據(jù)本算法識(shí)別的附面層參數(shù)計(jì)算位移厚度對(duì)彎曲壓縮型面進(jìn)行黏性修正后，激波高度與無黏設(shè)計(jì)結(jié)果之間的相對(duì)誤差從修正前的2.0%降低至0.3%，壁面末端壓力相對(duì)誤差從6.6%降低至2.3%，表明該方法可取得良好的效果。