• <tr id="yyy80"></tr>
  • <sup id="yyy80"></sup>
  • <tfoot id="yyy80"><noscript id="yyy80"></noscript></tfoot>
  • 99热精品在线国产_美女午夜性视频免费_国产精品国产高清国产av_av欧美777_自拍偷自拍亚洲精品老妇_亚洲熟女精品中文字幕_www日本黄色视频网_国产精品野战在线观看 ?

    寬速域RANS-LES混合方法的發(fā)展及應(yīng)用

    2017-07-03 16:08:49肖志祥羅堃宇
    空氣動力學(xué)學(xué)報 2017年3期
    關(guān)鍵詞:邊界層湍流尺度

    肖志祥, 羅堃宇, 劉 健

    (清華大學(xué) 航天航空學(xué)院, 北京 100086)

    ?

    寬速域RANS-LES混合方法的發(fā)展及應(yīng)用

    肖志祥*, 羅堃宇, 劉 健

    (清華大學(xué) 航天航空學(xué)院, 北京 100086)

    傳統(tǒng)的雷諾平均方法(RANS)已經(jīng)不能滿足大范圍分離、激波振蕩、壓力脈動、動載荷等極端工況下的流動預(yù)測需求;大渦模擬(LES)、直接數(shù)值模擬(DNS)等方法資源耗費(fèi)多、效率低,離工程湍流問題仍較為遙遠(yuǎn)。RANS-LES混合方法結(jié)合了RANS高效率和LES高精度的特點(diǎn),近期有望大規(guī)模應(yīng)用到工程湍流問題中。首先對現(xiàn)有的RANS-LES混合方法進(jìn)行了歸類,對各自的構(gòu)造思想、特點(diǎn)進(jìn)行了分析。然后報告了脫體渦模擬(DES)類方法的發(fā)展歷程和現(xiàn)狀,討論了使用DES類方法計算分離流動時,對流項(xiàng)離散格式對分離特性、小尺度結(jié)構(gòu)及頻譜特性等的影響,并構(gòu)造了自適應(yīng)耗散函數(shù)。最后介紹了近年來國內(nèi)外RANS-LES混合方法在寬馬赫數(shù)范圍(馬赫數(shù)從0.1到20)內(nèi)的機(jī)理研究和工程應(yīng)用?,F(xiàn)有的以DES類方法為代表的RANS-LES混合方法能夠較為精細(xì)地模擬非定常大分離流動中的復(fù)雜現(xiàn)象,但在計算效率等方面還有較大的改進(jìn)空間;植入式DES方法在模擬全機(jī)帶部件流動上具有較高的效率和模擬精度,是重要的發(fā)展方向。RANS-LES混合方法在動態(tài)失速、燃燒、氣動彈性、氣動噪聲、氣動光學(xué)等與非定常流動密切相關(guān)的方面也有廣闊的應(yīng)用前景。

    RANS-LES混合方法;脫體渦模擬;計算流體力學(xué);非定常流動;流動分離;高雷諾數(shù)湍流

    0 引 言

    計算流體力學(xué)(Computational Fluid Dynamics,CFD)在航空航天領(lǐng)域的應(yīng)用在近幾十年得到了迅猛發(fā)展,大大縮短了飛行器的設(shè)計周期,縮減了研發(fā)成本,同時幫助設(shè)計人員深入地理解流動機(jī)理、探索氣動設(shè)計的新思路[1-2]。

    在實(shí)際工程問題中,湍流常常主導(dǎo)著整個流動,決定飛行器的氣動特性。因此,湍流模擬方法是飛行器氣動設(shè)計、復(fù)雜湍流模擬中的關(guān)鍵技術(shù)。近年來,基于傳統(tǒng)渦粘性假設(shè)的湍流模型在追求設(shè)計點(diǎn)性能的飛行器氣動設(shè)計中取得了巨大的成功[3-5]。然而實(shí)踐表明,僅僅關(guān)注設(shè)計點(diǎn)是遠(yuǎn)遠(yuǎn)不夠的。在非設(shè)計點(diǎn)甚至極端工況下的復(fù)雜流動現(xiàn)象,例如邊界層非定常分離、旋渦破裂、激波/邊界層干擾、流動轉(zhuǎn)捩、壓力脈動和動載荷等,以及這些復(fù)雜非定常流動導(dǎo)致的強(qiáng)度、振動、噪聲、氣動熱及熱防護(hù)等問題,可能嚴(yán)重影響飛行性能,甚至危及飛行安全。因此,研究、發(fā)展并應(yīng)用高精度非定常湍流預(yù)測方法就成為了學(xué)術(shù)界和工業(yè)界的共同目標(biāo)。

    1 高精度湍流預(yù)測的現(xiàn)狀與前景

    直接求解Navier-Stokes(N-S)方程的方法即直接數(shù)值模擬(Direct Numerical Simulation,DNS)不需要引入任何額外的假設(shè)和模型,所有尺度的湍流信息都用足夠密的網(wǎng)格和時間步長直接解析出來。隨著雷諾數(shù)的增加,DNS對于網(wǎng)格量的要求非常高(所需網(wǎng)格量為Re9/4),對于實(shí)際的大雷諾數(shù)工程問題(Re≈1×107),計算所需網(wǎng)格總量在1×1014以上。在目前的計算能力限制下,采用DNS解決工程湍流問題,幾乎是不可能的。

    迄今為止,求解雷諾平均的N-S(Reynolds-Averaged N-S,RANS)方程組仍是計算工程湍流問題最常用的方法,它將湍流信息分解為時間無關(guān)的定常部分和時間相關(guān)的脈動部分,即:

    (1)

    上波浪線表示密度加權(quán)的Favre平均。將上式代入N-S方程后即可得到RANS方程,其動量方程寫作:

    (2)

    (3)

    其認(rèn)為雷諾應(yīng)力可類比于分子粘性應(yīng)力,并與平均場的應(yīng)變率成正比;通過引入等效的渦粘性系數(shù)將其封閉,將核心問題變?yōu)榍蠼鉁u粘系數(shù)μt,此后便出現(xiàn)了大量基于線性渦粘假設(shè)的湍流模型。在工業(yè)界得到最廣泛應(yīng)用的是一方程Spalart-Allmaras(S-A)模式[6]、由Menter提出的兩方程k-ωShear-Stress Transport(SST)模式[7]等,然而它們的局限在于適用面窄、經(jīng)驗(yàn)性強(qiáng),對大分離流動的預(yù)測能力差,無法準(zhǔn)確獲得頻譜特征,難以精細(xì)預(yù)測復(fù)雜流動。

    大渦模擬(Large-Eddy Simulation,LES)方法采取部分?;乃悸?,通過設(shè)置某種過濾器,將流場信息分解成隨時間變化的大周期低頻波動和高頻脈動;大尺度湍流的能量及行為直接解析得到;小尺度湍流脈動(小于慣性子區(qū)的湍流結(jié)構(gòu))由于遵守更為簡單的規(guī)律,可通過統(tǒng)一的方法來?;?,如經(jīng)典的Smagorinsky模型[8]。但是,在邊界層內(nèi)大部分湍流能量集中在高頻小尺度渦,采用LES方法求解需要的計算量與DNS相當(dāng),故實(shí)際應(yīng)用中不得不引入經(jīng)驗(yàn)性較強(qiáng)的壁面模型以模擬邊界層內(nèi)的LES附加應(yīng)力分布。

    Spalart估計了不同湍流模擬方法計算非定常工程湍流問題(雷諾數(shù)為106~107量級)所需的計算資源和運(yùn)用到工程實(shí)踐中的時間[9],見表1。RANS-LES混合方法在計算量上相對URANS方法略有增加,但相對與DNS和LES仍然具有相當(dāng)大的優(yōu)勢。近年來芯片行業(yè)的摩爾定律正在逐步失效,計算機(jī)芯片計算能力的提升速度可能會越來越慢,DNS和LES的實(shí)用化時間還可能進(jìn)一步推遲。美國國家航空航天局(NASA)在《CFD Vision 2030 Study》報告中指出[10],由于無法準(zhǔn)確可靠地預(yù)測大分離湍流流動,CFD技術(shù)在航空航天設(shè)計中的應(yīng)用已經(jīng)受到很大的局限;RANS方法的發(fā)展難以取得突破,DNS或LES方法的工程應(yīng)用在可預(yù)見的將來也不現(xiàn)實(shí);RANS-LES混合方法和壁面?;腖ES(Wall-modeled LES,WM-LES)[11]最有希望打破目前的困局。因此,在未來相當(dāng)長的一段時間內(nèi),RANS-LES混合方法可能是極少數(shù)能夠滿足工程應(yīng)用需求的高精度湍流預(yù)測手段,值得我們進(jìn)行深入的研究。

    表1 不同湍流預(yù)測方法達(dá)到工程實(shí)用化 所需的資源和預(yù)計時間[9]Table 1 Resources needed and estimated ready time of application for different turbulence simulation strategies

    2 RANS-LES混合方法的分類及國內(nèi)外研究現(xiàn)狀

    廣義而言,RANS-LES混合方法包括了所有在單個數(shù)值模擬中同時采用了RANS和LES方法來獲得解析湍流的計算方法。經(jīng)過近20年的研究,RANS-LES混合方法已經(jīng)發(fā)展出了非常多的分支。下面簡要介紹一下目前在科學(xué)研究和工程應(yīng)用中常見的幾種混合方法類型。

    2.1 湍流量混合模型

    由于時間平均的RANS方程和空間過濾的LES方程在形式上具有相似性,故可以在不同的流動區(qū)域引入不同的權(quán)重,對RANS和LES計算的湍流量進(jìn)行加權(quán),如對?;瘧?yīng)力加權(quán):

    (4)

    也可以直接基于渦粘性進(jìn)行加權(quán)[12]:

    (5)

    另一種思路是直接在雷諾應(yīng)力上引入一個貢獻(xiàn)函數(shù)。在不同的流動區(qū)域,RANS湍流模型對?;瘧?yīng)力的貢獻(xiàn)應(yīng)當(dāng)不一樣。例如Speziale提出的一種貢獻(xiàn)函數(shù)的形式[13]:

    (6)

    貢獻(xiàn)函數(shù)的本質(zhì)也是引入網(wǎng)格過濾尺度,使得在網(wǎng)格尺度足夠小的地方(例如小于Kolmogorov長度尺度η=(ν3/ε)1/4)減少雷諾應(yīng)力的?;?,保證足夠多的湍流信息能夠解析出來。針對非線性粘性問題有更復(fù)雜的改進(jìn)形式,如限制數(shù)值尺度(Limited Numerical Scale,LNS)方法[14]。

    這類方法的關(guān)鍵問題是對權(quán)函數(shù)的構(gòu)造,通常與網(wǎng)格尺度、壁面距離、湍流尺度等有關(guān);但權(quán)函數(shù)通常以代數(shù)形式出現(xiàn),經(jīng)驗(yàn)性過強(qiáng)。

    2.2 RANS-LES界面模型

    RANS-LES界面模型指的是在實(shí)際的計算流場中,近壁區(qū)域采用RANS計算,分離區(qū)域采用LES。其核心問題是RANS-LES的界面如何設(shè)置以及界面兩側(cè)如何連續(xù)過渡。具有代表性的方法包括脫體渦模擬(DES)類方法、分區(qū)混合方法和RANS限制的LES方法。

    2.2.1 DES類方法

    脫體渦模擬(Detached-Eddy Simulation,DES)類方法,是目前應(yīng)用最廣泛的RANS-LES混合方法。DES類方法構(gòu)造形式極為簡單、對復(fù)雜外形的適應(yīng)能力強(qiáng)。DES類方法通過引入網(wǎng)格尺度過濾,不同的計算區(qū)域采用不同的湍流模型,分界面動態(tài)變化,在壁面附近用RANS來模擬小尺度的湍流,分離區(qū)域采用類似Smagorinsky的亞格子應(yīng)力模型。Spalart在2009年對DES類方法進(jìn)行了比較全面深入的綜述[15]。后文將詳細(xì)闡述DES類方法的發(fā)展。

    2.2.2 分區(qū)混合方法

    分區(qū)混合方法的RANS-LES界面是人為指定的,典型代表是Deck等提出的分區(qū)DES(Zonal DES,ZDES)[16]。ZDES的基本框架建立在對分離流動的三種區(qū)分上(如圖1):

    (1) 幾何誘導(dǎo)的流動分離。這類流動分離由幾何形狀的突變導(dǎo)致,分離點(diǎn)的邊界層厚度δ遠(yuǎn)小于分離的幾何尺度H,如深腔流動。

    (2) 壓力梯度導(dǎo)致的流動分離。這類流動的分離主要由當(dāng)?shù)氐哪鎵禾荻葘?dǎo)致,在分離點(diǎn)的邊界層厚度δ仍然遠(yuǎn)小于分離的幾何尺度H,如翼型吸力面后緣的分離。

    (3) 與來流湍流邊界層性質(zhì)密切相關(guān)的流動分離。這類流動在分離點(diǎn)的邊界層厚度δ與分離的幾何尺度H相當(dāng),如淺腔流動。

    圖1 ZDES方法對分離流動的分類[17]Fig.1 Classification of separation by ZDES[17]

    ZDES方法針對這三種分離流動,分別定義不同的混合長度尺度和網(wǎng)格尺度:

    (7)

    圖2 三段翼型的ZDES分區(qū)[17]Fig.2 ZDES zones for the three-element airfoil[17]

    分區(qū)方法要求使用者在計算前即對流動情況有非常清晰的認(rèn)識,經(jīng)驗(yàn)性比較強(qiáng)。對于類型III的流動,需要在計算域邊界生成邊界層內(nèi)的湍流脈動信息才能正確預(yù)測分離,于是發(fā)展出了后文將要討論的植入式方法。

    2.2.3 RANS限制的LES方法

    RANS限制的LES方法(RANS-Limited LES)從LES出發(fā),用RANS的?;瘧?yīng)力限制LES在壁面附近的?;瘧?yīng)力。Delanghe[18-19]渦粘性假設(shè)引入亞格子湍動能kτ和耗散率ετ,構(gòu)造了適用于RANS和LES的統(tǒng)一渦粘系數(shù)形式:

    (8)

    (9)

    (10)

    2.3 第二代URANS方法

    所謂的第二代URANS方法在形式上與傳統(tǒng)的URANS一致,但分離區(qū)的雷諾應(yīng)力或多或少表現(xiàn)出LES的特征,因此也可以歸為RANS-LES混合方法的一類。具有代表性的包括由Girimaji[21]提出的部分平均的N-S方程(Partially-Averaged N-S,PANS)和Menter等[22-23]發(fā)展的尺度自適應(yīng)模擬(Scale-Adaptive Simulation,SAS)。

    PANS方法從k-ε湍流模式導(dǎo)出,采用兩個參數(shù)來控制計算所采用的模型:

    fk=ku/k

    (11)

    其中fk和fε分別代表未解析的(用下標(biāo)u表示)湍動能和耗散與總湍動能與耗散的比值。當(dāng)fk=1時,所有湍動能都未解析,即RANS;當(dāng)fk=0時,所有湍動能都解析,即DNS;當(dāng)0

    SAS方法的核心是在SST模式中引入與速度二階導(dǎo)數(shù)相關(guān)的von Karman長度尺度LvK:

    (12)

    對于現(xiàn)有的大部分兩方程模型而言,除k方程外,第二個方程大都參照k方程的形式進(jìn)行構(gòu)造,并非準(zhǔn)確的輸運(yùn)方程,很可能會遺漏必要的項(xiàng)或損失重要的物理機(jī)制。SAS方法將LvK整合到ω方程中得到額外的生成項(xiàng):

    (13)

    SAS方法能根據(jù)解析出的渦尺度動態(tài)調(diào)整?;耐牧鞒叨?,在分離區(qū)減小渦粘系數(shù),避免了傳統(tǒng)URANS方法只能得到單一模態(tài)渦結(jié)構(gòu)的局限;同時沒有顯式地引入網(wǎng)格過濾尺度,在較粗的網(wǎng)格上也能得到比以往更加精細(xì)的流動結(jié)構(gòu)。國內(nèi)外近年來對SAS方法的驗(yàn)證與應(yīng)用研究比較活躍,如[26-30]等。

    2.4 植入式方法

    植入式方法是指在RANS計算域中植入一個使用高精度非定常湍流預(yù)測方法的子計算域;其中非定常湍流預(yù)測方法的選取則很靈活,例如純LES、RANS-LES混合方法、DNS等。與傳統(tǒng)RANS-LES混合方法降低非定常預(yù)測成本的目標(biāo)不同,植入式方法的目標(biāo)是提高對局部流動物理機(jī)制的解析和預(yù)測能力。因此,在分區(qū)混合方法基礎(chǔ)上,植入式方法在RANS-LES分區(qū)邊界要傳遞湍流脈動信息。

    植入式方法是RANS-LES混合方法研究的熱門方向,有望解決傳統(tǒng)RANS-LES類混合方法存在的灰區(qū)問題。目前的研究集中于如何在RANS和LES的邊界上進(jìn)行湍流數(shù)據(jù)的重構(gòu)和傳遞。一方面,LES計算得到的數(shù)據(jù)要經(jīng)過統(tǒng)計平均得到URANS所需的邊界信息;另一方面,RANS計算的流場缺乏非定常信息,需要從湍流統(tǒng)計量中重構(gòu)出LES所需的湍流脈動量,而這樣的重構(gòu)既要保證與真實(shí)湍流足夠相似,又不能付出過大的計算代價。歐盟Go4Hybrid項(xiàng)目[31]對植入式方法進(jìn)行了較為系統(tǒng)的研究;Shur等開發(fā)了植入式方法所需的循環(huán)[32]及合成湍流技術(shù)[33]。陳海昕等[34]基于k-kL模式開發(fā)了窗口嵌入式(Window-Embedded)RANS-LES方法;甕哲等[35]對植入式IDDES(EIDDES)界面合成湍流(圖3)和循環(huán)方法進(jìn)行了初步探索。李棟等[36]則基于SST-DES采用合成渦方法生成界面的湍流信息。

    3 DES類混合方法的發(fā)展

    3.1 原始DES方法

    Spalart等[37]在1997年提出了DES方法的思想(DES-97),即在湍動能方程耗散項(xiàng)中引入DES長度尺度LDES=min(dw,CDESΔ),其中dw為當(dāng)前網(wǎng)格單元中心到最近壁面的法向距離,網(wǎng)格尺度Δ=max(Δx,Δy,Δz)。其后Strelets[39]在兩方程模式基礎(chǔ)上引入RANS模型中的湍流長度尺度LRANS和LES長度尺度LLES=CDESΔ,提出了DES方法的一般形式(DES-2001):

    (14)

    由LES長度尺度的形式可以看出,系數(shù)CDES的作用類似于Smagorinsky模式中的系數(shù)Cs,需要用各向同性湍流衰減或槽道流動算例進(jìn)行標(biāo)定。

    圖3 各向異性合成湍流Fig.3 Anisotropic synthetic turbulence

    由于原始DES方法中的RANS-LES界面完全由幾何外形和網(wǎng)格尺寸決定,當(dāng)加密邊界層內(nèi)流向或展向網(wǎng)格使得LRANS>Δ?LLES時,RANS-LES界面位于邊界層內(nèi)部,LES在邊界層內(nèi)啟動但又無法在較粗的網(wǎng)格上正確解析雷諾應(yīng)力,導(dǎo)致總雷諾應(yīng)力不足、產(chǎn)生非物理的網(wǎng)格誘導(dǎo)分離現(xiàn)象(Grid Induced Separation,GIS)。

    除此以外,DES類方法在RANS-LES界面上沒有向LES區(qū)提供充足的湍流脈動信息。在附著邊界層RANS計算到分離剪切層LES計算的轉(zhuǎn)換初期,部分區(qū)域處于RANS和LES之間的“灰區(qū)”,解析和?;膽?yīng)力都不足?;覅^(qū)問題對于大分離流動的影響并不是特別顯著;但對于淺臺階分離、射流自由混合層等依賴于剪切層內(nèi)湍流脈動特性的算例,其影響比較明顯。王翔宇等[38]用標(biāo)準(zhǔn)k方程亞格子模式將SST-DES重新構(gòu)造,得到了一種改進(jìn)的SST-DES方法,在AS239翼型小分離模擬中減輕了“灰區(qū)”問題。

    3.2 延遲的DES方法(DDES)

    針對原始DES方法產(chǎn)生網(wǎng)格誘導(dǎo)分離的問題,Spalart等[40]提出了延遲DES方法(Delayed DES,DDES),引入物理相關(guān)的延遲函數(shù)將LES區(qū)域阻隔到邊界層外部。其混合長度尺度寫作:

    LDDES=LRANS-fdmax(0,LRANS-LLES)

    (15)

    其中的fd即為延遲函數(shù),取值在0~1之間。原始的DES中,RANS到LES的轉(zhuǎn)換完全依靠網(wǎng)格尺度大??;而DDES在轉(zhuǎn)換判據(jù)中引入了當(dāng)?shù)亓鲃幼兞康挠绊?,fd函數(shù)在邊界層內(nèi)部接近1,邊界層外為0,從而可以推遲RANS向LES轉(zhuǎn)換。以上函數(shù)的具體構(gòu)造較為冗長,可參見文獻(xiàn)[40]。

    基于SST模式的DDES方法(DDES-2003)也有類似的形式[41]:

    (16)

    其中FSST可取為SST模式中的混合函數(shù)F1或F2。FSST在壁面附近接近1,在邊界層外為0,同樣可以在邊界層內(nèi)延遲RANS到LES的轉(zhuǎn)換。

    3.3 改進(jìn)的DDES方法(IDDES)

    (17)

    其形式與DDES接近,但修改了混合加權(quán)函數(shù)的定義:

    (18)

    其中fdt反映當(dāng)?shù)亓鲃佑绊?,fB反映幾何影響,fe=0時IDDES回退到DDES。fdt、fB及fe的具體形式在文獻(xiàn)[42]中可找到,本文不再贅述。

    另外IDDES中的RANS長度尺度LRANS和LES長度尺度LLES與DDES定義相同,但對網(wǎng)格尺度作了修改:

    (19)Δmin和Δmax分別為x/y/z方向網(wǎng)格尺度的最小和最大值,系數(shù)Cw=0.15。若流向網(wǎng)格尺寸為法向第一層的100~1000倍,壁面y+~1,法向網(wǎng)格拉伸率1.15,則新的網(wǎng)格尺度減小了y+=10~100處(即對數(shù)區(qū))的網(wǎng)格尺度,減小了亞格子粘性,消除了LLM現(xiàn)象。

    (20)

    WMLES機(jī)制的引入使得IDDES能夠響應(yīng)來流湍流信息,加快附著邊界層從RANS到分離區(qū)LES的轉(zhuǎn)換,緩解DES類方法的“灰區(qū)”問題。當(dāng)上游缺乏湍流脈動信息時IDDES則退回到DDES。圖4對比了間距L=3.7D的串列雙圓柱(TandemCylinders,TC-3.7)算例中兩個圓柱表面的壓力系數(shù)均方根分布。相比DDES-2003,IDDES預(yù)測結(jié)果與實(shí)驗(yàn)吻合較好,尤其對于后圓柱,IDDES結(jié)果幾乎與有轉(zhuǎn)捩帶的實(shí)驗(yàn)值完全一致。詳細(xì)比較可參見文獻(xiàn)[43]。由于IDDES方法最初是基于SA模式的,Gritskevich等[44]對基于SST模式的DDES和IDDES方法進(jìn)行了標(biāo)定和簡化,討論了在SA和SST等不同湍流模式中IDDES的形式和系數(shù)細(xì)節(jié)問題,但對結(jié)果并沒有本質(zhì)上的重大影響,有興趣的讀者可以參考。

    (a) Front cylinder

    (b) Rear cylinder圖4 DDES-2003和IDDES預(yù)測串列 雙圓柱表面壓力脈動對比Fig.4 Comparison between DDES-2003 and IDDESpredictions of surface pressure fluctuations for TC-3.7

    4 與RANS-LES混合方法匹配的數(shù)值方法

    作為預(yù)測非定常分離流動的手段,RANS-LES方法需要與合適的時間推進(jìn)和空間離散格式相匹配。實(shí)際工程湍流問題中飛行器外形和流動較為復(fù)雜,因此要求數(shù)值格式適用性強(qiáng)、魯棒性好;具有一定的時間精度和較高的并行效率;在物面、遠(yuǎn)場邊界、激波等處穩(wěn)定、在分離區(qū)能夠盡可能多解析流動結(jié)構(gòu)。

    4.1 時間推進(jìn)方法

    目前常見的是時間方向和空間方向分開的半離散格式,包括顯式推進(jìn)和隱式推進(jìn)。顯式時間推進(jìn)有多步龍格庫塔法[45-47]等,其形式簡單易于實(shí)現(xiàn)、可達(dá)到較高精度且易于并行計算;但受數(shù)值穩(wěn)定性的限制,統(tǒng)一時間步長必須很小,影響計算效率。隱式時間推進(jìn)有交錯對角迭代(ADI)[48]、近似因子(AF)[49]、上下-高斯賽德爾迭代(LU-SGS)[50-51]、數(shù)據(jù)并行-上下松弛(DP-LUR)[52-53]等,穩(wěn)定性要求相對寬松,可采用較大的時間步長,效率較高、工程實(shí)用性更強(qiáng);但是,隱式方法需要借助雙時間步子迭代等方法減小近似分解誤差,一般只有二階精度。對于傳統(tǒng)LES方法,隱式和顯式兩種時間推進(jìn)方法都可滿足預(yù)測精度的要求[54-55]。作為面向工程實(shí)際應(yīng)用的湍流數(shù)值預(yù)測手段,RANS-LES方法一般采用含子迭代的隱式時間推進(jìn)就可以滿足效率和精度的需求[56]。

    4.2 空間離散格式

    基于N-S方程組的特點(diǎn),其擴(kuò)散項(xiàng)一般采用中心格式離散,而對流項(xiàng)則有多種離散格式,其精度直接影響到數(shù)值解的特性。目前對流項(xiàng)空間離散格式基本上分為中心型和迎風(fēng)型兩大類。

    中心型格式最大的問題是其本身不具備格式耗散,容易造成計算不穩(wěn)定,甚至發(fā)散。以Jameson[47]為代表的中心型格式,在離散單元體界面上采用二階中心差分離散并引入人工粘性:

    (21)

    D(2)和D(4)分別為二階和四階人工粘性項(xiàng),前者為了捕捉激波,后者則在光滑區(qū)提供背景耗散,抑制數(shù)值振蕩。對于低速流動,可舍去二階人工粘性,將通量的差分精度提高至四階(記作C4格式)。中心型格式本身適合低速流動的模擬,但是人工粘性系數(shù)的經(jīng)驗(yàn)性較強(qiáng)。

    迎風(fēng)型格式考慮了特征線在界面上的傳播方向,物理意義相對明確,對流場間斷具有較好的捕捉能力,在高速流動中應(yīng)用廣泛,常見的有矢通量微分的Roe格式[57]、矢通量分裂的van Leer格式[58]等。對于LES區(qū)域內(nèi)需要解析的小尺度結(jié)構(gòu)而言,常規(guī)的迎風(fēng)型格式耗散過大[59],會嚴(yán)重抑制對小尺度結(jié)構(gòu)的解析。

    Bui[55]等在Roe格式耗散項(xiàng)前乘以一個小于1的常數(shù)φ降低迎風(fēng)格式的耗散:

    (22)

    基于對稱TVD格式[60]的思想,可將色散和耗散獨(dú)立出來單獨(dú)控制,寫出對稱通量和迎風(fēng)耗散混合格式的一般形式:

    Fi+1/2=Fsymmetric,i+1/2-

    (23)

    系數(shù)φ小于1,用于控制格式的迎風(fēng)耗散大小。對稱通量可以采用Roe格式或六階中心格式計算,界面左右變量qL和qR可由MUSCL插值、五階WENO插值[61]或MDCD插值[62]得到,例如可以用6階中心格式和5階WENO插值構(gòu)造適用于低速流動的S6WENO5格式。

    4.3 自適應(yīng)格式耗散

    雖然降低格式耗散有利于解析小尺度流動結(jié)構(gòu),然而整個流場統(tǒng)一降低格式耗散常引起數(shù)值振蕩。事實(shí)上,在流動的不同區(qū)域需要不同的數(shù)值耗散,即上文的常數(shù)φ應(yīng)該隨當(dāng)?shù)亓鲃幼赃m應(yīng)變化。一個合理的思路是發(fā)揮中心型和迎風(fēng)型格式各自的優(yōu)勢,引入加權(quán)混合的概念,使得在流場的不同區(qū)域采用不同的離散格式[63]:

    (24)

    其中權(quán)函數(shù)(亦稱為自適應(yīng)耗散函數(shù))φ在LES工作的分離區(qū)域趨近于0,使得格式耗散降低至接近中心格式的水平,盡可能解析小尺度湍流結(jié)構(gòu),排除數(shù)值耗散對物理耗散的影響;在壁面及遠(yuǎn)場無旋區(qū)域趨近于1,格式還原為穩(wěn)定的迎風(fēng)型格式,盡可能消除物面邊界和遠(yuǎn)場邊界帶來的數(shù)值誤差和振蕩,保證數(shù)值格式的穩(wěn)定性。公式(22)和(23)中的耗散系數(shù)φ都可以替換為這種自適應(yīng)耗散函數(shù)。

    借鑒SST模式中的混合函數(shù)可以構(gòu)造雙曲正切形式(圖5)的自適應(yīng)耗散函數(shù)[64]:

    φ=φmaxtanh(A3)

    A=max{[(CDESΔ/lturb)/g-0.5], 0}

    g=tanhB4

    (25)

    圖5 雙曲正切函數(shù)y=tanh(x)圖像Fig.5 Hyperbolic tangent function y=tanh(x)

    在遠(yuǎn)場無旋區(qū)域,CDESΔ>lturb,S>0,Ω?1,則B為小值,g接近0,因此A為較大值,函數(shù)φ接近于1;對于分離區(qū)域,S與Ω較為接近,B的值趨近于2,g接近1,CDESΔ

    仍以TC-3.7為例[65],渦粘系數(shù)和自適應(yīng)耗散函數(shù)的分布如圖6所示,在大分離區(qū)自適應(yīng)耗散函數(shù)值很小,對應(yīng)于較大的模化渦粘系數(shù), 小尺度結(jié)構(gòu)不會被耗散抹平。對比不同的格式(圖7),未降低耗散的原始Roe格式預(yù)測的前圓柱表面壓力脈動頻譜主頻大于實(shí)驗(yàn),幅值也有明顯差異;降低耗散后的S6WENO5格式與實(shí)驗(yàn)值的吻合度有顯著提高,后圓柱高頻部分衰減也較小。張揚(yáng)等采用基于Roe格式的二階精度耗散自適應(yīng)混合格式在單圓柱繞流和NACA0021大攻角分離流動中得到了較好的預(yù)測效果,計算結(jié)果受基準(zhǔn)湍流模式的影響較小。

    圖6 TC-3.7算例的渦粘系數(shù)μt和自適應(yīng)耗散φ分布Fig.6 Distribution of the eddy viscosity and adaptive function of TC-3.7

    對于低速流動,自適應(yīng)耗散函數(shù)最小值取0<φmin<0.1較好,實(shí)際計算中取0.05~0.1即可得到較理想的結(jié)果。

    對于存在流場間斷(如激波)的流動,在間斷附近需要足夠大的格式耗散來抑制數(shù)值振蕩,可將自適應(yīng)耗散函數(shù)φ改寫為:

    (26)

    其中φ1為上文所述的自適應(yīng)耗散,φ2為激波探測器[67]:

    (27)

    ‖Ω‖為渦量的模。在當(dāng)?shù)厮俣忍荻群艽筇帵?=1,提供足夠的格式耗散(圖8)。

    需要注意的是,格式耗散并非越低越好;過低的數(shù)值耗散會帶來明顯的非物理數(shù)值噪聲(如圖9中后圓柱剪切層),應(yīng)根據(jù)具體情況調(diào)整自適應(yīng)函數(shù)的下限值,詳細(xì)討論可參見文獻(xiàn)[65]。

    (a) Front cylinder, 135°

    (b) Rear cylinder, 45°圖7 不同格式得到的TC-3.7表面壓力脈動頻譜對比Fig.7 Comparisons of pressure fluctuation spectra by different schemes on the surface of TC-3.7

    (a) p

    (b) φ1

    (c) φ2

    (d) φ圖8 激波與對應(yīng)的激波探測器函數(shù)分布Fig.8 Shockwaves and the shockwave detector distribution

    (a) Expriment

    (b) φmin=0.1

    (c) φmin=0.03圖9 不同格式耗散下的串列雙圓柱瞬時展向渦量Fig.9 Instantaneous spanwise vorticity for Tandem cylinders

    5 RANS-LES混合方法在寬速域范圍內(nèi)的應(yīng)用

    RANS-LES混合方法經(jīng)過多年的發(fā)展,已經(jīng)在航天航空、航海、汽車等領(lǐng)域得到了較為廣泛的應(yīng)用(如F-15E的全機(jī)大攻角模擬[68]),可以勝任從完全不可壓到高超聲速的寬速域范圍實(shí)際工程流動的精細(xì)預(yù)測。早期的發(fā)展和應(yīng)用可參考Fr?hlich等[69]和Spalart[15]的綜述文章。此處僅介紹近幾年來國內(nèi)外在機(jī)理和工程應(yīng)用方面一部分有代表性的工作。

    5.1 低速流動

    在低速流動應(yīng)用中,采用RANS-LES方法耦合Ffowcs Williams-Hawkings(FW-H)方程[70]對流致噪聲進(jìn)行預(yù)測是較為熱門的研究方向。

    串列圓柱繞流是非常典型的低速大分離流動,根據(jù)圓柱間距不同具有差異很大的流動特征。Lockard等[71]研究了網(wǎng)格密度、展向?qū)挾葘Υ须p圓柱遠(yuǎn)場噪聲的影響并嘗試解釋了實(shí)驗(yàn)中出現(xiàn)的流動不對稱現(xiàn)象;之后的BANC-I Workshop對湍流模式、混合方法、網(wǎng)格密度等因素進(jìn)行了全方位的評估[72]。Weinmann等[73]比較了基于不同湍流模式的IDDES、SAS等多種RANS-LES混合方法對串列雙圓柱遠(yuǎn)場噪聲的影響;Xiao等[74]采用IDDES方法研究了窄間距串列三圓柱的流動特征。

    起落架是大型客機(jī)在起降過程中的重要噪聲源,也是RANS-LES混合方法研究的一個重要對象。胡寧等[75]比較了URANS、DES和DDES對RLG流動預(yù)測結(jié)果的影響,DDES結(jié)果較優(yōu);Xiao等[76]分析了波音基本起落架(Rudimentary Landing Gear,RLG)近場噪聲源的產(chǎn)生機(jī)理及其與表面流動拓?fù)涞膶?yīng)關(guān)系,IDDES結(jié)果比DDES更接近實(shí)驗(yàn)值;Spalart等[77]對比了不同積分面選取對遠(yuǎn)場噪聲的影響。關(guān)于LAGOON前起落架[78-79]和PDCC前起落架[80-81]的流動和噪聲也有相當(dāng)多的RANS-LES混合方法研究。 采用RANS-LES混合方法研究翼型和增升裝置在縫翼和襟翼處的渦干擾、尾緣分離等也很常見,包括DU96-W180翼型[82]、帶縫翼或襟翼的兩段翼型[83-84]、30P30N增升裝置構(gòu)型[85]等。Escobar等[86]完成了客機(jī)翼身組合體帶增升裝置的DES數(shù)值模擬。 背負(fù)式進(jìn)氣道加V形尾翼是長航時低可探測無人機(jī)的典型布局。若進(jìn)氣道堵塞,其溢流可能對尾翼結(jié)構(gòu)造成影響。采用添加背壓的方法模擬單側(cè)進(jìn)氣道堵塞后的情況,IDDES方法成功捕捉到了停車進(jìn)氣道內(nèi)部的非定常流動以及進(jìn)氣道口向外的溢流(圖10)。

    圖10 單側(cè)進(jìn)氣道堵塞溢流及其對V形尾翼的影響Fig.10 Inlet overflow and the interaction with V-shape tail

    從圖11中可看到,在出現(xiàn)溢流的右側(cè)進(jìn)氣道表面及其后的垂尾翼根處,表面壓力脈動要顯著高于未發(fā)生溢流的一側(cè)。

    圖11 左右進(jìn)氣道和垂尾表面動載荷頻譜對比Fig.11 Comparison of power spectral density (PSD) at surface sample points of the inlet and V-shape tail

    此外,IDDES方法在螺旋槳短艙相互干擾[87]、后臺階流動[88-89]、翼型大攻角分離和強(qiáng)迫振蕩[90]、雙三角翼渦破裂[91-92]等低速算例中也已經(jīng)得到了應(yīng)用。

    5.2 跨/超聲速流動

    采用RANS-LES混合方法模擬跨/超聲速流動可以較好地再現(xiàn)激波/邊界層干擾導(dǎo)致的分離和激波抖振、剪切層失穩(wěn)破碎等非定常流動特征。Gaitonde[93]和張偉偉等[94]在綜述中指出RANS-LES混合方法在激波/邊界層干擾流動中具有很好的應(yīng)用前景。

    Huang等[95]采用IDDES方法得到了在馬赫數(shù)0.73、攻角3.5°時OAT15A翼型表面激波前后抖振的形態(tài),用渦流發(fā)生器進(jìn)行控制取得了較為顯著的效果。 超聲速進(jìn)氣道內(nèi)部也有大量的激波/邊界層干擾和激波相互干擾現(xiàn)象。采用IDDES方法模擬某進(jìn)氣道典型不起動流場,來流馬赫數(shù)3.0,兩個燃料噴流柱附近的渦系分布如圖12所示。第一噴流柱附近的渦結(jié)構(gòu)較完整,非定?,F(xiàn)象較弱;第二噴流柱附近的渦結(jié)構(gòu)復(fù)雜,小尺度結(jié)構(gòu)較多,非定常效應(yīng)顯著,增混效果較明顯。這是由于第一噴流柱附近激波較弱,相應(yīng)的激波/激波、激波/邊界層干擾并不明顯;第二噴流柱附近存在復(fù)雜的激波/激波、激波/邊界層干擾(圖13),且處于第一噴流柱的尾跡影響范圍內(nèi),故小尺度結(jié)構(gòu)較多,利于油氣混合。

    (a) 第一噴流柱

    (b) 第二噴流柱圖12 噴流柱附近渦系分布Fig.12 Vortex structures around the flow injectors

    隨著軍用隱身飛機(jī)的發(fā)展,內(nèi)埋彈艙的非定常流動情況成為了工業(yè)界關(guān)注的焦點(diǎn);簡化的空腔流動具有彈艙流動的主要特征,可以反映流動的主要機(jī)理。關(guān)于空腔和彈艙的研究可以參看Lawson等[96]的綜述文章。M219空腔標(biāo)模[97]是一個被廣泛用于RANS-LES混合方法研究和標(biāo)定的算例,如歐盟DESider項(xiàng)目[98]。Barakos和Lawson[99-100]、Liggett等[101]、Temmerman等[102]和Babu等[103]分別利用M219空腔進(jìn)行了DES、SAS等方法和不同基準(zhǔn)湍流模式的評估。Luo等[104]對M219空腔的風(fēng)洞外形進(jìn)行了IDDES預(yù)測,在局部加密的網(wǎng)格上空腔底面壓力脈動均方根分布和頻譜與實(shí)驗(yàn)的吻合度有較大改善;劉瑜等[105]利用DDES方法研究了M219空腔前緣鋸齒擾流片控制方案;劉俊等[106]采用DDES方法研究了來流邊界層厚度對M219空腔非定常流動的影響,發(fā)現(xiàn)隨著邊界層厚度增加,空腔底部壓力脈動減弱,各階Rossiter模態(tài)頻率向低頻移動。

    (a) 第一噴流柱

    (b) 第二噴流柱圖13 噴流柱附近的波系對比Fig.13 Wave series around the flow injectors

    進(jìn)一步地,針對真實(shí)彈艙甚至全機(jī)帶彈艙外形的RANS-LES數(shù)值預(yù)測近年來也逐漸增多,如Lawson等[107]采用DES方法對1303 UCAV無人機(jī)彈艙模型進(jìn)行了多種馬赫數(shù)和艙門狀態(tài)的數(shù)值預(yù)測,考慮了艙門鋸齒邊緣、鉸鏈等部件,同時進(jìn)行了粗略的聲場計算。Kannepalli等[108]對F-35戰(zhàn)斗機(jī)彈艙風(fēng)洞模型進(jìn)行了Zonal RANS-LES模擬,研究了艙門開閉等不同狀態(tài)下彈艙內(nèi)部壓力脈動特性,并對前緣噴氣、隔音板等降噪措施效果進(jìn)行了評估。Kim等[109]把DES方法應(yīng)用于彈艙投彈的數(shù)值模擬。這些RANS-LES混合方法研究都顯示出了非常強(qiáng)的工程實(shí)用性,表明RANS-LES混合方法已經(jīng)逐漸成熟,已成為工程非定常流動預(yù)測的必備手段。

    Luo等開展了不規(guī)則外形彈艙流動機(jī)理研究[110]、鋸齒擾流片控制方案幾何參數(shù)優(yōu)選[111]等工作;也嘗試在全機(jī)復(fù)雜外形(包括考慮導(dǎo)彈及掛架等附屬物)上用IDDES方法評估鋸齒控制裝置的效果,總網(wǎng)格量超過4000萬。IDDES方法成功捕捉到了流動控制裝置對彈艙前緣剪切層的抬升作用(圖14)。當(dāng)導(dǎo)彈處于待發(fā)狀態(tài)時,剪切層與彈體頭部存在相互作用,可能會對彈體振動和發(fā)射產(chǎn)生影響。這是混合方法在極端復(fù)雜外形流動預(yù)測中非常成功的案例。

    (a) 不帶彈狀態(tài) (b) 懸掛狀態(tài) (c) 發(fā)射狀態(tài)

    圖14 復(fù)雜外形彈艙內(nèi)的瞬時展向渦量分布

    Fig.14 Instantaneous spanwise vorticity distributions inside the complex weapon bay

    突起物流動或翼身組合體流動中的馬蹄渦及背風(fēng)面大分離也適合RANS-LES混合方法的特點(diǎn)。應(yīng)用IDDES方法對馬赫數(shù)0.74來流中的突起物及加裝整流罩后的非定常流動、壓力脈動特性等進(jìn)行評估分析;加裝整流罩對抑制突起物及其附近的表面壓力脈動有比較顯著的效果(圖15)。整流前突起物后方大分離區(qū)有大片的高脈動區(qū)域,整流后顯著縮小,且離突起物本身比較遠(yuǎn),突起物和整流罩表面的壓力脈動并不大。

    (a) 無整流罩 (b) 有整流罩圖15 突起物及整流罩附近的表面壓力脈動均方根分布Fig.15 Surface pressure fluctuations near around the hump and fairing

    類似的工程應(yīng)用研究還有亞聲速段火箭助推器表面壓力脈動和非定常流動預(yù)測,如Pain等[112]和Heinrich等[113]分別將ZDES和IDDES應(yīng)用于對阿麗亞娜火箭助推段在跨聲速條件下的底部大分離流動進(jìn)行了預(yù)測并對比了不同箭體構(gòu)型的非定常氣動載荷抖振情況。

    5.3 高超聲速流動

    載人航天和高速全球投送的現(xiàn)實(shí)需求使得工業(yè)界對高超聲速流動越來越關(guān)注。載人返回艙、高速飛行器等都在高超聲速劇烈氣動加熱環(huán)境下工作。Brock等[114]對Ma=6.4的Orion返回艙進(jìn)行了DES模擬,比較了多種網(wǎng)格密度和空間離散格式的影響;結(jié)果表明DES方法得到的結(jié)果與實(shí)驗(yàn)?zāi)軌蛭呛希^密的網(wǎng)格和耗散較低的格式有利于計算結(jié)果的改善。Salazar等[115]也利用Orion返回艙算例對基于Menter的BSL模式的RANS-LES混合方法進(jìn)行了驗(yàn)證,與實(shí)驗(yàn)和DES計算結(jié)果相近。劉佳等[116]在DES方法中加入了可壓縮性修正,對Orion返回艙氣動力和氣動熱進(jìn)行了預(yù)測。董祥瑞等[117]采用DES方法對Ma=7流場中激波/邊界層干擾進(jìn)行了模擬并采用不同安裝位置的單/雙楔形單元進(jìn)行控制。Barnhardt等[118]使用DES方法對馬赫數(shù)高達(dá)20的大氣層再入裝置的底部流動進(jìn)行了模擬,計算與實(shí)驗(yàn)?zāi)軌蛭呛?,證明DES類方法可以解決高超聲速再入條件下的大分離流動和氣動熱問題。

    高超聲速條件下湍流帶來的氣動加熱比層流大很多。在風(fēng)洞實(shí)驗(yàn)中也需要在縮比模型上加裝粗糙單元模擬真實(shí)飛行中的流動轉(zhuǎn)捩現(xiàn)象。因此轉(zhuǎn)捩也是影響高超聲速飛行器表面流動和氣動力/熱的重要問題。由于轉(zhuǎn)捩模式的相對不成熟,這方面采用混合方法進(jìn)行研究的不多,一般需要在轉(zhuǎn)捩位置前人為將流動設(shè)置為層流;Duan等[119]對不同形狀和數(shù)目的粗糙單元強(qiáng)制轉(zhuǎn)捩進(jìn)行了IDDES模擬;Xiao等[120]采用IDDES方法研究了不同攻角下表面凹腔誘導(dǎo)流動轉(zhuǎn)捩的機(jī)理?;谵D(zhuǎn)捩模式的RANS-LES混合方法正在發(fā)展中。

    6 未來展望

    RANS-LES方法已被證明適用于全速域各種復(fù)雜外形的流動精細(xì)預(yù)測。未來的研究方向主要是完善現(xiàn)有方法的不足,進(jìn)一步提升計算效率;目前DDES、IDDES等先進(jìn)DES類方法還沒有在工程應(yīng)用中大面積普及,還需要進(jìn)一步推動,推廣工程應(yīng)用。

    完善現(xiàn)有RANS-LES混合方法是當(dāng)前的研究重點(diǎn)之一,如在非線性模式或轉(zhuǎn)捩模式基礎(chǔ)上開發(fā)RANS-LES混合方法;或改進(jìn)網(wǎng)格尺度定義使其適應(yīng)準(zhǔn)二維流動展向網(wǎng)格分布[121]:

    (28)

    N=(Nx,Ny,Nz)是沿當(dāng)?shù)販u量方向的單位向量。

    第二代URANS方法如SAS等相對于傳統(tǒng)的RANS方法具備更強(qiáng)的非定常預(yù)測能力,對網(wǎng)格依賴更小、計算效率更高,國外研究中已經(jīng)將其應(yīng)用于戰(zhàn)斗機(jī)全機(jī)復(fù)雜流動的快速預(yù)測[28]。

    植入式方法在RANS-LES界面上添加構(gòu)造的湍流擾動信息(圖16),加速剪切層從RANS到LES的轉(zhuǎn)換,有望解決存在多年的“灰區(qū)”問題,擴(kuò)大RANS-LES混合方法的適用范圍、提高預(yù)測精度。植入式方法可將非定常預(yù)測區(qū)域限制在更小的關(guān)注范圍內(nèi),為全機(jī)與局部細(xì)節(jié)之間差異巨大的尺度提供了統(tǒng)一處理的途徑,可用于渦流發(fā)生器等流動控制裝置的設(shè)計和機(jī)理研究。人工構(gòu)造的湍流擾動發(fā)展為真實(shí)湍流需要一定的過渡區(qū)域,如何縮短過渡區(qū)域長度是目前重點(diǎn)研究的方向。另一方面也有必要將自由來流湍流度的影響考慮進(jìn)來,解決鈍體迎風(fēng)面附著流區(qū)域脈動量偏低以及流動轉(zhuǎn)捩等問題。

    RANS-LES混合方法在學(xué)科交叉耦合中將得到越來越多地應(yīng)用,如燃燒和化學(xué)反應(yīng)[122-124]、氣動彈性[125]、氣動噪聲[126]、氣動光學(xué)[127]等方面的非定常預(yù)測。

    圖16 燃燒室流動在入口添加合成湍流[35]Fig.16 Synthetic turbulence before a combustor[35]

    7 結(jié) 論

    本文介紹了適用于復(fù)雜外形非定常大分離湍流流動預(yù)測的RANS-LES混合方法、與之匹配的數(shù)值格式和耗散,及近年來國內(nèi)外在RANS-LES混合方法研究及工程應(yīng)用中的部分代表性成果。

    1) RANS-LES混合方法的核心思想是綜合RANS和LES兩種方法的優(yōu)點(diǎn),在流動復(fù)雜區(qū)域采用LES方法盡可能解析豐富的流動結(jié)構(gòu),在邊界層內(nèi)部或不關(guān)心的流動區(qū)域采用RANS方法降低資源需求、提高計算效率。

    2) RANS-LES混合方法要求時間和空間離散格式都具備較高的精度。時間方向可以采用統(tǒng)一時間步長的二階精度隱式時間推進(jìn)??臻g離散格式除了要具備較高的離散精度以外,還要具備與流動相適應(yīng)的自適應(yīng)數(shù)值耗散。在壁面、遠(yuǎn)場無旋區(qū)及流動間斷處,數(shù)值耗散要足夠大以抑制數(shù)值振蕩、維持計算穩(wěn)定;在大分離區(qū),要減小數(shù)值耗散以便解析小尺度流動結(jié)構(gòu),壁面數(shù)值耗散對物理耗散的影響。

    3) DES類RANS-LES混合方法構(gòu)造簡單、分區(qū)邊界無需人工干預(yù),對寬馬赫數(shù)范圍內(nèi)的多種類型復(fù)雜流動都能得到精細(xì)預(yù)測結(jié)果,具備預(yù)測大分離、剪切層失穩(wěn)、激波/邊界層干擾等典型非定常流動特征的能力,具有廣闊的應(yīng)用空間。目前DDES、IDDES等先進(jìn)DES類方法還沒有在工程應(yīng)用中大面積普及,還需要進(jìn)一步推廣。

    4) 以DES類方法為代表的RANS-LES混合方法仍然有較大的改進(jìn)空間。第二代URANS方法具備更強(qiáng)的非定常預(yù)測能力,有望逐步應(yīng)用于全機(jī)復(fù)雜流動的快速預(yù)測。對于飛機(jī)局部部件流動,植入式方法具備更高的物理解析能力和計算效率,在流動控制裝置設(shè)計、驗(yàn)證和機(jī)理研究上能發(fā)揮重要作用;自由來流湍流度對流場精細(xì)預(yù)測的影響也是需要加以考慮的內(nèi)容。在動態(tài)失速、燃燒、氣動彈性、氣動噪聲、氣動光學(xué)等與非定常流動密切耦合的方面,不論是機(jī)理分析還是工程應(yīng)用,RANS-LES混合方法也都有很廣闊的發(fā)展空間。

    致謝:本文部分內(nèi)容承蒙國家自然科學(xué)基金(11372159)、國家重點(diǎn)研發(fā)計劃(2016YFA0401200)和歐盟地平線2020研究創(chuàng)新項(xiàng)目IMAGE(688971)的資助。感謝清華大學(xué)信息科學(xué)與技術(shù)國家實(shí)驗(yàn)室提供計算資源。感謝黃靜波博士、肖良華博士、段志偉博士和甕哲工程師在UNITs軟件發(fā)展中的貢獻(xiàn)。

    [1]Jameson A.Re-engineering the design process through computation[J].Journal of Aircraft, 1999, 36(1): 36-50.

    [2]Kraft E M.Integrating computational science and engineering with testing to re-engineer the aeronautical development process.AIAA 2010-139[R].Orlando: AIAA, 2010.

    [3]YanC, Yu J, Xu J, et al.On the achievements and prospects for the methods of computational fluid dynamics[J].Advances in Mechanics, 2011, 41(5): 562-589.(in Chinese)閻超, 于劍, 徐晶磊, 等.CFD模擬方法的發(fā)展成就與展望[J].力學(xué)進(jìn)展, 2011, 41(5): 562-589.

    [4]Menter F R.Review of the shear-stress transport turbulence model experience from an industrial perspective[J].Interna-tional Journal of Computational Fluid Dynamics, 2009, 23(4): 305-316.

    [5]Zhang Y.Aerodynamic optimization of civil aircraft design based on advanced computational fluid dynamics[D].Beijing: Tsinghua University, 2010.(in Chinese)張宇飛.基于先進(jìn)CFD方法的民用客機(jī)氣動優(yōu)化設(shè)計[D].北京: 清華大學(xué), 2010.

    [6]Spalart P R, Allmaras S R.A one-equation turbulence model for aerodynamic flows.AIAA 92-0439[R].Reno: AIAA, 1992.

    [7]Menter F R.Two-equation eddy-viscosity turbulence models for engineering applications[J].AIAA Journal, 1994, 32(8): 1598-1605.

    [8]Smagorinsky J.General circulation experiments with the primitive equations[J].Monthly Weather Review, 1963, 91: 99-164.

    [9]Spalart P R.Strategies for turbulence modelling and simulations[J].International Journal of Heat and Fluid Flow, 2000, 21(3): 252-263.

    [10]Slotnick J, Khodadoust A, Alonso J, et al.CFD vision 2030 study: a path to revolutionary computational aerosciences.NASA/CR 2014-218178[R].Hampton: NASA, 2014.

    [11]Piomelli U, Balaras E.Wall-layer models for large-eddy simulations[J].Annual Review of Fluid Mechanics, 2002, 34: 349-374.

    [12]FanT C, Tian M, Edwards J R, et al.Validation of a hybrid Reynolds-averaged/large-eddy simulation method for simulating cavity flameholder configurations.AIAA 2001-2929[R].Anaheim: AIAA, 2001.

    [13]Speziale C G.Turbulence modeling for time-dependent RANS and VLES: a review[J].AIAA Journal, 1998, 36(2): 173-84.

    [14]Batten P, Goldberg U, Chakravarthy S.Interfacing statistical turbulence closures with large-eddy simulation[J].AIAA Journal, 2004, 42(3): 485-492.

    [15]Spalart P R.Detached-Eddy Simulation[J].Annual Review of Fluid Mechanics, 2009, 41(1): 181-202.

    [16]Deck S.Zonal-detached-eddy simulation of the flow around a high-lift configuration[J].AIAA Journal, 2005, 43(11): 2372-2384.

    [17]Deck S.Recent improvements in the Zonal Detached Eddy Simulation (ZDES) formulation[J].Theoretical and Computational Fluid Dynamics, 2012, 26(6): 523-550.

    [18]Delanghe C, Merci B, Dick E.Very large eddy simulation and RNG turbulence models.AIAA 2001-3041[R].Anaheim: AIAA, 2001.

    [19]Delanghe C, Merci B, Dick E.Hybrid RANS/LES modelling with an approximate renormalization group.I: model development[J].Journal of Turbulence, 2005, 6(13): 1-18.

    [20]Chen S Y, Xia Z H, Pei S Y, et al.Reynolds-stress-constrained large-eddy simulation of wall-bounded turbulent flows[J].Journal of Fluid Mechanics, 2012, 703: 1-28.

    [21]Girimaji S S.Partially-averaged Navier-Stokes model for turbulence: a Reynolds-averaged Navier-Stokes to direct numerical simulation bridging method[J].Journal of Applied Mechanics-Transactions of the ASME, 2006, 73(3): 413-421.

    [22]Menter F, Egorov Y.A scale adaptive simulation model using two-equation models.AIAA 2005-1095[R].Reno: AIAA, 2005.

    [23]Menter F R, Egorov Y.The scale-adaptive simulation method for unsteady turbulent flow predictions.Part 1: theory and model description[J].Flow, Turbulence and Combustion, 2010, 85(1): 113-138.

    [24]Du R F, Yan C, Luo D H.Property assessment of PANS method for numerical simulation of flow around tandem cylinders[J].Journal of Beijing University of Aeronautics and Astronautics, 2015, 41(8): 1374-1380.(in Chinese)杜若凡, 閻超, 羅大海.PANS方法在雙圓柱繞流數(shù)值模擬中的性能分析[J].北京航空航天大學(xué)學(xué)報, 2015, 41(8): 1374-1380.

    [25]Luo D H, Yan C, Wang X Y.Partially averaged Navier-Stokes method for simulation of supersonic flow over a ramped-cavity[J].Journal of Aerospace Power, 2015(09): 2167-2173.(in Chinese)羅大海, 閻超, 王小永.部分平均Navier-Stokes方法模擬超聲速斜面空腔流動[J].航空動力學(xué)報, 2015(09): 2167-2173.

    [26]Egorov Y, Menter F.Development and application of SST-SAS turbulence model in the DESIDER project[C]//Peng S, Haase W.Advances in hybrid RANS-LES modelling.Berlin: Springer Berlin Heidelberg, 2008, 261-270.

    [27]Zhao R, Xu J L, Yan C, et al.Scale-adaptive simulation of flow past wavy cylinders at a subcritical Reynolds number[J].Acta Mechanica Sinica, 2011, 27(5): 660-667.

    [28]Egorov Y, Menter F R, Lechner R, et al.The scale-adaptive simulation method for unsteady turbulent flow predictions.Part 2: application to complex flows[J].Flow, Turbulence and Combustion, 2010, 85(1): 139-165.

    [29]Du L, Ning F F.Scale adaptive simulation of flows around a circular cylinder at high sub-critical Reynolds number[J].Chinese Journal of Theoretical and Applied Mechanics, 2014, 4(46): 487-496.(in Chinese)杜磊, 寧方飛.高亞臨界雷諾數(shù)下圓柱繞流的尺度自適應(yīng)模擬[J].力學(xué)學(xué)報, 2014, 4(46): 487-496.

    [30]Babu S V, Zografakis G, Barakos G N.Evaluation of scale-adaptive simulations for transonic cavity flows[C]//Girimaji S, Haase W, Peng S, et al.Progress in hybrid RANS-LES Modelling.Switzerland: Springer International Publishing, 2015: 433-444.

    [31]Mockett C, Haase W, Thiele F.Go4Hybrid: A European initiative for improved hybrid RANS-LES modelling[C]//Girimaji S, Haase W, Peng S, et al.Progress in hybrid RANS-LES modelling.Switzerland: Springer International Publishing, 2015, 299-303.

    [32]Shur M, Spalart P R, Strelets M, et al.A rapid and accurate switch from RANS to LES in boundary layers using an overlap region[J].Flow, Turbulence and Combustion, 2011, 86(2): 179-206.

    [33]Shur M, Spalart P R, Strelets M, et al.Synthetic turbulence generators for RANS-LES interfaces in zonal simulations of aerodynamic and aeroacoustic problems[J].Flow, Turbulence and Combustion, 2014, 93(1): 63-92.

    [34]Chen H, Li Z, Zhang Y.A window-embedded RANS/LES method[C]//The Chinese Congress of Theoretical and Applied Mechanics.Shanghai: The Chinese Society of Theoretical and Applied Mechanics, 2015.(in Chinese)陳海昕, 李釗, 張宇飛.一種RANS/LES嵌入耦合方法[C]//中國力學(xué)大會, 上海: 中國力學(xué)學(xué)會, 2015.

    [35]Weng Z.Generation of synthetic turbulence for the embedded-IDDES model[D].Beijing: Tsinghua University, 2015.(in Chinses).甕哲.植入式混合方法之合成湍流生成方法研究[D].北京:清華大學(xué), 2015.

    [36]Li D, Wang X Y.Applying embedded hybrid RANS/LES method to channel flow[J].Journal of Northwestern Poly-technical University, 2015, 33(5): 799-803.(in Chinese)李棟, 王翔宇.嵌入式RANS/LES混合方法在槽道流動數(shù)值模擬中的應(yīng)用[J].西北工業(yè)大學(xué)學(xué)報, 2015, 33(5): 799-803.

    [37]Spalart P R, Jou W, Strelets M, et al.Comments on the feasibility of LES for wings, and on a hybrid RANS/LES approach[C]//Liu C, Liu Z.Advances in DNS/LES.Columbus: Greyden Press, 1997.

    [38]Wang X Y, Li D.Improved SST-DES in numerical simulation of mild separation[J].Journal of Beijing University of Aeronautics and Astronautics, 2014, 40(9): 1245-1249.(in Chinese)王翔宇, 李棟.SST-DES在小分離流動數(shù)值模擬中的改進(jìn)[J].北京航空航天大學(xué)學(xué)報, 2014, 40(9): 1245-1249.

    [39]Strelets M.Detached eddy simulation of massively separated flows.AIAA 2001-0879[R].Reno: AIAA, 2001.

    [40]Spalart P R, Deck S, Shur M L, et al.A new version of detached-eddy simulation, resistant to ambiguous grid densities[J].Theoretical and Computational Fluid Dynamics, 2006, 20(3): 181-195.

    [41]Menter F R, Kuntz M, Langtry R.Ten years of industrial experience with the SST turbulence model[C]//Hanjalic K, Nagano Y, Tummers M.Turbulence, heat and mass transfer 4.Antalya: Begell House, 2003.

    [42]Shur M L, Spalart P R, Strelets M K, et al.A hybrid RANS-LES approach with delayed-DES and wall-modelled LES capabilities[J].International Journal of Heat and Fluid Flow, 2008, 29(6): 1638-1649.

    [43]Xiao Z X, Liu J, Huang J B, et al.Comparisons of three improved DES methods on unsteady flows past tandem cylinders[C]//Fu S, Haase W, Peng S, et al.Progress in Hybrid RANS-LES modelling.Berlin: Springer-Verlag Berlin Heidelberg, 2012.

    [44]Gritskevich M S, Garbaruk A V, Schütze J, et al.Development of DDES and IDDES formulations for thek-ωshear stress transport model[J].Flow, Turbulence and Combustion, 2012, 88(3): 431-449.

    [45]Swanson R C, Turkel E.A multistage time-stepping scheme for the Navier-Stokes equations.AIAA 85-0035[R].Reno: AIAA, 1985.

    [46]Ramboer J, Broeckhoven T, Lacor C, et al.Optimized Runge-Kutta schemes for use with compact central and upwind schemes.AIAA 2005-4992[R].Toronto: AIAA, 2005.

    [47]Jameson A, Schmidt W, Turkel E.Numerical solutions of euler equations by finite volume methods with Runge-Kutta time stepping schemes.AIAA 81-1259[R].Palo Alto: AIAA, 1981.

    [48]Beam R, Warming R F.An implicit factored scheme for the compressible Navier-Stokes equations[J].AIAA Journal, 1978, 16(4): 393-402.

    [49]Pulliam T H, Chaussee D S.A diagonal form of an implicit approximate factorization algorithm[J].Journal of Computa-tional Physics, 1981, 39(2): 347-363.

    [50]Jameson A, Yoon S.Lower-upper implicit scheme with multiple grids for the Euler equations[J].AIAA Journal, 1987, 25(7): 929-935.

    [51]YoonS, Kwak D.Three-dimensional incompressible Navier-Stokes equations solver using lower-upper symmetric-Gauss-Siedel algorithm[J].AIAA Journal, 1991, 29(6): 874-875.

    [52]Candler G V, Wright M J.Data-parallel lower-upper relaxation method for reacting flows[J].AIAA Journal, 1994, 32(12): 2380-2386.

    [53]Wright M J, Candler G V, Prampolini M.Data-parallel lower-upper relaxation method for the Navier-Stokes equations[J].AIAA Journal, 1996, 34(7): 1371-1377.

    [54]Shen Y Q, Zha G C.Comparison of high order schemes for large eddy simulation of circular cylinder flow.AIAA 2009-945[R].Orlando: AIAA, 2009.

    [55]BuiT T.A parallel, finite-volume algorithm for large-eddy simulation of turbulent flow[J].Computers & Fluids, 2000, 29(8): 877-915.

    [56]Fu S, Xiao Z X, Chen H X, et al.Simulation of wing-body junction flows with hybrid RANS/LES methods[J].International Journal of Heat and Fluid Flow, 2007, 28(6): 1379-1390.

    [57]Roe P L.Approximate Riemann solver, parameter vectors, and difference schemes[J].Journal of Computational Physics, 1981, 43: 357-372.

    [58]Anderson W K, Thomas J L, Van Leer B.A comparison of finite volume flux vector splitting for the Euler equations.AIAA 85-0122[R].Reno: AIAA, 1985.

    [59]Moin P.Numerical and physical issues in large eddy simulation of turbulent flows[J].JSME International Journal Series B, 1998, 41(2): 454-463.

    [60]Yee H C, Sandham N D, Djomehri M J.Low dissipative high order shock-capturing methods using characteristic-based filters[J].Journal of Computational Physics, 1999, 150(1): 199-238.

    [61]Jiang G, Shu C W.Efficient implementation of weighted ENO schemes[J].Journal of Computational Physics, 1996, 126(1): 202-228.

    [62]Sun Z, Ren Y, Larricq C, et al.A class of finite difference schemes with low dispersion and controllable dissipation for DNS of compressible turbulence[J].Journal of Computational Physics, 2011, 12(230): 4616-4635.

    [63]Strelets M.Detached eddy simulation of massively separated flows.AIAA 2001-0879[R].Reno: AIAA, 2001.

    [64]Mockett C.A comprehensive study of detached-eddy simulation[D].Berlin: Technical University of Berlin, 2009.

    [65]Xiao Z X, Liu J, Huang J B, et al.Numerical dissipation effects on massive separation around tandem cylinders[J].AIAA Journal, 2012, 50(5): 1119-1136.

    [66]Zhang Y, Zhang L P, He X, et al.Detached-eddy simulation based on unstructured and hybrid grid[J].Acta Aeronautica et Astronautica Sinica, 2015, 36(9): 2900-2910.(in Chinese)張揚(yáng), 張來平, 赫新, 等.基于非結(jié)構(gòu)/混合網(wǎng)格的脫體渦模擬算法[J].航空學(xué)報, 2015, 36(9): 2900-2910.

    [67]Ducros F, Ferrand V, Nicoud F, et al.Large-eddy simulation of the shock/turbulence interaction[J].Journal of Computa-tional Physics, 1999, 152(2): 517-549.

    [68]Forsythe J R, Squires K D, Wurtzler K E, et al.Detached-eddy simulation of the F-15E at high alpha[J].Journal of Aircraft, 2004, 41(2): 193-200.

    [69]Fr?hlich J, Von Terzi D.Hybrid LES/RANS methods for the simulation of turbulent flows[J].Progress in Aerospace Sciences, 2008, 44(5): 349-377.

    [70]Brentner K S, Farassat F.Analytical comparison of the acoustic analogy and Kirchhoff formulation for moving surfaces[J].AIAA Journal, 1998, 36(8): 1379-1386.

    [71]Lockard D P, Choudhari M M, Khorrami M R, et al.Aeroacoustic simulations of tandem cylinders with subcritical spacing.AIAA 2008-2862[R].Vancouver: AIAA, 2008.

    [72]Lockard D P.Summary of the tandem cylinder solutions from the benchmark problems for airframe noise computations-I workshop.AIAA 2011-353[R].Orlando: AIAA, 2011.

    [73]Weinmann M, Sandberg R D, Doolan C.Tandem cylinder flow and noise predictions using a hybrid RANS/LES approach[J].International Journal of Heat and Fluid Flow, 2014, 50: 263-278.

    [74]Xiao Z X, Luo K Y.Improved delayed detached-eddy simulation of massive separation around triple cylinders[J].Acta Mechanica Sinica, 2015, 31(6): 799-816.

    [75]Hu N, Hao X, Su C, et al.Aeroacoustic study of landing gear by detached eddy simulation[J].Acta Aerodynamica Sinica, 2015, 33(1): 99-106.(in Chinese)胡寧, 郝璇, 蘇誠, 等.基于分離渦模擬的起落架氣動噪聲研究[J].空氣動力學(xué)學(xué)報, 2015, 33(1): 99-106.

    [76]Xiao Z X, Liu J, Luo K Y, et al.Investigation of flows around a rudimentary landing gear with advanced detached-eddy-simulation approaches[J].AIAA Journal, 2013, 51(1): 107-125.

    [77]Spalart P R, Shur M L, Strelets M K, et al.Sensitivity of landing-gear noise predictions by large-eddy simulation to numerics and resolution.AIAA 2012-1174[R].Nashville: AIAA, 2012.

    [78]Sanders L, Manoha E, Khelil S B, et al.LAGOON: CFD/CAA coupling for landing gear noise and comparison with experimental database.AIAA 2011-2822[R].Portland: AIAA, 2011.

    [79]Sanders L, Manoha E, Khelil S B, et al.LAGOON: New mach landing gear noise computation and further analysis of the CAA process.AIAA 2012-2281[R].Colorado Springs: AIAA, 2012.

    [80]Vuillot F, Lupoglazoff N, Luquent D, et al.Hybrid CAA solutions for nose landing gear noise.AIAA 2012-2283[R].Colorado Springs: AIAA, 2012.

    [81]Vatsa V N, Khorrami M R, Park M A, et al.Aeroacoustic simulation of nose landing gear on adaptive unstructured grids with FUN3D.AIAA 2013-2071[R].Berlin: AIAA, 2013.

    [82]Herr M, Ewert R, Rautmann C, et al.Broadband trailing-edge noise predictions: Overview of BANC-III results.AIAA 2015-2847[R].Dallas: AIAA, 2015.

    [83]Xu J M, Song W B.Zonal LES/DES method and its application in slat flow field simulation[J].Acta Aerodynamica Sinica.2014, 32(5): 668-674.(in Chinese)徐佳敏, 宋文濱.分區(qū)LES/DES混合方法及縫翼三維流場模擬[J].空氣動力學(xué)學(xué)報, 2014, 32(5): 668-674.

    [84]Xia M, Zhang X D, Zhong B W.Numerical simulation of separated flow over Gurney flap based on DDES[J].Testing and Simulation, 2014, 16: 114-117.(in Chinese)夏明, 張曉東, 鐘伯文.基于延遲脫體渦模擬方法對Gurney襟翼分離流動的數(shù)值模擬[J].航空制造技術(shù), 2014, 16: 114-117.

    [85]Choudhari M M, Lockard D P.Assessment of slat noise predictions for 30P30N high-lift configuration from BANC-III workshop.AIAA 2015-2844[R].Dallas: AIAA, 2015.

    [86]Escobar J A, Suarez C A, Silva C, et al.Detached-eddy simulation of a wide-body commercial aircraft in high-lift configuration[J].Journal of Aircraft, 2015, 52(4): 1112-1121.

    [87]Chen R Q, Wans X, You Y C.Numerical investigation of nacelle′s effects on propeller slipstream based on IDDES model[J].Acta Aeronautica et Astronautica Sinica, 2016: 1-11.(in Chinese)陳榮錢, 王旭, 尤延鋮.短艙對螺旋槳滑流影響的IDDES數(shù)值模擬研究[J].航空學(xué)報, 2016: 1-11.

    [88]Wu J F, Ning F F.Hybrid RANS-LES method applied to backward facing step flow[J].Journal of Beijing University of Aeronautics and Astronautics, 2011(06): 701-704.(in Chinese)吳晶峰, 寧方飛.后臺階流動的Hybrid RANS/LES模擬[J].北京航空航天大學(xué)學(xué)報, 2011(06): 701-704.

    [89]Hu R Y, Wang L, Fu S.Improved delayed detached eddy simulation of flow structures behind a backward-facing step.AIAA 2016-1104[R].San Diego: AIAA, 2016.

    [90]Liu Z, Yang Y J, Zhou W J, et al.Study of unsteady separation flow around airfoil at high angle of attack using hybrid RANS-LES method[J].Acta Aeronautica et Astronautica Sinica, 2014, 35(2): 372-380.(in Chinese)劉周, 楊云軍, 周偉江, 等.基于RANS-LES混合方法的翼型大迎角非定常分離流動研究[J].航空學(xué)報, 2014, 35(2): 372-380.

    [91]Li Q, Sun D, Zhang H X.Detached-eddy simulations and analyses on new vortical flows over a 76°/40° double delta wing[J].Science China Physics, Mechanics and Astronomy, 2013, 56(6): 1062-1073.

    [92]Liu J, Sun H S, Liu Z T, et al.Numerical investigation of unsteady vortex breakdown past 80°/65° double-delta wing[J].Chinese Journal of Aeronautics, 2014, 27(3): 521-530.

    [93]Gaitonde D V.Progress in shock wave/boundary layer interactions[J].Progress in Aerospace Sciences.2015, 72: 80-99.

    [94]Zhang W W, Gao C Q, Ye Z Y.Research advances of win/airfoil transonic buffet[J].Acta Aeronautica et Astronautica Sinica, 2015(04): 1056-1075.(in Chinese)張偉偉, 高傳強(qiáng), 葉正寅.機(jī)翼跨聲速抖振研究進(jìn)展[J].航空學(xué)報, 2015(04): 1056-1075.

    [95]Huang J B, Xiao Z X, Liu J, et al.Simulation of shock wave buffet and its suppression on an OAT15A supercritical airfoil by IDDES[J].Science China-Physics Mechanics & Astronomy, 2012, 55(2): 260-271.

    [96]Lawson S J, Barakos G N.Review of numerical simulations for high-speed, turbulent cavity flows[J].Progress in Aerospace Sciences, 2011, 47(3): 186-216.

    [97]Henshaw M J.M219 cavity case.ADPO10729[R].Yorkshire: Defense Technical Information Center, 2000.

    [98]Haase W, Braza M, Revell A.DESider-A European effort on hybrid RANS-LES modelling[M].Springer Berlin Heidelberg, 2009.

    [99]Barakos G N, Lawson S J, Steijl R, et al.Numerical simulations of high-speed turbulent cavity flows[J].Flow, Turbulence and Combustion, 2009, 83(4): 569-585.

    [100]Lawson S J, Barakos G N.Computational fluid dynamics analyses of flow over weapons-bay geometries[J].Journal of Aircraft, 2010, 47(5): 1605-1623.

    [101]Liggett N D, Smith M J.Cavity flow assessment using advanced turbulence methods[J].Journal of Aircraft, 2011, 48(1): 141-156.

    [102]Temmerman L, Tartinville B, Hirsch C.URANS investigation of the transonic M219 cavity[C]//Fu S, Haase W, Peng S, et al.Progress in hybrid RANS-LES modelling.Germany: Springer Berlin Heidelberg, 2012, 471-481.

    [103]BabuS V, Zografakis G, Barakos G N.Evaluation of scale-adaptive simulations for transonic cavity flows[C]//Girimaji S, Haase W, Peng S, et al.Progress in hybrid RANS-LES modelling.Switzerland: Springer International Publishing, 2015: 433-444.

    [104]Luo K Y, Xiao Z X.Improved delayed detached-eddy simulation of transonic and supersonic cavity flows[C]//Girimaji S, Haase W, Peng S, et al.Progress in hybrid RANS-LES modelling.Switzerland: Springer International Publish-ing, 2015: 163-174.

    [105]Liu Y, Tong M B, Hu Z W.DDES of aeroacoustic over an open cavity with and without a spoiler[J].Acta Aerodynamica Sinica, 2015, 33(5): 643-648.(in Chinese)劉瑜, 童明波, Hu Z W.基于DDES算法的有擾流片腔體氣動噪聲分析[J].空氣動力學(xué)學(xué)報, 2015, 33(5): 643-648.

    [106]Liu J,Yang D G, Wang X S, et al.Effect of the turbulent boundary layer thickness on a three-dimensional cavity flow[J].Acta Aeronautica et Astronautica Sinica, 2016, 37(2): 475-483.劉俊, 楊黨國, 王顯圣, 等.湍流邊界層厚度對三維空腔流動的影響研究[J].航空學(xué)報, 2016, 37(2): 475-483.

    [107]Lawson S J, Barakos G N.Evaluation of DES for weapons bays in UCAVs[J].Aerospace Science and Technology, 2010, 14(6): 397-414.

    [108]Kannepalli C, Chartrand C, Birkbeck R, et al.Computational modeling of geometrically complex weapons bays.AIAA 2011-2774[R].Portland: AIAA , 2011.

    [109]Kim D H, Choi J H, Kwon O J.Detached eddy simulation of weapons bay flows and store separation[J].Computers & Fluids, 2015, 121: 1-10.

    [110]Luo K Y, Xiao Z X.Application of the improved delayed detached-eddy simulation to transonic and supersonic cavity flows[C]// 3rd Symposium on Fluid-Structure-Sound Inter-actions and Control.Perth: Curtin University, 2015.

    [111]Weng Z, Luo K Y, Fu S, et al.Improved-delayed-detached-eddy simulations of saw-tooth spoilers control before a supersonic weapon bay[C]//Turbulence, Heat and Mass Transfer 8, Sarajevo: Begell House, 2015.

    [112]Pain R, Weiss P E, Deck S.Zonal detached eddy simulation of the flow around a simplified launcher afterbody[J].AIAA Journal, 2014, 52(9): 1967-1979.

    [113]Heinrich L, Daniel M J, Klaus H.Launch vehicle base flow analysis using improved delayed detached-eddy simulation[J].AIAA Journal, 2015, 53(9): 2454-2471.

    [114]Brock J M, Subbareddy P K, Candler G V.Detached-eddy simulations of hypersonic capsule wake flow[J].AIAA Journal, 2015, 53(1): 70-80.

    [115]Salazar G, Edwards J R.Mach 6 wake flow simulations using a large-eddy simulation/reynolds-averaged Navier-Stokes model[J].Journal of Spacecraft and Rockets, 2014, 51(4): 1329-1348.

    [116]Liu J, Yan C, Zhao R, et al.Simulation of re-entry capsule

    thermodynamics environment by DES method[J].Journal of Beijing University of Aeronautics and Astronautics, 2013, 39(5): 590-594.(in Chinese)劉佳, 閻超, 趙瑞, 等.DES方法對返回艙氣動熱環(huán)境數(shù)值模擬[J].北京航空航天大學(xué)學(xué)報, 2013(05): 590-594.

    [117]Dong X R, Chen Y H, Dong G, et al.DES numerical study of shock wave/boundary layer interactions in hypersonic flows controlled by double micro-ramps[J].Acta Aeronautica et Astronautica Sinica, 2016, 37(6): 1771-1780.(in Chinese) 董祥瑞, 陳耀慧, 董剛, 等.基于DES方法的高超聲速激波/邊界層干擾的雙微楔控制數(shù)值研究[J].航空學(xué)報, 2016, 37(6): 1771-1780.

    [118]Barnhardt M, Candler G V.Detached-eddy simulation of the Reentry-F flight experiment[J].Journal of Spacecraft and Rockets, 2012, 49(4): 691-699.

    [119]Duan Z W, Xiao Z X, Fu S.Simulation of transition triggered by isolated roughness in hypersonic boundary layer.AIAA 2012-3076[R].New Orleans: AIAA, 2012.

    [120]Xiao L H, Xiao Z X, Duan Z W, et al.Improved-delayed-detached-eddy simulation of cavity-induced transition in hypersonic boundary layer[J].International Journal of Heat and Fluid Flow, 2015, 51: 138-150.

    [121]Chauvet N, Deck S, Jacquin L.Zonal detached eddy simulation of a controlled propulsive jet[J].AIAA Journal, 2007, 45(10): 2458-2473.

    [122]Won S, Jeung I, Parent B, et al.Numerical investigation of transverse hydrogen jet into supersonic crossflow using detached-eddy simulation[J].AIAA Journal, 2010, 48(6): 1047-1058.

    [123]Dong H, Wen X Y, Li P, et al.DES of spray combustion in can-annular combustor[J].Journal of Naval University of Engineering, 2014, (5): 57-61.(in Chinese)董紅, 聞雪友, 李鵬, 等.環(huán)管型燃燒室噴霧燃燒DES模擬[J].海軍工程大學(xué)學(xué)報, 2014, (5): 57-61.

    [124]Wang H S, Shan F L, Niu J, et al.IDDES simulation of supersonic combustion based on flamelet/progress variable model[J].Journal of Propulsion Technology, 2015, 36(3): 321-327.(in Chinese)王浩蘇, 單繁立, 牛健, 等.基于火焰面/進(jìn)度變量模型的超聲速燃燒IDDES模擬[J].推進(jìn)技術(shù), 2015, 36(3): 321-327.

    [125]Hee G J, Chwalowski P, SchusterD M, et al.Plans and example results for the 2ndAIAA aeroelastic prediction workshop.AIAA 2015-0437[R].Kissimmee: AIAA, 2015.

    [126]Spalart P R, Shur M L, Strelets M K, et al.Sensitivity of landing-gear noise predictions by large-eddy simulation to numerics and resolution.AIAA 2012-1174[R].Nashville: AIAA, 2012.

    [127]Wang M, Mani A, Gordeyev S.Physics and computation of aero-optics[J].Annual Review of Fluid Mechanics, 2012, 44: 299-321.

    Developments and applications of hybrid RANS-LES methods for wide-speed-range flows

    XIAO Zhixiang*, LUO Kunyu, LIU Jian

    (SchoolofAerospaceEngineering,TsinghuaUniversity,Beijing100086,China)

    The growing demands of aerospace industry require accurate prediction of unsteady flow details.Current Reynolds averaged Navier-Stokes (RANS) methods are unable to provide dynamic loads for unsteady turbulent flows at high Reynolds numbers, such as massively separated flows past complex geometries.Large eddy simulation (LES) and direct numerical simulation (DNS) are still too expensive for engineering applications.Hybrid RANS-LES methods, combining near-wall RANS regions and outer LES regions, are the most promising techniques for unsteady turbulent flows in engineering.First of all, a general introduction to several categories of hybrid RANS-LES methods was given to discuss their basic ideas and characteristics.Then, the development history of detached-eddy simulation (DES) type methods was presented together with the influences of high-accuracy time-marching methods, spatial discretization schemes and proper numerical dissipation.For both mechanism studies and engineering applications, multiple cases at Mach numbers varying from 0.1 to 20 show that hybrid RANS-LES methods are an ideal choice for the prediction of turbulent flow in engineering.Hybrid RANS-LES methods are capable of predicting complex features involved in massive separation.Future improvement includes the promotion in computation efficiency and embedded RANS/LES strategy for detailed flow past small local component.Big potential also exists for

    hybrid RANS-LES methods in areas closely related to unsteady flow, such as dynamic stall, combustion and aero-elastics/acoustics/optics.

    hybrid RANS-LES methods; detached-eddy simulation; computational fluid dynamics; unsteady flow; separation; turbulence at high reynolds number

    0258-1825(2017)03-0338-16

    2017-03-22;

    2017-04-10

    國家自然科學(xué)基金面上項(xiàng)目(11372159); 國家重點(diǎn)研發(fā)計劃(2016YFA0401200); 歐盟地平線2020研究創(chuàng)新項(xiàng)目IMAGE(688971)

    肖志祥*(1974-),男,四川資陽人,特別研究員,博士,研究方向:RANS-LES混合方法、高精度湍流預(yù)測、可壓縮轉(zhuǎn)捩/湍流模式、高超聲速氣動力和氣動熱、計算氣動聲學(xué).E-mail:xiaotigerzhx@tsinghua.edu.cn

    肖志祥, 羅堃宇, 劉健.寬速域RANS-LES混合方法的發(fā)展及應(yīng)用[J].空氣動力學(xué)學(xué)報, 2017, 35(3): 338-353.

    10.7638/kqdlxxb-2017.0048 XIAO Z X, LUO K Y, LIU J.Developments and applications of hybrid RANS-LES methods for wide-speed-range flows[J].Acta Aerodynamica Sinica, 2017, 35(3): 338-353.

    V211.3

    A doi: 10.7638/kqdlxxb-2017.0048

    猜你喜歡
    邊界層湍流尺度
    財產(chǎn)的五大尺度和五重應(yīng)對
    基于HIFiRE-2超燃發(fā)動機(jī)內(nèi)流道的激波邊界層干擾分析
    重氣瞬時泄漏擴(kuò)散的湍流模型驗(yàn)證
    宇宙的尺度
    太空探索(2016年5期)2016-07-12 15:17:55
    一類具有邊界層性質(zhì)的二次奇攝動邊值問題
    非特征邊界的MHD方程的邊界層
    9
    “青春期”湍流中的智慧引渡(三)
    “青春期”湍流中的智慧引渡(二)
    弱分層湍流輸運(yùn)特性的統(tǒng)計分析
    日本wwww免费看| 国产精品国产三级国产专区5o| 日韩欧美 国产精品| 成人亚洲欧美一区二区av| 亚洲精华国产精华液的使用体验| 晚上一个人看的免费电影| 一级av片app| 日本黄色片子视频| 免费看不卡的av| 成人亚洲精品一区在线观看 | 亚洲av成人精品一区久久| 国产精品麻豆人妻色哟哟久久| 欧美+日韩+精品| 久久鲁丝午夜福利片| 亚洲精品久久久久久婷婷小说| 一本色道久久久久久精品综合| 国产免费视频播放在线视频| 欧美成人午夜免费资源| 三级经典国产精品| 少妇人妻 视频| av在线天堂中文字幕| 最后的刺客免费高清国语| 在线亚洲精品国产二区图片欧美 | 91在线精品国自产拍蜜月| 亚洲精品成人av观看孕妇| 黄色欧美视频在线观看| 国产高清不卡午夜福利| 免费av毛片视频| 久久99热这里只频精品6学生| 麻豆乱淫一区二区| 高清视频免费观看一区二区| 国产精品无大码| 国产精品爽爽va在线观看网站| 热99国产精品久久久久久7| 亚洲图色成人| 伊人久久精品亚洲午夜| 美女主播在线视频| 日日摸夜夜添夜夜爱| 亚洲人与动物交配视频| 干丝袜人妻中文字幕| 黄色视频在线播放观看不卡| 少妇熟女欧美另类| 超碰97精品在线观看| 在线观看人妻少妇| 最近最新中文字幕大全电影3| 国产精品99久久久久久久久| 网址你懂的国产日韩在线| 一个人看的www免费观看视频| 天天一区二区日本电影三级| 成人高潮视频无遮挡免费网站| 免费观看av网站的网址| 永久网站在线| 三级男女做爰猛烈吃奶摸视频| 大陆偷拍与自拍| 91久久精品电影网| 欧美亚洲 丝袜 人妻 在线| 国产精品福利在线免费观看| 五月伊人婷婷丁香| 国产综合精华液| 亚洲欧美一区二区三区黑人 | 免费av不卡在线播放| 久久鲁丝午夜福利片| 热99国产精品久久久久久7| 国产免费视频播放在线视频| 国模一区二区三区四区视频| 少妇人妻精品综合一区二区| 国产综合懂色| 一个人看的www免费观看视频| 麻豆久久精品国产亚洲av| 99视频精品全部免费 在线| 99视频精品全部免费 在线| 国产有黄有色有爽视频| av网站免费在线观看视频| av.在线天堂| 九九在线视频观看精品| 寂寞人妻少妇视频99o| 岛国毛片在线播放| 黄色日韩在线| 国产精品一区二区性色av| 免费大片黄手机在线观看| 99re6热这里在线精品视频| 久久鲁丝午夜福利片| 国产有黄有色有爽视频| 制服丝袜香蕉在线| 九九在线视频观看精品| 中文资源天堂在线| 亚洲在线观看片| 亚洲三级黄色毛片| 又黄又爽又刺激的免费视频.| 国产亚洲av片在线观看秒播厂| 亚洲精品乱码久久久v下载方式| 中文字幕亚洲精品专区| 高清在线视频一区二区三区| 水蜜桃什么品种好| 亚洲欧美中文字幕日韩二区| 亚洲精品成人av观看孕妇| 三级国产精品片| 亚洲人成网站在线播| 欧美日韩亚洲高清精品| 精品人妻熟女av久视频| 国产男女内射视频| 国产男人的电影天堂91| 成人黄色视频免费在线看| 亚洲综合色惰| 久久久久久九九精品二区国产| 国产成人精品福利久久| 在线精品无人区一区二区三 | 精品一区二区三区视频在线| 亚洲国产最新在线播放| 2021少妇久久久久久久久久久| 狂野欧美白嫩少妇大欣赏| 啦啦啦啦在线视频资源| 女的被弄到高潮叫床怎么办| 别揉我奶头 嗯啊视频| 青春草视频在线免费观看| 亚洲精品456在线播放app| 青春草国产在线视频| 人妻少妇偷人精品九色| 久久6这里有精品| 日本av手机在线免费观看| 中文精品一卡2卡3卡4更新| 欧美国产精品一级二级三级 | 国产熟女欧美一区二区| 男的添女的下面高潮视频| 日韩欧美精品v在线| 欧美bdsm另类| 日韩大片免费观看网站| 国产成年人精品一区二区| 国产国拍精品亚洲av在线观看| 视频中文字幕在线观看| 一区二区三区免费毛片| 久久久久性生活片| 国产久久久一区二区三区| 久久久a久久爽久久v久久| 青春草视频在线免费观看| 熟女av电影| h日本视频在线播放| 人人妻人人澡人人爽人人夜夜| 成人高潮视频无遮挡免费网站| 三级国产精品欧美在线观看| 亚洲精品自拍成人| 亚洲av国产av综合av卡| 国产乱人视频| 国产黄频视频在线观看| 少妇人妻 视频| 午夜激情久久久久久久| 卡戴珊不雅视频在线播放| 秋霞在线观看毛片| 97在线人人人人妻| 亚洲精华国产精华液的使用体验| 国产黄a三级三级三级人| 一区二区三区免费毛片| 亚洲国产色片| 国产精品99久久99久久久不卡 | 天堂中文最新版在线下载 | 欧美性猛交╳xxx乱大交人| 噜噜噜噜噜久久久久久91| 国产综合精华液| 王馨瑶露胸无遮挡在线观看| 最近2019中文字幕mv第一页| 在线免费观看不下载黄p国产| 日日摸夜夜添夜夜爱| 热99国产精品久久久久久7| 中文字幕制服av| 中国美白少妇内射xxxbb| 亚洲第一区二区三区不卡| 午夜福利在线在线| 国产精品人妻久久久影院| 久久久精品欧美日韩精品| 天天躁日日操中文字幕| 晚上一个人看的免费电影| 国产精品不卡视频一区二区| 精品熟女少妇av免费看| 欧美日韩在线观看h| 一区二区三区精品91| 国产一区亚洲一区在线观看| 少妇裸体淫交视频免费看高清| 久久久久九九精品影院| 久久精品国产亚洲av天美| 99久久精品一区二区三区| 国产成人精品福利久久| 三级国产精品欧美在线观看| 午夜福利在线观看免费完整高清在| 2021天堂中文幕一二区在线观| 少妇人妻一区二区三区视频| 国产国拍精品亚洲av在线观看| 欧美zozozo另类| 青青草视频在线视频观看| 狂野欧美激情性bbbbbb| 欧美极品一区二区三区四区| 国产午夜精品一二区理论片| 国产免费福利视频在线观看| 日韩大片免费观看网站| 六月丁香七月| 亚洲精品国产av蜜桃| 免费看光身美女| av国产免费在线观看| 视频中文字幕在线观看| 中文乱码字字幕精品一区二区三区| 99热全是精品| 男人添女人高潮全过程视频| 国产av码专区亚洲av| 中文字幕av成人在线电影| 国产av不卡久久| 欧美一区二区亚洲| 91在线精品国自产拍蜜月| 简卡轻食公司| 免费少妇av软件| 成人美女网站在线观看视频| 国产精品一二三区在线看| 麻豆久久精品国产亚洲av| 大片免费播放器 马上看| 国产毛片a区久久久久| 国产精品久久久久久久电影| 成人黄色视频免费在线看| 嫩草影院精品99| 国语对白做爰xxxⅹ性视频网站| 国产视频内射| 精品一区二区三卡| 99精国产麻豆久久婷婷| 国产精品秋霞免费鲁丝片| 在线观看国产h片| 日韩亚洲欧美综合| 美女视频免费永久观看网站| 日韩一区二区视频免费看| 国产高清不卡午夜福利| 精华霜和精华液先用哪个| 深爱激情五月婷婷| 亚洲va在线va天堂va国产| 噜噜噜噜噜久久久久久91| www.色视频.com| 中国国产av一级| 日韩av不卡免费在线播放| 日日啪夜夜撸| 视频区图区小说| 一级二级三级毛片免费看| 日本一二三区视频观看| 国产爽快片一区二区三区| 久热这里只有精品99| 三级男女做爰猛烈吃奶摸视频| 五月伊人婷婷丁香| 国产精品精品国产色婷婷| 王馨瑶露胸无遮挡在线观看| 国产久久久一区二区三区| 国产精品一区二区三区四区免费观看| 国产精品福利在线免费观看| 91久久精品国产一区二区三区| 精品人妻熟女av久视频| 又黄又爽又刺激的免费视频.| 熟女电影av网| 国产成人91sexporn| 亚洲国产精品成人综合色| 女人被狂操c到高潮| 亚洲,欧美,日韩| 亚洲av中文av极速乱| 国产国拍精品亚洲av在线观看| 精品人妻熟女av久视频| 亚洲最大成人手机在线| 亚洲最大成人手机在线| 国产成人a区在线观看| 欧美区成人在线视频| 哪个播放器可以免费观看大片| 一边亲一边摸免费视频| 18禁裸乳无遮挡免费网站照片| 特大巨黑吊av在线直播| 99久久精品一区二区三区| 亚洲av日韩在线播放| 我的女老师完整版在线观看| 一级毛片久久久久久久久女| 久久精品国产亚洲网站| 国产成人午夜福利电影在线观看| 国内精品美女久久久久久| 亚州av有码| 美女cb高潮喷水在线观看| 涩涩av久久男人的天堂| 人体艺术视频欧美日本| 亚州av有码| 在现免费观看毛片| 亚洲在线观看片| 制服丝袜香蕉在线| 在线观看三级黄色| 亚洲av免费高清在线观看| 中文在线观看免费www的网站| 高清午夜精品一区二区三区| 中国三级夫妇交换| 亚洲综合精品二区| 久久精品国产鲁丝片午夜精品| 成人免费观看视频高清| 91aial.com中文字幕在线观看| 美女xxoo啪啪120秒动态图| 麻豆成人午夜福利视频| 欧美潮喷喷水| 亚洲av国产av综合av卡| 美女xxoo啪啪120秒动态图| 香蕉精品网在线| 少妇的逼水好多| 日本午夜av视频| 日韩免费高清中文字幕av| 久久精品熟女亚洲av麻豆精品| av一本久久久久| 性色av一级| 91精品一卡2卡3卡4卡| 最新中文字幕久久久久| 国产亚洲av嫩草精品影院| 大又大粗又爽又黄少妇毛片口| 久久99蜜桃精品久久| 欧美zozozo另类| 国产综合精华液| 亚洲国产成人一精品久久久| 九九在线视频观看精品| 日本猛色少妇xxxxx猛交久久| 国产成人aa在线观看| 18禁动态无遮挡网站| 久久久精品免费免费高清| 男的添女的下面高潮视频| av又黄又爽大尺度在线免费看| av国产精品久久久久影院| 国产精品熟女久久久久浪| 久久6这里有精品| 亚洲欧美中文字幕日韩二区| 亚洲av成人精品一区久久| 99精国产麻豆久久婷婷| 国产精品人妻久久久影院| 成年人午夜在线观看视频| 亚洲国产精品成人久久小说| 一本一本综合久久| 99re6热这里在线精品视频| 国产男人的电影天堂91| 麻豆精品久久久久久蜜桃| 成年女人看的毛片在线观看| 大码成人一级视频| 乱系列少妇在线播放| 网址你懂的国产日韩在线| 精品久久久久久久久av| 嘟嘟电影网在线观看| av天堂中文字幕网| 一本色道久久久久久精品综合| 国产成人免费无遮挡视频| 国产精品人妻久久久久久| 亚洲成人精品中文字幕电影| 丝袜美腿在线中文| 国产人妻一区二区三区在| 精品酒店卫生间| 2018国产大陆天天弄谢| 女人被狂操c到高潮| 三级经典国产精品| 永久免费av网站大全| 在线观看国产h片| 国模一区二区三区四区视频| 国产亚洲av片在线观看秒播厂| 成人毛片a级毛片在线播放| 777米奇影视久久| 男插女下体视频免费在线播放| 高清视频免费观看一区二区| 午夜福利在线在线| 亚洲aⅴ乱码一区二区在线播放| 久久久久久久国产电影| 午夜福利在线观看免费完整高清在| 日韩三级伦理在线观看| 国产精品久久久久久精品电影小说 | 午夜免费鲁丝| 建设人人有责人人尽责人人享有的 | 久久韩国三级中文字幕| 高清在线视频一区二区三区| 亚洲精品久久久久久婷婷小说| 亚洲国产精品成人综合色| 日韩电影二区| 国内揄拍国产精品人妻在线| 六月丁香七月| 一本久久精品| 亚洲自拍偷在线| 精品一区二区免费观看| 亚洲成人一二三区av| 日韩一本色道免费dvd| 好男人视频免费观看在线| 赤兔流量卡办理| 免费av毛片视频| 三级国产精品片| 最近中文字幕2019免费版| 国产精品精品国产色婷婷| 日本一本二区三区精品| av国产免费在线观看| 日韩免费高清中文字幕av| 一个人看视频在线观看www免费| 色哟哟·www| 丝袜喷水一区| 1000部很黄的大片| 午夜亚洲福利在线播放| 老女人水多毛片| 国产黄色免费在线视频| 男插女下体视频免费在线播放| 亚洲,欧美,日韩| 97在线视频观看| 99热这里只有是精品50| 久久99热这里只有精品18| 韩国高清视频一区二区三区| 国产高潮美女av| av国产精品久久久久影院| 赤兔流量卡办理| 国产国拍精品亚洲av在线观看| 国产乱人偷精品视频| 91aial.com中文字幕在线观看| 久久精品国产亚洲网站| 亚洲经典国产精华液单| 校园人妻丝袜中文字幕| tube8黄色片| 高清日韩中文字幕在线| 亚洲美女搞黄在线观看| 亚洲欧美成人精品一区二区| 日本色播在线视频| 日本一二三区视频观看| 最近中文字幕2019免费版| 毛片一级片免费看久久久久| 国产欧美亚洲国产| 亚洲精品国产色婷婷电影| 国产成人a∨麻豆精品| 国产探花极品一区二区| 五月伊人婷婷丁香| 热re99久久精品国产66热6| 日本黄大片高清| 亚洲,一卡二卡三卡| 亚洲国产精品999| 亚洲国产欧美人成| 色视频www国产| 美女脱内裤让男人舔精品视频| 深爱激情五月婷婷| 欧美老熟妇乱子伦牲交| 精品久久国产蜜桃| tube8黄色片| 一级毛片 在线播放| 亚洲图色成人| 亚洲欧美精品专区久久| 18禁裸乳无遮挡免费网站照片| 99久久九九国产精品国产免费| 午夜亚洲福利在线播放| 亚洲国产精品999| 舔av片在线| 男女无遮挡免费网站观看| 少妇人妻久久综合中文| 国产欧美日韩一区二区三区在线 | 日日啪夜夜爽| 国内揄拍国产精品人妻在线| 婷婷色av中文字幕| 18禁裸乳无遮挡免费网站照片| 熟女人妻精品中文字幕| 精品国产露脸久久av麻豆| 亚洲国产精品成人综合色| 天天一区二区日本电影三级| 国产精品爽爽va在线观看网站| 在线亚洲精品国产二区图片欧美 | 国产白丝娇喘喷水9色精品| 成人国产av品久久久| 黄色日韩在线| 网址你懂的国产日韩在线| 国产精品无大码| 三级国产精品片| 亚洲av日韩在线播放| 欧美激情在线99| 国产色爽女视频免费观看| 26uuu在线亚洲综合色| 欧美极品一区二区三区四区| 97人妻精品一区二区三区麻豆| 一级黄片播放器| 日韩欧美 国产精品| 18禁在线无遮挡免费观看视频| 日韩一区二区视频免费看| 97超视频在线观看视频| 国产美女午夜福利| 精品熟女少妇av免费看| 日日摸夜夜添夜夜爱| 国产伦理片在线播放av一区| 国产视频首页在线观看| 亚洲av中文字字幕乱码综合| 国产亚洲午夜精品一区二区久久 | 亚洲欧美精品专区久久| 欧美3d第一页| 久久99蜜桃精品久久| 色视频www国产| 亚洲一区二区三区欧美精品 | 久久热精品热| 男人爽女人下面视频在线观看| 全区人妻精品视频| 国产精品久久久久久av不卡| 18+在线观看网站| 黄片wwwwww| 成年av动漫网址| 成人国产麻豆网| 亚洲av不卡在线观看| 国产午夜精品久久久久久一区二区三区| 在线免费十八禁| 精品熟女少妇av免费看| 亚洲精品影视一区二区三区av| 天天一区二区日本电影三级| 偷拍熟女少妇极品色| 天天一区二区日本电影三级| 丝袜美腿在线中文| 日日摸夜夜添夜夜添av毛片| 国产成人福利小说| 99热这里只有是精品50| 亚洲av.av天堂| 成人鲁丝片一二三区免费| 久久久午夜欧美精品| 嘟嘟电影网在线观看| 国内精品美女久久久久久| 国产av国产精品国产| 人体艺术视频欧美日本| 最近中文字幕高清免费大全6| 十八禁网站网址无遮挡 | 日本黄色片子视频| 国产精品无大码| 日韩亚洲欧美综合| 在线天堂最新版资源| 老司机影院成人| 1000部很黄的大片| 日韩免费高清中文字幕av| 欧美xxⅹ黑人| 春色校园在线视频观看| 亚洲欧美成人精品一区二区| 美女xxoo啪啪120秒动态图| 亚洲欧美日韩东京热| 久久久久久久精品精品| 97超碰精品成人国产| 国产成年人精品一区二区| 国产一区有黄有色的免费视频| 国产一区亚洲一区在线观看| 小蜜桃在线观看免费完整版高清| 久久久久久国产a免费观看| 青春草视频在线免费观看| 国产国拍精品亚洲av在线观看| 免费人成在线观看视频色| 一级片'在线观看视频| 国产探花在线观看一区二区| 久久久国产一区二区| 日本三级黄在线观看| 秋霞伦理黄片| 我的女老师完整版在线观看| 精品人妻一区二区三区麻豆| 亚洲av福利一区| 成人综合一区亚洲| 国产午夜精品久久久久久一区二区三区| 超碰97精品在线观看| 成人国产麻豆网| 国产成人a区在线观看| 又粗又硬又长又爽又黄的视频| 青春草视频在线免费观看| av在线天堂中文字幕| 国产淫语在线视频| a级毛片免费高清观看在线播放| 97超碰精品成人国产| 欧美zozozo另类| 精品熟女少妇av免费看| 男人爽女人下面视频在线观看| 久久99精品国语久久久| 狠狠精品人妻久久久久久综合| 亚洲国产av新网站| 欧美97在线视频| 一级毛片电影观看| 欧美+日韩+精品| 久久99热这里只频精品6学生| 国产精品一区二区在线观看99| 久久久精品94久久精品| 高清视频免费观看一区二区| 精品少妇久久久久久888优播| 亚洲av二区三区四区| 最近手机中文字幕大全| 亚洲av欧美aⅴ国产| 午夜激情福利司机影院| 2022亚洲国产成人精品| 久久99热6这里只有精品| 国产极品天堂在线| 建设人人有责人人尽责人人享有的 | 国产色婷婷99| 亚洲怡红院男人天堂| 青春草国产在线视频| 男人舔奶头视频| 99热6这里只有精品| 亚洲精品久久久久久婷婷小说| 韩国av在线不卡| 久久精品国产自在天天线| 99热这里只有是精品50| 中文精品一卡2卡3卡4更新| 日本一本二区三区精品| 国产亚洲91精品色在线| 久久精品熟女亚洲av麻豆精品| 丰满乱子伦码专区| 欧美老熟妇乱子伦牲交| 亚洲婷婷狠狠爱综合网| 国产乱人视频| 少妇人妻精品综合一区二区| 国产精品国产三级专区第一集| 性色av一级| 欧美bdsm另类| 熟妇人妻不卡中文字幕| 国产成年人精品一区二区| 97超碰精品成人国产| 好男人视频免费观看在线| 91午夜精品亚洲一区二区三区| 日韩一区二区三区影片| 婷婷色综合大香蕉| 18禁动态无遮挡网站| 日韩大片免费观看网站| 欧美一级a爱片免费观看看| 综合色av麻豆| 日韩国内少妇激情av| 一个人观看的视频www高清免费观看| 免费观看在线日韩| 黄色视频在线播放观看不卡| 亚洲成人精品中文字幕电影| 亚洲精品乱码久久久v下载方式| 亚洲精品色激情综合| 久久精品熟女亚洲av麻豆精品| 国产亚洲最大av| 久久久久久久久久成人| 最近中文字幕高清免费大全6| 三级男女做爰猛烈吃奶摸视频|