湍流入流條件在強(qiáng)迫轉(zhuǎn)捩研究中的應(yīng)用

逯學(xué)志，黃章峰，2，羅紀(jì)生

（1．天津大學(xué)力學(xué)系，天津300072；2．空氣動(dòng)力學(xué)國(guó)家重點(diǎn)實(shí)驗(yàn)室，四川 綿陽(yáng)621000）

0 引 言

直接數(shù)值模擬（Direct Numerical Simulation，DNS）可提供流場(chǎng)的高分辨率的全部信息，是研究湍流機(jī)制的主要工具之一，尤其是對(duì)超聲速和高超聲速邊界層流動(dòng)。然而采用空間模式的直接數(shù)值模擬（Spatial Direct Numerical Simulation，SDNS），其入流條件應(yīng)能反映上游的特性，對(duì)于湍流計(jì)算而言入流條件如何給是一大難題。本文介紹了一種新的提法，即將由時(shí)間模式直接數(shù)值模擬（Temporal Direct Numerical Simulation，TDNS）得到的充分發(fā)展湍流流場(chǎng)（只需要一個(gè)時(shí)刻的流場(chǎng)）以適當(dāng)方式轉(zhuǎn)換成SDNS的入口流場(chǎng)。計(jì)算結(jié)果證實(shí)該方法可行且有效。

對(duì)于高超聲速飛行器，下表面往往有兩到三個(gè)壓縮折角，飛行器下部的氣流經(jīng)過(guò)這些壓縮折角后進(jìn)入發(fā)動(dòng)機(jī)。如果折角處流動(dòng)狀態(tài)為層流狀態(tài)，往往會(huì)在折角處發(fā)生流動(dòng)分離，從而導(dǎo)致氣流偏離下表面，只有較少部分的氣體流入了發(fā)動(dòng)機(jī)，嚴(yán)重影響了發(fā)動(dòng)機(jī)的性能。然而在邊界層的前部，流動(dòng)一般是層流，而且自然狀態(tài)下層流段可以很長(zhǎng)。這時(shí)為使轉(zhuǎn)捩提前，需要采用加人工粗糙度的方法，使之在第一個(gè)壓縮折角前就發(fā)生強(qiáng)迫轉(zhuǎn)捩。

美國(guó)主要采用地面試驗(yàn)和飛行試驗(yàn)的驗(yàn)證方法來(lái)研究各種轉(zhuǎn)捩裝置的效果［1－3］，并通過(guò)傳熱系數(shù)是否快速上升來(lái)判斷人工粗糙度促使轉(zhuǎn)捩是否成功。但實(shí)際上，實(shí)驗(yàn)和數(shù)值模擬結(jié)果表明，由于壁面粗糙引起的流向渦也會(huì)導(dǎo)致St上升［4－5］，如圖1所示，其趨勢(shì)和由于轉(zhuǎn)捩導(dǎo)致的很相似。其中St的定義為：因此傳熱系數(shù)上升既可能是由于轉(zhuǎn)捩導(dǎo)致的，也可能是由于較強(qiáng)的流向渦引起的，不是強(qiáng)迫轉(zhuǎn)捩后達(dá)到湍流狀態(tài)的可靠判據(jù)。有必要對(duì)添加湍流發(fā)生器后下游流場(chǎng)是否真的達(dá)到了湍流狀態(tài)進(jìn)行研究。

圖1 由于壁面粗糙引起的流向渦導(dǎo)致傳熱系數(shù)St曲線(xiàn)上升Fig．1 Increasing of St led by stream－wise vortex caused by the wall roughness

本文利用充分發(fā)展湍流的相似性特性和可維持特性來(lái)研究強(qiáng)迫轉(zhuǎn)捩發(fā)生的可能性。采用空間模式直接數(shù)值模擬（SDNS）的方法，在入口處添加包含湍流信息的擾動(dòng)，如果在下游流場(chǎng)能夠恢復(fù)湍流狀態(tài)，說(shuō)明強(qiáng)迫轉(zhuǎn)捩有可能發(fā)生，此時(shí)入口擾動(dòng)與湍流信息存在一定的關(guān)系；如果在下游流場(chǎng)趨近層流狀態(tài)，說(shuō)明即使是在入口添加了包含湍流信息的擾動(dòng)也無(wú)法達(dá)到湍流狀態(tài)，強(qiáng)迫轉(zhuǎn)捩很有可能不會(huì)發(fā)生，此時(shí)需要考慮更下游的位置。

