關(guān)于約束大渦模擬方法的一些思考

夏振華, 史一蓬

(1. 浙江大學 航空航天學院工程力學系, 杭州 310027; 2. 北京大學工學院 湍流與復雜系統(tǒng)國家重點實驗室, 北京 100871)

0 引 言

高雷諾數(shù)壁湍流是國防和航空航天領(lǐng)域關(guān)鍵的問題之一。對高雷諾數(shù)壁湍流的精確模擬和控制是指導飛行器航天器的設(shè)計和改進的基礎(chǔ)，是國家安全、經(jīng)濟建設(shè)和科學技術(shù)中的重點研究對象之一，也是制約我國航空航天事業(yè)發(fā)展的基礎(chǔ)科學難題之一。近些年，隨著計算機技術(shù)突飛猛進的發(fā)展以及數(shù)值計算方法的進一步完善，數(shù)值模擬方法正逐漸成為湍流研究的重要手段之一。

經(jīng)典的湍流數(shù)值模擬方法有三種：直接數(shù)值模擬方法(Direct Numerical Simulation, DNS)，大渦模擬方法(Large-eddy Simulation，LES)和雷諾平均方法(Raynolds-Averaged Navier-Stokes，RANS)。由于網(wǎng)格需求大，格式精度要求高，現(xiàn)階段DNS方法只停留在簡單流動的物理機理研究上[1]。現(xiàn)在工程上最常用的方法還是更易于實現(xiàn)的RANS[2-5]。但是雷諾應(yīng)力模型普適性較差，特別是在大分離、非定常流動時效果不好。LES方法的計算量介于RANS和DNS之間，它可以用來模擬大分離、非定常流動，一般模型參數(shù)也相對較少[6]。但是對于高雷諾數(shù)壁湍流而言，在壁面解析的情況下，LES的網(wǎng)格需求仍然很大[7-8]，因此在工程應(yīng)用中，壁面解析的LES使用得并不多。

為了更好地模擬大分離、非定常流動，人們提出了各種改進的方法。例如，壁面?；腖ES方法(Wall-modelled LES， WMLES)[9-10]，分離渦模擬方法(Detached Eddy Simulation, DES)等RANS-LES混合方法(Hybrid RANS-LES Method)[11-14]。RANS-LES混合方法結(jié)合了RANS方法和LES方法的優(yōu)點，同時又避免了這兩種方法各自的缺點，對分離流動的模擬效果大大優(yōu)于RANS方法。RANS-LES混合方法被認為是未來數(shù)十年湍流模擬最有希望的方法。雖然以DES為代表的RANS-LES混合方法被越來越多地用來模擬復雜流動問題，但是此類方法有一個最致命的問題就是在模擬附著壁湍流時，存在對數(shù)區(qū)不匹配現(xiàn)象(Log-layer Mismatch, LLM)。例如，如果用1997年原始的DES方法[15]來模擬槽道湍流，會出現(xiàn)兩段不同的對數(shù)律區(qū)間，其中底層的RANS區(qū)域有一段對數(shù)律，而外層的LES區(qū)域也有一段對數(shù)律。兩段對數(shù)律的斜率非常接近，但是截距差別明顯[16]。此外，如果DES方法使用的網(wǎng)格不是很好，使得DES方法在很粗的網(wǎng)格里切換到LES分支，這樣會使得模型化部分的應(yīng)力偏小，而總體雷諾應(yīng)力又沒有能力彌補這個缺失，從而出現(xiàn)所謂的模型應(yīng)力虧損(Modeled Stress Depletion, MSD)[17]，甚至有可能出現(xiàn)網(wǎng)格引起的虛假分離(Grid-induced Separation，GIS)[17-18]。為了克服這些缺點，研究人員相繼提出了DDES(Delayed DES)[17]、IDDES(improved DDES)[14,18]等修正方法。雖然IDDES方法在槽道湍流的驗證中消除了LLM現(xiàn)象，但是其引入了新的尺度定義、新的混合函數(shù)和數(shù)個模型參數(shù)，從而使得其適用性大為降低。Spalart悲觀地認為IDDES是DES方法發(fā)展的一個不幸的趨勢[12]。

為了消除LLM，本文作者及其合作者提出了約束大渦模擬方法(Constrained Large Eddy Simulation, CLES)[19-20]。CLES方法同之前的DES類方法不同，它全場作LES計算，但是在靠近壁面的區(qū)域采用約束了的SGS模型。CLES方法解決了RANS-LES混合方法最主要的缺陷——對近壁區(qū)的模擬缺少小尺度脈動，從而消除了LLM。本文將首先詳細介紹CLES方法的基本原理，而后對CLES方法應(yīng)用過程中的一些問題展開討論分析，最后對CLES方法未來的一些發(fā)展方向做了一點概括。對于CLES方法的數(shù)值驗證，我們將簡略帶過，讀者有興趣可以參考一些已發(fā)表的文獻。