1 數(shù)值方法

1．1 時(shí)間模式直接數(shù)值模擬

DNS的控制方程為無(wú)量綱的守恒型的N－S方程，不做任何假設(shè)。在結(jié)構(gòu)化的網(wǎng)格上采用有限差分法進(jìn)行空間離散，并采用高精度的緊致差分格式計(jì)算；時(shí)間上采用滿(mǎn)足TVD（Total Variation Diminished）特性的高階Runge－Kutta方法離散；為了使用迎風(fēng)差分格式，對(duì)非線(xiàn)性項(xiàng)還根據(jù)Jacobian矩陣特征值的正負(fù)進(jìn)行通量分裂。除了給出合理的初始條件外，還需要給定邊界條件。對(duì)于邊界層流而言，在流向和展向均為周期性邊界條件，而在法向，一個(gè)是固壁一個(gè)是出流條件。

采用時(shí)間模式直接數(shù)值模擬（TDNS）的好處是計(jì)算量小，可以一直從層流計(jì)算到充分發(fā)展湍流；其缺點(diǎn)是流向存在局部平行流假設(shè)。研究結(jié)果表明，對(duì)于零壓力梯度的高速平板邊界層，其雷諾數(shù)很大，空間模式直接數(shù)值模擬（SDNS）和TDNS的計(jì)算結(jié)果相比，在定性和定量上的差別都很?。?］。

1．2 空間模式直接數(shù)值模擬

與TDNS不同，空間模式直接數(shù)值模擬（SDNS）在流向還需要入流條件，而入流條件應(yīng)能反映上游的流動(dòng)特性。對(duì)于湍流計(jì)算而言，入流條件如何給是一大難題。

1．3 采用TDNS結(jié)果作為SDNS的入流條件

采用TDNS在某一時(shí)刻的結(jié)果，這時(shí)有確定的邊界層厚度，有相應(yīng)的Reynolds數(shù)。利用充分發(fā)展湍流的相似性特性，對(duì)TDNS在某一時(shí)刻的結(jié)果在法向進(jìn)行相似變換，就可以得到指定Reynolds數(shù)下的流場(chǎng)，這一流場(chǎng)的Reynolds數(shù)也就是SDNS入口處的Reynolds數(shù)。

轉(zhuǎn)換的步驟如下［7－8］：（1）采用 TDNS計(jì)算得到零壓力梯度的平板邊界層的充分發(fā)展湍流流場(chǎng)，取某一時(shí)刻的結(jié)果XT（x，y，z），其湍流邊界層厚度為δT；（2）采用數(shù)值方法求出具體問(wèn)題的層流解Xc（x，y，z），其計(jì)算域入口的層流邊界層厚度為δC；（3）假定流體在計(jì)算域入口處之前剛剛完成轉(zhuǎn)捩，則在計(jì)算域入口已是湍流，其邊界層厚度變?yōu)閗δδC，其中kδ取2到3．5；（4）將TDNS的結(jié)果在法向做相似變換XT（x，z），其中＝y(tǒng)kδδC／δT，并分解為平均量與脈動(dòng)量之和；（5）在湍流邊界層內(nèi)采用TDNS的平均流場(chǎng)，在湍流邊界層外采用SDNS的層流解 Xc（x0，y，z），中間采用光滑函數(shù)進(jìn)行連接，從而得到SDNS入口處的平均流場(chǎng)；（6）以適當(dāng)方式將TDNS的脈動(dòng)量X′T空間分布轉(zhuǎn)換成SDNS所需的入口流場(chǎng)的時(shí)間序列，從而得到SDNS入口處的脈動(dòng)流場(chǎng)。