1 CLES方法的基本原理

CLES方法是湍流約束思想[21-22]在LES建模中的應(yīng)用，本質(zhì)上還是一種LES方法。與傳統(tǒng)的LES方法不同的是，它在靠近壁面的內(nèi)區(qū)需要使用約束的亞網(wǎng)格模型(Constrained Subgrid-scale Model， CSGS)。接下來，我們將首先推導CSGS模型的基本約束等式。不失一般性，我們接下來以不可壓縮問題為例。對于三維非定常不可壓縮常密度流體，其控制方程如下：

(1)

可壓縮問題情況類似，具體推導可參考文獻[20]。

1.1 RANS控制方程

對上述方程(1)做系綜平均運算〈〉，便可得到RANS方程組：

(2)

這里Rij=〈uiuj〉-〈ui〉〈uj〉為雷諾應(yīng)力張量。

1.2 LES控制方程

類似地，如果我們對方程(1)做濾波運算～，就可以得到LES方程組：

(3)

從方程(2)、(3)可見，LES和RANS方程在形式上完全一致。這也是RANS/LES混合方法的最基本出發(fā)點。

1.3 CLES的約束等式

如果我們將方程(3)做一個系綜平均，那么其動量方程會變成：

(4)

(5)

(6)

即：總雷諾應(yīng)力可以分解為亞網(wǎng)格應(yīng)力的平均值和可解雷諾應(yīng)力兩部分。對于可壓縮問題，我們同樣推導出了類似的關(guān)于雷諾應(yīng)力的約束等式以及關(guān)于熱流的約束等式[20]。

1.4 CLES方法

有了上面的雷諾應(yīng)力約束等式，我們就可以構(gòu)建CSGS模型。對于亞網(wǎng)格應(yīng)力張量，我們首先將其做一個分解：

(7)

它的兩部分可以分別給出：其平均部分可以由上面的雷諾應(yīng)力約束等式(6)給出:

而其脈動部分我們可以簡單地由已有的SGS模型來構(gòu)造，以Smagorinsky模型[23]為例：

(8)

在CLES方法的具體應(yīng)用的時候，在外區(qū)，我們可以使用傳統(tǒng)的SGS模型，而在靠近壁面的內(nèi)區(qū)，我們則可以使用CSGS模型。通過上面的推導，我們也發(fā)現(xiàn)，雖然我們?nèi)珗龆荚谧鯨ES計算，但是我們在內(nèi)區(qū)實際上做了一個隱含的RANS計算。由于CLES方法在內(nèi)區(qū)和外區(qū)都開展LES計算，因此小尺度脈動能夠在全場存在，不會出現(xiàn)RANS/LES混合方法那樣由于缺乏近壁小尺度脈動而導致的LLM現(xiàn)象。

CLES方法在不可壓縮/可壓縮槽道湍流[19-20]、周期丘陵繞流[26]、不可壓縮/可壓縮圓柱繞流[19,27-28]等算例中得到了驗證，并被用來模擬了大型客機外型在迎角14°時的流動狀態(tài)[29-31]。圖1給出了在不可壓縮槽道湍流中驗證結(jié)果。由圖可見，當CLES方法用來模擬不可壓縮槽道湍流時，能得到與DNS結(jié)果基本吻合的平均速度型。相反，DES方法出現(xiàn)了LLM，而傳統(tǒng)的LES-DSM則會出現(xiàn)平均速度在對數(shù)區(qū)的整體上移。而且我們還發(fā)現(xiàn)，在槽道湍流中，使用代數(shù)渦黏模型和SA模型來預(yù)估雷諾應(yīng)力，對結(jié)果影響不大。

圖1Reτ=1000時幾種不同方法計算得到的平均速度型[19](RSC-A和RSC-SA分別指代使用代數(shù)渦黏模型和SA模型來預(yù)估雷諾應(yīng)力的CLES模擬結(jié)果)

Fig.1 Mean velocity profiles atReτ=1000 by using different simulation methods[19](RSC-A and RSC-SA denote the CLES simulations using the algebraic eddy viscosity model and the SA model to estimate the prescribed Reynolds stresses, respectively)

2 關(guān)于CLES方法的一些討論

CLES方法中，CSGS模型至關(guān)重要。

2.1 關(guān)于CSGS模型的思考

CSGS建模的想法可追溯到湍流中的約束關(guān)系的提出[21-22]。圖2給出了一些常見的SGS模型的物理/數(shù)學考量。由圖我們可以看出，大多數(shù)SGS模型都是有或多或少的物理/數(shù)學考量的。當我們將更多的物理/數(shù)學考量加入到SGS模型中，我們就能得到更好的模型。例如當將Germano等式這個數(shù)學約束加入到常系數(shù)Smagorinsky模型之后，我們就可以得到比常系數(shù)模型更準確的動態(tài)Smagorinsky模型(Dynamic Smagorinsky Model, DSM)。類似地，我們?nèi)绻趧討B(tài)混合相似性模型(dynamic mixed similarity model)基礎(chǔ)上加入能流約束，我們就可以得到精確度更高的能流約束模型[32]；在壁湍流的模擬中，如果我們在DSM基礎(chǔ)上加入雷諾應(yīng)力約束，就能得到雷諾應(yīng)力約束模型[19-20]。

需要特別指出的是，約束等式(6)在數(shù)學的意義上是一個必要條件，而不是一個充要條件，這和Germano等式類似。事實上，Meneveau[33]曾經(jīng)推導過亞網(wǎng)格應(yīng)力模型需要滿足的一些必要條件。在他推導的必要條件里，可以允許兩者之間相差一個無散的常張量Cij。在這里，我們考慮到在壁面上的相容條件，總雷諾應(yīng)力、可解雷諾應(yīng)力和亞網(wǎng)格雷諾應(yīng)力都需要為0，因此Cij=0。對于其他的約束，我們可以允許差一個常數(shù)張量，但是它們在實際計算中并沒有貢獻，因為雷諾應(yīng)力是以散度的形式出現(xiàn)在動量方程中的。

此外，我們還需要說明的是基于雷諾應(yīng)力約束的CSGS模型更多的是一種方法論，采用不同基準的SGS模型，可以得到不同的CSGS模型。例如，在文獻[34]中，我們就考慮了基于Vreman模型[35]的CSGS模型。

圖2 SGS模型的物理/數(shù)學考量Fig.2 Some basic considerations on the SGS models

2.2 總雷諾應(yīng)力的影響

2.3 CLES方法與傳統(tǒng)LES和DES的聯(lián)系

如果我們以常系數(shù)Smagorinsky模型為基礎(chǔ)，此時CSGS模型(8)可以被改寫為：

(9)

類似地，我們可以將方程式(8)作另一種改寫：

(10)

如果內(nèi)區(qū)的總雷諾應(yīng)力是由求解一方程或者兩方程的DES型(即原始的RANS模型方程中的特征長度尺度按照DES方法做了修改)RANS模型方程來預(yù)估的，而外區(qū)的SGS模型也是使用DES的LES分支(也可以繼續(xù)使用Smagorinsky模型)，那么此時CSGS模型可以看成是在DES方法的基礎(chǔ)上，在內(nèi)區(qū)加了一個附加力，該力為張量的散度。值得提出的是，Keating和Piomelli[36]曾經(jīng)提出過在界面處加入隨機里的方法來減弱RANS/LES混合方法中的LLM現(xiàn)象，但這種附加的隨機力本身與流場沒有相關(guān)性，需要根據(jù)實際需要調(diào)節(jié)隨機力的強度。CSGS模型事實上在整個內(nèi)區(qū)都增加一個與流場相關(guān)的附加力。該附加力使得小尺度脈動在內(nèi)區(qū)一致存在，并且在界面附近得到充分發(fā)展，從而可以消除LLM。事實上，在文獻[19]中，我們分析了DES方法中的模型應(yīng)力和可解應(yīng)力的法向變化規(guī)律。同傳統(tǒng)的LES和CLES不同，DES方法在內(nèi)區(qū)是RANS分支，此時模型應(yīng)力占主導，而可解部分的應(yīng)力其實是隨著壁面距離的增加而緩慢增加的。這種外區(qū)長度尺度的切換，從外部也改變了內(nèi)部的雷諾應(yīng)力分布，從而使得總體雷諾應(yīng)力偏小。而CLES方法也可以視為在內(nèi)區(qū)增加了一個力，這個力可以使得內(nèi)區(qū)的脈動迅速增長起來，從而改變了內(nèi)區(qū)的可解雷諾應(yīng)力。

2.4 約束界面的影響