轉(zhuǎn)換的示意圖如圖2所示。將轉(zhuǎn)換關(guān)系寫(xiě)成x＝ct，其中x的原點(diǎn)在TDNS中流場(chǎng)的出口處，方向朝向上游。這樣做的基本概念是流場(chǎng)中x處的擾動(dòng)以速度c向下游傳播，在t時(shí)刻將到達(dá)出口處，而這正好可以作為SDNS的入口擾動(dòng)。當(dāng)x的值大于TDNS的流向計(jì)算域長(zhǎng)度后，可利用流向周期性邊界條件的性質(zhì)，減去流向計(jì)算域長(zhǎng)度即可。這里取c＝0．92［7］。

該方法已經(jīng)被證實(shí)是解決空間模式直接數(shù)值模擬（SDNS）入口條件的一種有效方法，并且對(duì)不同Mach數(shù)、Reynolds數(shù)、壁面條件、典型邊界層類(lèi)型均可行［7－8］。

圖2 由TDNS得到的某時(shí)刻充分發(fā)展湍流流場(chǎng)轉(zhuǎn)換成SDNS的入口流場(chǎng)的示意圖Fig．2 Schematic diagram of the transformation from TDNS to SDNS

2 物理模型和計(jì)算工況

2．1 二維層流解的計(jì)算

本文選取某高超聲速飛行器的對(duì)稱(chēng)面作為研究對(duì)象，物理模型的頭部為圓弧，下表面有三個(gè)坡道。計(jì)算域如圖3所示。以來(lái)流速度U0＝1789m／s和頭部半徑r＝2．5mm為基本參量進(jìn)行無(wú)量綱化，無(wú)量綱參數(shù)Ma＝6，單位Re＝5．0×106，攻角為0°。采用SDNS的方法計(jì)算得到二維層流解，其網(wǎng)格參數(shù)為：流向2300個(gè)點(diǎn)；法向201個(gè)點(diǎn)，頭部激波內(nèi)部的點(diǎn)數(shù)不少于170個(gè)，靠近壁面的1個(gè)邊界層名義厚度內(nèi)的點(diǎn)數(shù)不少于75個(gè)，距離壁面最近的網(wǎng)格尺度不超過(guò)0．002個(gè)無(wú)量綱長(zhǎng)度。

2．2 三維流場(chǎng)的計(jì)算

截取二維SDNS中下表面某坡道內(nèi)的層流流場(chǎng)，并在展向做擴(kuò)展，作為三維SDNS的層流解，其中x方向沿著壁面向右，y方向垂直壁面向下，z方向垂直xy平面向外。三維SDNS中的入口位置是指入口壁面處在二維SDNS計(jì)算域的x坐標(biāo)值（圖3（b）所示），流向長(zhǎng)度指x方向的長(zhǎng)度，展向長(zhǎng)度指z方向的長(zhǎng)度，為21個(gè)無(wú)量綱長(zhǎng)度，有128個(gè)網(wǎng)格點(diǎn)。本文選取文獻(xiàn)［7］中的平板湍流邊界層TDNS的結(jié)果，結(jié)合SDNS入口處層流邊界層厚度，按照1．3節(jié)給出三維SDNS的入流條件。

圖3 計(jì)算域示意圖Fig．3 Schematic diagram of computational domain

2．3 計(jì)算工況

在下表面某位置添加湍流發(fā)生器后，此湍流發(fā)生器為下游流體提供了有限幅值的擾動(dòng)，有可能導(dǎo)致強(qiáng)迫轉(zhuǎn)捩。如果在下游某個(gè)位置發(fā)生了轉(zhuǎn)捩，那么在轉(zhuǎn)捩后的湍流流場(chǎng)應(yīng)該能夠維持湍流狀態(tài)。

圖4 第一個(gè)坡道的計(jì)算結(jié)果Fig．4 Results in the first ramp