在CLES方法的應(yīng)用中，其內(nèi)外區(qū)界面的位置決定了使用CSGS模型的區(qū)域，從而對結(jié)果會產(chǎn)生不同程度的影響?？梢灶A(yù)見的是，如果約束區(qū)域很小，那么其效果就會不好，當界面的位置設(shè)置于固壁上，CLES中CSGS模型的使用區(qū)域就會隨著消失，CLES方法也會退化為傳統(tǒng)的LES方法。而如果CLES中CSGS區(qū)域太大，一個問題是我們很難保證遠離壁面的外區(qū)我們還能得到準確的雷諾應(yīng)力，特別是針對復雜的分離流動而言。此外，約束區(qū)域的增大也會相應(yīng)增加計算量，畢竟CSGS模型的計算量比傳統(tǒng)的SGS模型要大。在槽道湍流等簡單流動中，我們研究過不同約束界面位置對結(jié)果的影響。如分析的那樣，如果約束界面位于緩沖區(qū)(y+=5.05)，那么平均速度會上移；如果界面位于對數(shù)區(qū)(y+=20.1和44.9),結(jié)果基本相同[19]。如果用CLES方法模擬復雜壁湍流，如何選取界面位置的確是一個需要研究的問題。目前大多數(shù)計算中，CLES的界面位置會參考DES的界面位置，并比DES的位置略微靠近壁面一點。

此外，我們還需要考慮內(nèi)外區(qū)之間的銜接問題。雖然CSGS模型和傳統(tǒng)的SGS模型在數(shù)學上不存在連續(xù)性問題，但是在實際計算中，兩種不同的模型在固定的界面處連接在一起，勢必會在界面處出現(xiàn)一些統(tǒng)計量不連續(xù)的情況。為了避免這個問題，我們可以參照一些RANS/LES混合方法的做法，可以在界面處使用一個連接函數(shù)將兩個區(qū)域連接起來。一個最簡單的嘗試就是界面處的亞網(wǎng)格應(yīng)力為兩種模型的平均值。

2.5 系綜平均的實現(xiàn)

在CSGS模型中，我們看到需要進行系綜平均的操作，這是一個耗費巨大計算量的工作。在簡單的流動里，系統(tǒng)本身存在統(tǒng)計均勻方向，此時我們可以用統(tǒng)計均勻方向的空間平均來近似系綜平均。例如在槽道湍流中，我們可以使用面平均來近似系綜平均。然而，在一般的復雜流動中，可能不存在統(tǒng)計均勻方向，這時我們無法采用空間平均來近似系綜平均。在這類問題中，我們可以仿照Meneveau等[37]構(gòu)造Lagrangian平均來近似，但是Lagrangian平均往往計算復雜，計算量也很大，在復雜流動的計算中很難實用。另一個選擇是采用Euler時間平均來近似系綜平均：

〈f〉n+1=εfn+1+(1-ε)〈f〉n

ε=Δt/(Δt+T)

(11)

這里Δt是計算的時間步長，T是Euler時間平均的特征間隔時間。在槽道湍流里，我們驗證了時間平均能夠得到和空間平均相近的結(jié)果[19]。

3 CLES方法的未來研究建議

CLES方法經(jīng)過近十年的發(fā)展，已經(jīng)在很多經(jīng)典問題中得到驗證，并已經(jīng)被用于模擬復雜的工程流動問題。未來，我們希望CLES方法能夠被更多地應(yīng)用到工程湍流問題的模擬中去。對于CLES方法的理論研究，我們認為一下幾個方向值得研究：

(2) 復雜問題模擬中，如何給出準確的總雷諾應(yīng)力的預(yù)估值？總雷諾應(yīng)力是CLES方法中的目標值，因此對CLES方法來說，總雷諾應(yīng)力的預(yù)估值至關(guān)重要，會直接影響CLES最終的模擬結(jié)果。在簡單的壁湍流問題中，大多數(shù)RANS模型都可以給出較好的預(yù)測結(jié)果；在一般的復雜流動中，如何在CSGS應(yīng)用區(qū)域中給出合理的總雷諾應(yīng)力的預(yù)估值是值得研究的問題。在這里，基于機器學習的RANS建?？赏峁┮环N新思路。

(3) 約束界面位置的選取方法。在過去的研究中，我們通常采用的固定的界面位置。我們知道的是，固定的界面位置不能太靠近壁面，否則約束效果不好。如何選取合理的界面位置是一個值得研究的問題。

(4) 約束大渦模擬在氣動聲學和氣動光學中的應(yīng)用。CLES能夠模擬出比傳統(tǒng)RANS-LES混合方法更加豐富而精細的小尺度流動結(jié)構(gòu)，這些小尺度流動結(jié)構(gòu)對光學傳播有重要的影響，同時也是氣動噪聲的主要的聲源。因此CLES方法在氣動聲學和氣動光學的模擬方面有廣闊的應(yīng)用前景。