入口處的不同位置代表湍流發(fā)生器的不同位置，入口處加擾動(dòng)的強(qiáng)度代表發(fā)生器不同的形狀；通過(guò)考察湍流的可維持特性來(lái)研究湍流發(fā)生器對(duì)流場(chǎng)的影響，如果在入口引入湍流信息后，流場(chǎng)很快能調(diào)整到湍流狀態(tài)，說(shuō)明在計(jì)算域內(nèi)通過(guò)加湍流發(fā)生器能夠達(dá)到湍流狀態(tài)，否則在計(jì)算域范圍內(nèi)難以達(dá)到湍流狀態(tài)。本文考察了入口處不同位置和不同湍流度的情況，共9種情況，詳情見(jiàn)表1，其中最大脈動(dòng)速度均方根是指三個(gè)方向脈動(dòng)速度的平方和的均方根沿法向剖面的最大值。

表1 計(jì)算工況Table 1 Computational cases

3 結(jié)果分析

3．1 第一個(gè)坡道的結(jié)果

圖4（a）和圖4（b）給出了第一個(gè)坡道內(nèi)3個(gè)算例的壁面摩擦系數(shù)沿流向的分布、流向平均速度剖面在靠近出口處沿法向的分布?？梢钥闯觯闆rA和B的壁面摩擦系數(shù)逐漸減小，并趨近于層流狀態(tài)的值，而且入口位置越靠前，越趨近層流解。從圖4（b）可以看出，三種情況的平均速度剖面不同程度地偏離了層流解，并且入口位置越靠前、入口湍流脈動(dòng)幅值越小，平均速度剖面越不飽滿(mǎn)。這表明在第一坡道內(nèi)引入較小幅值的擾動(dòng)，流場(chǎng)會(huì)趨近于層流狀態(tài)。

情況C中的壁面摩擦系數(shù)在x＝0．28m之前是快速減小，而在之后減小的幅度較小，并且其幅值接近湍流狀態(tài)的幅值。此外在情況A、B和C三種的流向平均速度剖面來(lái)看，情況C的流向平均速度剖面最飽滿(mǎn)，有可能是湍流狀態(tài)。

圖4（c）給出了第一個(gè)坡道內(nèi)情況C中不同位置的湍流Mach數(shù)沿法向的分布，可以看出在x＝0．28m之后的湍流Mach數(shù)幾乎重合。根據(jù)湍流的相似特性，說(shuō)明此時(shí)流場(chǎng)達(dá)到了湍流狀態(tài)。這表明在第一坡道內(nèi)引入較大幅值的擾動(dòng)，流場(chǎng)有可能趨近湍流狀態(tài)，強(qiáng)迫轉(zhuǎn)捩有可能發(fā)生。

3．2 第二個(gè)坡道的結(jié)果

圖5（a）和圖5（b）給出了第二個(gè)坡道內(nèi)3個(gè)算例的壁面摩擦系數(shù)沿流向的分布、流向平均速度剖面在靠近出口處沿法向的分布?？梢钥闯?，添加擾動(dòng)的入口位置越靠后，壁面摩擦系數(shù)偏離層流解越多，平均速度剖面越飽滿(mǎn)。當(dāng)擾動(dòng)的脈動(dòng)速度均方根值最大值不小于0．15時(shí)，下游的流場(chǎng)均趨近湍流狀態(tài)。

圖5（c）給出了第二個(gè)坡道內(nèi)情況C中不同位置的湍流Mach數(shù)沿法向的分布，可以看出在x＝0．54m之后的湍流Mach數(shù)也同樣具有相似性。這一點(diǎn)與圖5中情況C中x＞0．54m時(shí)壁面摩擦系數(shù)較平緩減小的結(jié)論是一致的，說(shuō)明此時(shí)流場(chǎng)達(dá)到了湍流狀態(tài)。

圖5 第二個(gè)坡道的計(jì)算結(jié)果Fig．5 Results in the second ramp

圖6 第三個(gè)坡道的計(jì)算結(jié)果Fig．6 Results in the third ramp

3．3 第三個(gè)坡道的結(jié)果

圖6 （a）和圖6（b）給出了第三個(gè)坡道內(nèi)3個(gè)算例的壁面摩擦系數(shù)沿流向的分布、流向平均速度剖在靠近出口處沿法向方向的分布。可以看出在入口處添加湍流較大幅值擾動(dòng)后，經(jīng)過(guò)一個(gè)快速的調(diào)整段，壁面摩擦系數(shù)在x＝0．67m開(kāi)始趨近湍流狀態(tài)的值，但如果擾動(dòng)的幅值較小，壁面摩擦系數(shù)任然會(huì)趨近于層流狀態(tài)的值。從流向平均速度剖面同樣也可以看出，情況A和B的剖面很飽滿(mǎn)，而情況C的剖面在靠近壁面處趨近于層流解。

圖6（c）給出了第三個(gè)坡道內(nèi)情況B中經(jīng)過(guò)Van Direst變換的平均速度剖面，可以看出在x＝0．68m到x＝0．70m的平均速度剖面幾乎重合，符合平均速度的相似性特性。并且平均速度剖面與線(xiàn)性律和對(duì)數(shù)律都吻合的很好，說(shuō)明情況B確實(shí)達(dá)到了湍流狀態(tài)。

4 結(jié) 論

采用空間模式直接數(shù)值模擬的方法，仿真了帶有三個(gè)壓縮折角的高超聲速飛行器下表面的湍流發(fā)生器，分析了其引起強(qiáng)迫轉(zhuǎn)捩的可能性，發(fā)現(xiàn)：

（1）在第一個(gè)坡道內(nèi)不會(huì)發(fā)生自然轉(zhuǎn)捩。

（2）在第一個(gè)坡道上，擾動(dòng)脈動(dòng)速度均方根值小于0．15時(shí)，流場(chǎng)趨近于層流狀態(tài)，不會(huì)發(fā)生強(qiáng)迫轉(zhuǎn)捩；如果大于0．225時(shí)，流場(chǎng)有趨近湍流狀態(tài)的趨勢(shì)，有可能發(fā)生強(qiáng)迫轉(zhuǎn)捩。

（3）在第二個(gè)和第三個(gè)坡道上，擾動(dòng)脈動(dòng)速度均方根值為0．15時(shí)，流場(chǎng)趨近湍流狀態(tài)，有可能發(fā)生強(qiáng)迫轉(zhuǎn)捩，但如果小于0．03時(shí)，流場(chǎng)有趨近層流狀態(tài)的趨勢(shì)，不會(huì)發(fā)生強(qiáng)迫轉(zhuǎn)捩。

［1］ TIRTEY S C，CHAZOT O．Characterization of hypersonic roughness induced transition for the EXPERT flight experiment［R］．AIAA 2009－7215．

［2］ BERRY S，DARYABEIGI K．Preliminary analysis of boundary layer transition on X－43Aflight 2［A］．52nd JANNAF Propulsion Meeting［C］，2004．

［3］ BERRY S，DARYABEIGI K，WURSTER K．Boundary layer transition on X－43A［R］．AIAA 2008－3736．

［4］ LIU J T．Nonlinear instability of developing stream－wise vortices with applications to boundary layer heat transfer intensification through an extended Reynolds analogy［J］．Phil．Trans．R．Soc．A．，2008，366：2699－2716．

［5］ MOMAYEZ L，DUPONT P，PEERHOSSAINI H．Some unexpected effects of wavelength and perturbation strength on heat transfer enhancement by Go¨rtler instability［J］．International Journal of Heat and Mass Transfer，2004，47：3783－3795．

［6］ 黃章峰，曹偉，周恒．超聲速平板邊界層轉(zhuǎn)捩中層流突變?yōu)橥牧鞯臋C(jī)理－－時(shí)間模式［J］．中國(guó)科學(xué)（G 輯），2005，35（5）：537－547．

［7］ 黃章峰，周恒．湍流邊界層的空間模式DNS的入口條件［J］．中國(guó)科學(xué)（G輯），2008，38（3）：310－318．

［8］ 董明，周恒．超聲速鈍錐湍流邊界層DNS入口邊界條件的研究［J］．應(yīng)用數(shù)學(xué)和力學(xué)，2008，29（8）：893－904．