趙耀民 *, 徐曉偉
* (北京大學(xué)應(yīng)用物理與技術(shù)研究中心,高能量密度物理數(shù)值模擬教育部重點(diǎn)實(shí)驗(yàn)室,北京 100871)
? (北京大學(xué)工學(xué)院力學(xué)與工程科學(xué)系,湍流與復(fù)雜系統(tǒng)國(guó)家重點(diǎn)實(shí)驗(yàn)室,北京 100871)
** (澳大利亞墨爾本大學(xué)機(jī)械工程系,澳大利亞墨爾本 3010)
湍流廣泛存在于自然界和工程應(yīng)用中,準(zhǔn)確認(rèn)識(shí)和預(yù)測(cè)湍流現(xiàn)象十分重要.湍流的研究有近百年的歷史,人們由理論、實(shí)驗(yàn)、和數(shù)值模擬方向出發(fā)取得了一系列的研究成果[1].然而,在以航空航天為代表的工程應(yīng)用中,流動(dòng)幾何邊界往往較為復(fù)雜,且流動(dòng)雷諾數(shù)較高,研究的難度很大,相應(yīng)流動(dòng)的基本規(guī)律和預(yù)測(cè)模型仍有待進(jìn)一步深入研究.因此,湍流的機(jī)理、建模及控制仍然是航空航天等領(lǐng)域的重要瓶頸之一,針對(duì)湍流的研究有著重要的基礎(chǔ)和應(yīng)用價(jià)值.
計(jì)算流體動(dòng)力學(xué)(computational fluid dynamics,CFD)是湍流的重要研究手段,常見(jiàn)的數(shù)值模擬方法可大致分為直接數(shù)值模擬(direct numerical simulation,DNS)、大渦模擬(large eddy simulation,LES)和雷諾平均模擬(Reynolds-averaged Navier-Stokes,RANS)等.其中,直接數(shù)值模擬和大渦模擬均需要較高的網(wǎng)格分辨精度.對(duì)于航空航天等工程應(yīng)用中的高雷諾數(shù)流動(dòng),直接數(shù)值模擬和大渦模擬所需要的網(wǎng)格量巨大,遠(yuǎn)超當(dāng)前超級(jí)計(jì)算機(jī)的算力極限[2].相對(duì)而言,RANS 雖然精度不及前兩種方法,但計(jì)算效率很高,是各種工程應(yīng)用中的主要研究手段.
雷諾平均模擬的結(jié)果往往受限于湍流模型的預(yù)測(cè)精度.以雷諾應(yīng)力輸運(yùn)模型[3]為代表的高階矩封閉模型雖然精度較高,但計(jì)算量較大且收斂性較差,其應(yīng)用范圍較為有限.實(shí)際工程中的常用湍流模型包括 SA 模型[4]、k-ε模型[5]、k-ω模型[6]、SST 模型[7]等.這些模型均基于Boussinesq 渦黏假設(shè),即假設(shè)雷諾應(yīng)力的各向異性部分與應(yīng)變率張量存在線性關(guān)系.然而,當(dāng)流動(dòng)較為復(fù)雜,特別是出現(xiàn)如邊界層分離、強(qiáng)烈混合等現(xiàn)象時(shí),雷諾應(yīng)力與應(yīng)變之間往往不再滿足簡(jiǎn)單的線性關(guān)系,因此基于渦黏假設(shè)的模型將會(huì)出現(xiàn)較大誤差[8-9].
近年來(lái),隨著相關(guān)算法的快速發(fā)展和可用數(shù)據(jù)的迅速增長(zhǎng),機(jī)器學(xué)習(xí)方法開(kāi)始在流體力學(xué)領(lǐng)域得到廣泛關(guān)注[10],其中一個(gè)重要的研究方向是利用機(jī)器學(xué)習(xí)方法開(kāi)發(fā)湍流模型[11-12].如圖1 所示,可以將數(shù)據(jù)驅(qū)動(dòng)的湍流建模大致分為3 個(gè)步驟,即數(shù)據(jù)預(yù)處理、模型訓(xùn)練以及模型驗(yàn)證和預(yù)測(cè).在數(shù)值模擬或者實(shí)驗(yàn)獲得的高可信度(high-fidelity,Hi-Fi)流場(chǎng)基礎(chǔ)上,首先通過(guò)數(shù)據(jù)預(yù)處理提取輸入特征變量(X表示)及輸出變量(Y表示);再由輸入及輸出變量出發(fā),采用不同的機(jī)器學(xué)習(xí)算法及策略訓(xùn)練模型;最后一般要對(duì)訓(xùn)練得到的模型(Y=f(X))進(jìn)行一定的測(cè)試,驗(yàn)證其預(yù)測(cè)能力.需要注意的是,最終模型的預(yù)測(cè)能力依賴于訓(xùn)練數(shù)據(jù)、特征變量的選擇、算法及策略等多種因素.
圖1 機(jī)器學(xué)習(xí)湍流建模示意圖Fig.1 Schematic for turbulence modelling with machine learning methods
根據(jù)輸出變量的不同,現(xiàn)有的數(shù)據(jù)驅(qū)動(dòng)湍流建模研究大致可分為兩個(gè)類(lèi)型.其中一部分研究關(guān)注于發(fā)展修正模型,即采用不同的機(jī)器學(xué)習(xí)方法修正現(xiàn)有的模型及其參數(shù),以改善原有模型的預(yù)測(cè)精度.Edeling 等[13]基于貝葉斯方法估計(jì)壓力梯度邊界層中模型誤差并嘗試修正模型參數(shù).Parish 和Duraisamy[14]提出了流場(chǎng)反演與機(jī)器學(xué)習(xí)(field inversion and machine learning,FIML)方法.此框架首先通過(guò)統(tǒng)計(jì)方法反演得到修正參數(shù)場(chǎng)作為訓(xùn)練目標(biāo),例如湍動(dòng)能方程中湍流生成項(xiàng)的修正系數(shù)的空間分布;進(jìn)一步地,采用機(jī)器學(xué)習(xí)方法重構(gòu)這一修正場(chǎng),建立輸入特征量與修正場(chǎng)之間的映射關(guān)系.FIML 框架的優(yōu)勢(shì)在于較為靈活,其中流場(chǎng)反演部分可以使用貝葉斯方法[15],也可基于伴隨方法求解[16];另外,機(jī)器學(xué)習(xí)部分也可以由高斯過(guò)程、神經(jīng)網(wǎng)絡(luò)等不同方法實(shí)現(xiàn)[17-18].值得注意的是,這種框架中反演得到的修正場(chǎng)并不總是“可訓(xùn)練的”[19],而訓(xùn)練過(guò)程的缺陷會(huì)導(dǎo)致反演場(chǎng)的誤差,進(jìn)而導(dǎo)致最終預(yù)測(cè)精度和泛化能力的不足.一種可能的改進(jìn)方向是將反演和學(xué)習(xí)兩個(gè)過(guò)程耦合,但是相關(guān)的研究還較為初步[19].
另外一個(gè)大類(lèi)的工作可歸類(lèi)為替代模型[12]的開(kāi)發(fā).研究者們不滿足于單純修正現(xiàn)有的模型參數(shù),而是探索構(gòu)建滿足一定物理約束的RANS 模型以替代原有模型.針對(duì)RANS 和高精度流場(chǎng)的雷諾應(yīng)力偏差,Wang 等[20]采用隨機(jī)森林方法構(gòu)建模型,在周期山測(cè)試算例中可以明顯改善雷諾應(yīng)力場(chǎng)的預(yù)測(cè)精度,但是模型在后驗(yàn)測(cè)試中對(duì)于平均流場(chǎng)的預(yù)測(cè)精度還需要進(jìn)一步驗(yàn)證.Ling 等[21]采用深度神經(jīng)網(wǎng)絡(luò)方法訓(xùn)練滿足伽利略不變性的雷諾應(yīng)力模型.訓(xùn)練基于速度梯度基張量和不變量[22],得到的模型在后驗(yàn)測(cè)試中可以改善對(duì)于方腔流的預(yù)測(cè)精度.Yin等[23]同樣基于神經(jīng)網(wǎng)絡(luò)方法,深入研究輸入特征的選擇標(biāo)準(zhǔn)以優(yōu)化模型訓(xùn)練框架.優(yōu)化后的方法應(yīng)用于周期山算例,顯著改善了流場(chǎng)的預(yù)測(cè)精度.張偉偉等[12,24-25]在機(jī)器學(xué)習(xí)湍流建模領(lǐng)域也開(kāi)展了一系列工作,例如利用神經(jīng)網(wǎng)絡(luò)方法直接構(gòu)建基于平均流場(chǎng)的渦黏場(chǎng)預(yù)測(cè)模型[26],針對(duì)近壁區(qū)、尾流區(qū)以及遠(yuǎn)場(chǎng)分別建模,成功地預(yù)測(cè)了機(jī)翼的升阻力系數(shù)、阻力分布等.Xie 等[27-29]在一系列工作中采用不同方法針對(duì)大渦模擬的亞格子應(yīng)力進(jìn)行建模,在后驗(yàn)測(cè)試中可以更為準(zhǔn)確地預(yù)測(cè)湍流能譜等.
在眾多數(shù)據(jù)驅(qū)動(dòng)湍流建模研究中,澳大利亞墨爾本大學(xué)Weatheritt 和Sandberg[30]采用基因表達(dá)式編程方法(gene-expression programming,GEP)[30],在最近的一些研究取得了一系列進(jìn)展.GEP 方法[31]是遺傳算法的一種,可以較為方便地實(shí)現(xiàn)符號(hào)回歸,近年來(lái)已在非穩(wěn)態(tài)RANS[32]、湍流傳熱[33]、大渦模擬亞格子應(yīng)力[34]以及燃燒[35]等多個(gè)領(lǐng)域的建模中得到了初步應(yīng)用.相較于神經(jīng)網(wǎng)絡(luò)等方法的“黑箱”模型,GEP 方法可以顯式給出模型方程,因此可以很方便地應(yīng)用于工程軟件中[36-37].
本文旨在對(duì)基于GEP 方法的湍流顯式模型相關(guān)工作進(jìn)行簡(jiǎn)要總結(jié).文章的結(jié)構(gòu)安排如下:第1 部分介紹基于GEP 的模型訓(xùn)練方法及其框架;第2 部分,介紹方法的幾種典型應(yīng)用;最后第3 部分是結(jié)論與展望.
GEP 方法[28]是遺傳算法的一種,基于達(dá)爾文自然選擇理論,通過(guò)種群的演化實(shí)現(xiàn)優(yōu)化過(guò)程.這里僅對(duì)GEP 方法進(jìn)行簡(jiǎn)要介紹,方法的具體細(xì)節(jié)可參考文獻(xiàn)[30].GEP 方法中的種群往往包含大量個(gè)體,而每個(gè)個(gè)體均代表一個(gè)候選模型.如圖2 所示,模型的基因型(genotype)是一個(gè)字符串,又被稱作染色體.該字符串由不同符號(hào)組成,包括運(yùn)算符(如 +,×等)、常數(shù)項(xiàng)及變量符號(hào)(x,y).在將字符串表達(dá)為方程形式的過(guò)程中,首先將不同符號(hào)按照樹(shù)結(jié)構(gòu)(即表達(dá)樹(shù),expression tree)編碼,將運(yùn)算符放置在非終端節(jié)點(diǎn),而常數(shù)項(xiàng)和變量符號(hào)放置在終端節(jié)點(diǎn).進(jìn)一步地,基于此樹(shù)結(jié)構(gòu)將其解碼為表現(xiàn)型(phenotype),該表達(dá)型就是一個(gè)候選模型方程.這一基于樹(shù)結(jié)構(gòu)的表達(dá)過(guò)程保證公式符合語(yǔ)法規(guī)則,可以被應(yīng)用于實(shí)際運(yùn)算.
圖2 GEP 方法的基因型及其表達(dá)過(guò)程Fig.2 Genotype and expression process in GEP
如圖3 所示,在訓(xùn)練的開(kāi)始,首先隨機(jī)生成一個(gè)包含大量個(gè)體的種群,種群中個(gè)體的優(yōu)劣可以由使用者定義的損失函數(shù)得到.類(lèi)似于自然界中的自然選擇過(guò)程,種群中損失函數(shù)較低的個(gè)體的“健壯度”較好,因而有更大的概率進(jìn)入下一代.進(jìn)一步地,通過(guò)種群繁殖和基因變異等操作,可以得到更新一代的種群.正是通過(guò)這樣的迭代過(guò)程,整個(gè)種群向著更適應(yīng)自然選擇標(biāo)準(zhǔn)的方向演化,演化最終得到的最優(yōu)個(gè)體即為訓(xùn)練得到的模型.
圖3 GEP 方法流程圖Fig.3 Schematic for GEP method
相較于其他機(jī)器學(xué)習(xí)方法,GEP 方法的一個(gè)顯著特點(diǎn)是可以通過(guò)符號(hào)回歸得到顯式的模型方程,因而具有模型可解釋、易應(yīng)用的特點(diǎn).另外,訓(xùn)練過(guò)程的尋優(yōu)過(guò)程主要由基因的遺傳和變異實(shí)現(xiàn).相較于梯度下降等方法,這種訓(xùn)練過(guò)程更接近全局搜索,不易收斂于局部最優(yōu)解;但同時(shí)收斂過(guò)程具有一定的隨機(jī)性,收斂速度無(wú)法得到保證,特別是可能很難得到模型中常系數(shù)的最優(yōu)解[31].
GEP 方法近年來(lái)被應(yīng)用于不同領(lǐng)域的湍流建模問(wèn)題.其中,既包括RANS 計(jì)算中的雷諾應(yīng)力封閉模型,也包括對(duì)湍流傳熱問(wèn)題的建模.本節(jié)分別介紹這兩種典型的建模過(guò)程.
1.2.1 顯式代數(shù)應(yīng)力模型
在雷諾平均模擬中,對(duì)Navier-Stokes 方程進(jìn)行雷諾平均,動(dòng)量方程中將產(chǎn)生一個(gè)非封閉項(xiàng),即雷諾應(yīng)力.這里ui是湍流脈動(dòng)速度分量.Boussinesq(1877) 假定雷諾剪切應(yīng)力與平均速度梯度間為線性相關(guān),該線性假定可推廣到各項(xiàng)異性雷諾應(yīng)力張量ai j與平均應(yīng)變率張量Sij之間,即
其中,k=為湍動(dòng)能,δij為 Kronecker 符號(hào),μt為湍流動(dòng)力學(xué)黏度系數(shù),νt=μt/ρ 為渦黏系數(shù).
Boussinesq 渦黏假設(shè)在大多數(shù)情況下并不適用.Pope[22]基于量綱分析和坐標(biāo)變換不變量,進(jìn)一步添加了非線性項(xiàng),將其擴(kuò)展為
式中,Ui,j為平均速度梯度,湍流時(shí)間尺度τ=1/ω=0.09k/ε ,ε 為湍動(dòng)能耗散率,ω 為比耗散率.在三維流場(chǎng)中,存在10 個(gè)線性無(wú)關(guān)的張量基與5 個(gè)獨(dú)立的標(biāo)量不變量,具體如下
式(2)中給出的顯式代數(shù)應(yīng)力模型是Boussinesq線性渦黏假設(shè)的一種修正,也是許多數(shù)據(jù)驅(qū)動(dòng)湍流建??蚣艿幕A(chǔ)[21,27-28].在GEP 湍流建模中,以標(biāo)量不變量Ik、常數(shù)項(xiàng)、以及數(shù)學(xué)運(yùn)算符號(hào)(如+,×)將作為輸入符號(hào),通過(guò)符號(hào)回歸的方式訓(xùn)練模型系數(shù) ζk的代數(shù)表達(dá)式.
1.2.2 湍流傳熱模型
與雷諾應(yīng)力相似,對(duì)能量方程進(jìn)行雷諾平均模擬時(shí)也將產(chǎn)生另一個(gè)非封閉項(xiàng),即湍流熱通量傳統(tǒng)湍流傳熱的建模一般根據(jù)標(biāo)準(zhǔn)梯度擴(kuò)散假設(shè)(standard gradient diffusion hypothesis,SGDH)
式中 Θ,j為平均溫度梯度,αt=νt/Prt為湍流擴(kuò)散系數(shù),Prt為湍流普朗特?cái)?shù).實(shí)際工程計(jì)算中往往將湍流普朗特?cái)?shù)取常數(shù).例如,空氣中Prt取值為0.90 左右.然而這種常數(shù)湍流普朗特?cái)?shù)的假設(shè)會(huì)引入較大的誤差,因此在機(jī)器學(xué)習(xí)湍流建模中,人們往往將Prt(或者 αt)作為訓(xùn)練目標(biāo).通過(guò)建立Prt隨平均溫度梯度等變量的分布函數(shù),可以提升湍流傳熱的預(yù)測(cè)精度.
另外,SGDH 中的擴(kuò)散系數(shù)為各向同性,且假設(shè)湍流熱通量與平均溫度梯度線性相關(guān),這也是模型誤差的另外一個(gè)可能來(lái)源.針對(duì)這一問(wèn)題,一種更為通用的模型可以將 αt替換為湍流擴(kuò)散系數(shù)張量Di j[38]
為進(jìn)一步提高湍流傳熱模型的計(jì)算精度,在機(jī)器學(xué)習(xí)建模中往往采用更為泛化的湍流擴(kuò)散系數(shù)張量模型.本文采用 Weatheritt 等[33]對(duì)Younis 等[40]提出的顯式湍流擴(kuò)散系數(shù)張量模型的簡(jiǎn)潔表述
與雷諾應(yīng)力相似,在湍流傳熱建模中,機(jī)器學(xué)習(xí)的目標(biāo)是建立模型方程(7)中g(shù)m的代數(shù)表達(dá)式.其中,輸入代數(shù)符號(hào)為標(biāo)量不變量Ik和Jγ、常數(shù)項(xiàng)以及數(shù)學(xué)運(yùn)算符號(hào)(如+,×).
機(jī)器學(xué)習(xí)方法訓(xùn)練得到的模型在應(yīng)用前必須進(jìn)行充分的測(cè)試.在湍流建模領(lǐng)域,一般將測(cè)試分為先驗(yàn)(a priori)和后驗(yàn)(a posteriori).所謂先驗(yàn)測(cè)試,即將訓(xùn)練或者測(cè)試的高精度流場(chǎng)數(shù)據(jù)作為輸入特征變量代入模型方程,進(jìn)而將模型預(yù)測(cè)的輸出變量與高精度數(shù)據(jù)作對(duì)比,以檢測(cè)模型對(duì)于雷諾應(yīng)力(或是湍流傳熱)的預(yù)測(cè)精度.與之相對(duì)應(yīng),后驗(yàn)測(cè)試則需要將訓(xùn)練得到的模型應(yīng)用于實(shí)際RANS 計(jì)算,因此模型的輸入特征變量由RANS 計(jì)算的流場(chǎng)提供,而模型的預(yù)測(cè)精度則根據(jù)RANS 計(jì)算的結(jié)果進(jìn)行判定.需要指出的是,一些數(shù)據(jù)驅(qū)動(dòng)的湍流模型盡管可以給出令人滿意的先驗(yàn)測(cè)試結(jié)果,但是在實(shí)際計(jì)算中可能導(dǎo)致較大誤差[41].
在機(jī)器學(xué)習(xí)建模中,先驗(yàn)測(cè)試評(píng)價(jià)模型本身對(duì)于高精度訓(xùn)練數(shù)據(jù)的擬合和預(yù)測(cè)能力,而后驗(yàn)測(cè)試則關(guān)注模型能否在實(shí)際計(jì)算中更好地預(yù)測(cè)流場(chǎng)信息.對(duì)于湍流建模而言,盡管先驗(yàn)測(cè)試十分重要,但是后驗(yàn)測(cè)試的CFD 計(jì)算中的預(yù)測(cè)效果往往才是最終的目標(biāo).
機(jī)器學(xué)習(xí)算法中損失函數(shù)的定義往往對(duì)于訓(xùn)練結(jié)果有決定性作用.對(duì)于湍流建模問(wèn)題而言,也同樣需要強(qiáng)調(diào)損失函數(shù)的重要性.湍流建模問(wèn)題中一個(gè)顯著的特點(diǎn)是訓(xùn)練數(shù)據(jù)與最終模型的應(yīng)用環(huán)境可能有很大區(qū)別.一方面,研究者一般認(rèn)為直接數(shù)值模擬或者大渦模擬的結(jié)果為“真實(shí)值”,作為訓(xùn)練和測(cè)試數(shù)據(jù);另一方面,模型的最終需要應(yīng)用于雷諾平均模擬,而雷諾平均模擬中的湍動(dòng)能、湍流耗散時(shí)間尺度等往往與訓(xùn)練數(shù)據(jù)存在較大區(qū)別.這會(huì)給損失函數(shù)的定義帶來(lái)一定的困難.本節(jié)介紹兩種不同的損失函數(shù)定義方法.
1.4.1 凍結(jié)訓(xùn)練
Akolekar 等[37]基于GEP 方法提出了凍結(jié)訓(xùn)練方法.圖4 以二方程雷諾平均模型為例介紹凍結(jié)訓(xùn)練的流程.由直接數(shù)值模擬等高可信度流場(chǎng),可以得到平均速度、平均速度梯度張量、雷諾應(yīng)力、湍動(dòng)能等數(shù)據(jù),并以此作為“真實(shí)值”.在此基礎(chǔ)上,將“凍結(jié)”的平均速度和湍動(dòng)能代入模型輸運(yùn)方程,如k-ω模型中的ω方程,得到適用于雷諾平均模擬的湍流耗散時(shí)間尺度.采用這種湍流特征時(shí)間尺度對(duì)速度梯度張量進(jìn)行無(wú)量綱化,可以得到1.2.1 小節(jié)中的張量基函數(shù)和標(biāo)量不變量作為模型的輸入變量.進(jìn)一步地,使用GEP 方法訓(xùn)練基于無(wú)量綱速度梯度張量基函數(shù)和標(biāo)量不變量的顯式代數(shù)應(yīng)力模型.
圖4 凍結(jié)訓(xùn)練流程示意圖Fig.4 Schematic for ‘frozen’ training method
在這一訓(xùn)練過(guò)程中,通過(guò)求解k-ω模型中的ω方程得到湍流耗散時(shí)間尺度,用于進(jìn)一步的速度梯度張量無(wú)量綱化.在求解過(guò)程中,平均速度場(chǎng)和湍動(dòng)能來(lái)自于“凍結(jié)的”高精度訓(xùn)練數(shù)據(jù),在迭代過(guò)程中不發(fā)生變化.這里僅對(duì)ω進(jìn)行迭代更新,直至收斂.因此該方法被稱為凍結(jié)訓(xùn)練.
凍結(jié)訓(xùn)練方法中采用求解模型方程得到的湍流耗散時(shí)間尺度進(jìn)行無(wú)量綱化,而不是直接使用高精度訓(xùn)練數(shù)據(jù)中的特征時(shí)間尺度.這是因?yàn)橹苯訑?shù)值模擬等方法得到的湍流耗散時(shí)間尺度與實(shí)際雷諾平均模擬運(yùn)算中的結(jié)果往往存在較大差異.這種差異可以認(rèn)為是高精度訓(xùn)練數(shù)據(jù)與雷諾平均模擬之間存在一定的相容性問(wèn)題.換而言之,完全依據(jù)直接數(shù)值模擬等高精度流場(chǎng)訓(xùn)練得到的雷諾應(yīng)力模型可能并不適用于雷諾平均模擬,這種情況下盡管模型可能在先驗(yàn)測(cè)試中準(zhǔn)確度較高,但在后驗(yàn)測(cè)試中仍然存在較大誤差.因此,必須在模型訓(xùn)練過(guò)程中考慮這種差異,通過(guò)凍結(jié)方法增強(qiáng)模型在實(shí)際雷諾平均模擬中的適用性,提升其后驗(yàn)精度.
凍結(jié)訓(xùn)練以修正模型預(yù)測(cè)的雷諾應(yīng)力為目標(biāo),對(duì)應(yīng)的損失函數(shù)是通過(guò)對(duì)比高精度訓(xùn)練數(shù)據(jù)中的雷諾應(yīng)力和模型給出的雷諾應(yīng)力來(lái)得到.通過(guò)求解模型輸運(yùn)方程,凍結(jié)訓(xùn)練引入適用于雷諾平均模擬的特征時(shí)間尺度進(jìn)行無(wú)量綱化,從而可以在一定程度上緩解高精度訓(xùn)練數(shù)據(jù)與RANS 模擬的相容性問(wèn)題.然而,最近一些學(xué)者研究發(fā)現(xiàn)[42-43],即便將直接數(shù)值模擬得到的雷諾應(yīng)力代入RANS 計(jì)算中,得到的平均流場(chǎng)也可能會(huì)產(chǎn)生較大的誤差.Wu 等[41]還發(fā)現(xiàn),這種預(yù)測(cè)誤差在雷諾數(shù)增高時(shí)會(huì)進(jìn)一步增大.這種現(xiàn)象導(dǎo)致數(shù)據(jù)驅(qū)動(dòng)的RANS 模型很難在后驗(yàn)驗(yàn)證中準(zhǔn)確預(yù)測(cè)平均流動(dòng).
1.4.2 雙向耦合方法
針對(duì)湍流建模中常見(jiàn)的高精度訓(xùn)練數(shù)據(jù)和RANS計(jì)算的相容性問(wèn)題,Zhao 等[44]首次提出了機(jī)器學(xué)習(xí)和CFD 雙向耦合的模型訓(xùn)練方法.雙向耦合方法的最大特點(diǎn)是訓(xùn)練過(guò)程中模型的優(yōu)劣通過(guò)實(shí)際RANS 計(jì)算的結(jié)果來(lái)評(píng)估,而不是僅僅基于高精度訓(xùn)練數(shù)據(jù)判定對(duì)于雷諾應(yīng)力的預(yù)測(cè)效果.如圖5 所示,對(duì)于每一代的候選模型,將一系列模型方程讀入到RANS 求解器中,并基于該模型實(shí)現(xiàn)RANS 計(jì)算;進(jìn)一步將計(jì)算的結(jié)果返回GEP 訓(xùn)練過(guò)程,以評(píng)估模型方程的損失函數(shù).由于在訓(xùn)練過(guò)程中,每一代模型的評(píng)估均需要GEP 算法和RANS 計(jì)算共同作用,因此稱之為機(jī)器學(xué)習(xí)和CFD 雙向耦合的訓(xùn)練方法.雙向耦合訓(xùn)練方法的步驟可以總結(jié)為以下幾個(gè)步驟:
圖5 雙向耦合方法訓(xùn)練流程示意圖Fig.5 Schematic for CFD-driven machine learning method
(1)初始化種群,得到一系列模型方程;
(2)種群中的一系列候選模型均由CFD 運(yùn)算進(jìn)行評(píng)估,這一過(guò)程可以并行操作,具體步驟如下:
①模型方程寫(xiě)入文件,再由RANS 求解器讀入此文件;
② RANS 求解器基于GEP 算法提供的候選模型運(yùn)行CFD 模擬;
③RANS 計(jì)算的結(jié)果分兩類(lèi).如果RANS 收斂,模型的損失函數(shù)可由收斂后流場(chǎng)計(jì)算得到;如果RANS 不收斂,模型損失函數(shù)設(shè)為無(wú)窮大;
④ 將CFD 計(jì)算得到的模型損失函數(shù)返回GEP 算法,以評(píng)估候選模型;
⑤ GEP 算法迭代,模型種群繁殖、變異,進(jìn)入下一代,直到得到期望的模型.
相比于凍結(jié)方法,雙向耦合方法具有以下特點(diǎn).首先,雙向耦合方法的訓(xùn)練目標(biāo)靈活,可以根據(jù)RANS 計(jì)算的結(jié)果設(shè)定不同的損失函數(shù)作為優(yōu)化目標(biāo);相對(duì)而言,凍結(jié)方法只能對(duì)比模型預(yù)測(cè)的雷諾應(yīng)力與高精度訓(xùn)練數(shù)據(jù).其次,雙向耦合方法對(duì)于訓(xùn)練數(shù)據(jù)的需求量較少,實(shí)驗(yàn)測(cè)量得到的少許平均流場(chǎng)數(shù)據(jù)即可作為訓(xùn)練數(shù)據(jù);相對(duì)而言,凍結(jié)訓(xùn)練往往需要高精度數(shù)值模擬得到的全場(chǎng)雷諾應(yīng)力等高階湍流統(tǒng)計(jì)量.最后,由于雙向耦合方法得到的模型在訓(xùn)練過(guò)程中就已經(jīng)經(jīng)過(guò)RANS 計(jì)算的檢驗(yàn),因此模型的收斂性和后驗(yàn)預(yù)測(cè)精度均有一定保證;相對(duì)而言,凍結(jié)方法得到的模型在應(yīng)用于真實(shí)RANS 計(jì)算時(shí),對(duì)于平均流場(chǎng)的預(yù)測(cè)精度可能存在較大誤差.雙向耦合方法面臨的主要挑戰(zhàn)在于計(jì)算效率方面.與凍結(jié)訓(xùn)練相比,盡管雙向耦合方法在GEP 訓(xùn)練中所需的種群代數(shù)相差不大(一般約為數(shù)百代),但是由于其迭代過(guò)程中每一代候選模型的評(píng)估均需運(yùn)行RANS 計(jì)算,因此訓(xùn)練所需的計(jì)算量往往遠(yuǎn)大于凍結(jié)方法.
第1 部分介紹了基因表達(dá)式編程方法以及如何基于此方法訓(xùn)練湍流和傳熱模型.接下來(lái),介紹GEP 方法在幾種湍流建模問(wèn)題中的具體應(yīng)用.
2.1.1 統(tǒng)計(jì)穩(wěn)態(tài)流動(dòng)的RANS 模型
在航空發(fā)動(dòng)機(jī)內(nèi)流中的葉片尾沿下游,葉片壓力側(cè)和吸力側(cè)的來(lái)流在強(qiáng)剪切層的主導(dǎo)下產(chǎn)生強(qiáng)烈混合,這是典型的尾流混合問(wèn)題.尾流混合往往導(dǎo)致顯著的流動(dòng)損失,因此準(zhǔn)確模擬尾流混合問(wèn)題對(duì)于內(nèi)流預(yù)測(cè)和葉片設(shè)計(jì)十分重要,然而基于Boussinesq假設(shè)的常用工程湍流模型對(duì)于尾流損失的預(yù)測(cè)精度往往不甚理想[45].最近幾年來(lái),GEP 方法在一系列的工作中[37,44,46-47]被應(yīng)用于尾流混合問(wèn)題,針對(duì)如圖6所示的高壓渦輪和低壓渦輪平面葉柵算例,目標(biāo)為發(fā)展具有較高預(yù)測(cè)精度和較強(qiáng)泛化能力的尾流混合模型.
圖6 航空發(fā)動(dòng)機(jī)內(nèi)流葉柵尾流混合算例示意圖Fig.6 Simulation setup for cases
表1 給出了一系列渦輪葉柵算例的參數(shù).基于直接數(shù)值模擬或是大渦模擬提供的高精度流場(chǎng)數(shù)據(jù),這些算例被用于模型的訓(xùn)練和測(cè)試.以典型航空發(fā)動(dòng)機(jī)工況下的高壓渦輪(表1 算例A) 為研究目標(biāo),Weatheritt 等[47]采用1.2.1 小節(jié)中介紹的凍結(jié)方法訓(xùn)練顯式代數(shù)應(yīng)力模型.訓(xùn)練中損失函數(shù)設(shè)置為
表1 渦輪葉柵尾流混合算例Table 1 Parameters of turbine wake mixing cases
即在N個(gè)網(wǎng)格點(diǎn)上,GEP 模型得到的雷諾應(yīng)力各向異性部分與高精度訓(xùn)練數(shù)據(jù)的相對(duì)誤差.訓(xùn)練得到的模型方程形式如下
對(duì)這一模型進(jìn)行先驗(yàn)測(cè)試,Weatheritt 等[47]發(fā)現(xiàn)相較于傳統(tǒng)Boussinesq 假設(shè),該模型可以明顯改善對(duì)于尾流區(qū)域雷諾應(yīng)力的預(yù)測(cè).但是,由于存在大量高階非線性項(xiàng),該模型在后驗(yàn)測(cè)試中很難得到收斂的流場(chǎng).
為進(jìn)一步發(fā)展魯棒性更強(qiáng)的尾流混合模型,Zhao 等[44]提出了1.2.2 小節(jié)中介紹的雙向耦合訓(xùn)練方法,并將其應(yīng)用于高壓渦輪算例A.不同于式(10)中給出的凍結(jié)訓(xùn)練損失函數(shù),雙向耦合方法訓(xùn)練過(guò)程中損失函數(shù)的設(shè)置較為靈活,不僅可以對(duì)比雷諾應(yīng)力等高階統(tǒng)計(jì)量,更可以選擇實(shí)驗(yàn)、直接數(shù)值模擬等高可信度數(shù)據(jù)中可提供的任意感興趣的平均流場(chǎng)信息.以渦輪尾流混合問(wèn)題為例,這里選擇根據(jù)RANS 模擬結(jié)果中的平均尾流損失設(shè)置損失函數(shù).具體形式如下
這里 Ω (y) 為尾流損失曲線,,po和pt(y) 分別代表來(lái)流總壓、出口靜壓和當(dāng)?shù)乜倝?Δx是GEP 模型預(yù)測(cè)的尾流損失曲線與高可信度流場(chǎng)數(shù)據(jù)的相對(duì)誤差,而最終損失函數(shù)JCFD則是兩個(gè)不同位置的相對(duì)誤差之和.值得注意的是,在這樣的雙向耦合訓(xùn)練過(guò)程中,僅從高可信度流場(chǎng)中提取少數(shù)平均尾流損失曲線作為訓(xùn)練數(shù)據(jù),而不需要凍結(jié)方法中必須的全場(chǎng)雷諾應(yīng)力信息,因此對(duì)于訓(xùn)練數(shù)據(jù)的要求大大降低.
采用雙向耦合方法訓(xùn)練并設(shè)置式(12)中的損失函數(shù),Zhao 等[44]得到的尾流混合模型如下
將此模型與凍結(jié)方法獲得的模型式(10) 作對(duì)比,發(fā)現(xiàn)雙向耦合方法得到的模型方程中高階非線性項(xiàng)大大減少.其原因在于這些高階項(xiàng)往往導(dǎo)致實(shí)際RANS 計(jì)算難以收斂,因此在雙向耦合訓(xùn)練過(guò)程中被自然選擇所淘汰.由此可見(jiàn),雙向耦合方法獲得的模型一般具有更強(qiáng)的收斂性和魯棒性.
將式(11)中的雙向耦合模型應(yīng)用于后驗(yàn)測(cè)試,并將其與去除高階項(xiàng)后的簡(jiǎn)化凍結(jié)模型的預(yù)測(cè)結(jié)果作對(duì)比.如圖7(a) 所示,相較于作為基準(zhǔn)模型的k-ωSST 模型和凍結(jié)模型,雙向耦合方法在后驗(yàn)測(cè)試中可以顯著提高對(duì)于式(12)中尾流損失 的預(yù)測(cè)精度.進(jìn)一步對(duì)模型進(jìn)行分析,Zhao 等[44]指出雙向耦合模型對(duì)于尾流損失預(yù)測(cè)的提升,其主要原因在于方程中占主導(dǎo)的項(xiàng)的系數(shù)幅值更大(2 .57 > 1 .334).考慮到,這說(shuō)明該模型在尾流區(qū)域?qū)е赂鼜?qiáng)的擴(kuò)散,因而雙向耦合模型預(yù)測(cè)的尾流損失的峰值和尾流寬度均更接近于大渦模擬的結(jié)果.由此可見(jiàn),GEP 方法訓(xùn)練得到的湍流模型是顯式方程,這樣可以分析、解釋模型預(yù)測(cè)精度提升的作用機(jī)制.
圖7 不同渦輪葉柵算例(見(jiàn)表1)中的尾流損失剖面Fig.7 Kinetic wake loss profiles from test cases in Table 1
為進(jìn)一步驗(yàn)證雙向耦合模型的泛化能力,還將其應(yīng)用于表1 所示的算例B,C 和D 中.這些算例與訓(xùn)練算例A 相比,其雷諾數(shù)、出流馬赫數(shù)、攻角以及幾何外形均有很大差別.如圖7(b)~圖7(d)所示,在這一系列后驗(yàn)測(cè)試中,這種雙向耦合方法訓(xùn)練得到的顯式代數(shù)應(yīng)力模型可以較為準(zhǔn)確地預(yù)測(cè)尾流混合導(dǎo)致的流動(dòng)損失.因此該模型具有較強(qiáng)的泛化能力.
2.1.2 非穩(wěn)態(tài)流動(dòng)的unsteady RANS 建模策略
上一節(jié)討論了針對(duì)統(tǒng)計(jì)穩(wěn)態(tài)流動(dòng)的建模過(guò)程,可以看到GEP 方法能有效提升RANS 對(duì)于尾流混合的預(yù)測(cè)精度.但是當(dāng)流動(dòng)處于非穩(wěn)態(tài)時(shí),往往應(yīng)采用非定常雷諾平均(unsteady RANS)模擬,而建模策略可能也需要相應(yīng)調(diào)整.
航空發(fā)動(dòng)機(jī)內(nèi)流領(lǐng)域一種常見(jiàn)的非穩(wěn)態(tài)現(xiàn)象是兩排葉片之間存在的相對(duì)運(yùn)動(dòng)導(dǎo)致的流動(dòng)周期性變化.Akolekar 等[48]針對(duì)來(lái)流為周期性尾流擾動(dòng)的低壓渦輪葉柵展開(kāi)研究,探討機(jī)器學(xué)習(xí)湍流建模的策略.圖8 給出了鎖相平均的葉柵湍動(dòng)能云圖.將一個(gè)來(lái)流擾動(dòng)周期分為20 個(gè)不同的相,由圖8 可發(fā)現(xiàn)周期性來(lái)流尾流對(duì)渦輪葉柵流場(chǎng)造成顯著影響,不同相的流動(dòng)之間呈現(xiàn)顯著差別.
圖8 來(lái)流尾流擾動(dòng)下低壓渦輪葉柵的鎖相平均湍動(dòng)能云圖[48]Fig.8 Phase-lock averaged TKE contours for LPT flow disturbed by incoming wakes[48]
針對(duì)這種非穩(wěn)態(tài)流動(dòng),湍流建模策略可有幾種:(1)針對(duì)不同相的流場(chǎng)單獨(dú)建模;(2)根據(jù)時(shí)間平均后的流場(chǎng)建模;(3)根據(jù)一定的流動(dòng)物理分組建模.顯然,針對(duì)不同相分別建模會(huì)給出適用于各相的一系列湍流模型,但是模型間的切換會(huì)給使用者的實(shí)際模擬帶來(lái)不便.Akolekar 等[46]進(jìn)一步發(fā)現(xiàn),針對(duì)時(shí)間平均的流場(chǎng)建模,雖然可以一定程度改善預(yù)測(cè)精度,但是由于忽略了各相間的差異導(dǎo)致其精度仍然有限.相對(duì)而言,最優(yōu)的策略是根據(jù)流動(dòng)物理將不同相分組,進(jìn)而針對(duì)各組具有相似物理的流動(dòng)分別建模.
在這種根據(jù)流動(dòng)物理對(duì)不同相進(jìn)行分組的過(guò)程中,GEP 方法可以起到關(guān)鍵作用.針對(duì)低頻來(lái)流尾流擾動(dòng)下的低壓渦輪葉柵流動(dòng),Akolekar 等[48]首先對(duì)20 相分別建模.得到這些單相模型后,一方面可以分別測(cè)試各模型的預(yù)測(cè)精度,發(fā)現(xiàn)單相訓(xùn)練得到的模型在不同相的測(cè)試中表現(xiàn)并不理想.另一方面,也可以根據(jù)GEP 方法顯式給出的方程形式,對(duì)不同相的模型進(jìn)行分析.在分析過(guò)程中,Akolekar 等[48]發(fā)現(xiàn)不同相的模型方程大致可分為兩組.對(duì)于其中一組(包含大約10 相),其模型方程中項(xiàng)的系數(shù)幅值遠(yuǎn)大于另外一組(包含剩余10 相),這說(shuō)明這些相中需要引入相對(duì)更強(qiáng)的湍流擴(kuò)散.進(jìn)一步分析流場(chǎng)特征,可以發(fā)現(xiàn)第一組的流動(dòng)主要特征是來(lái)流尾流與葉片邊界層相互作用,進(jìn)而導(dǎo)致葉柵尾流混合強(qiáng)度顯著增強(qiáng);而另外一組流動(dòng)中來(lái)流尾流與葉片則沒(méi)有明顯的相互作用,相應(yīng)需要的湍流擴(kuò)散修正較弱.這一分組過(guò)程首先依據(jù)GEP 訓(xùn)練的單相模型方程得到,進(jìn)一步得到了流場(chǎng)分析的驗(yàn)證.在此分組的基礎(chǔ)上,可以針對(duì)兩組流動(dòng)分別建模,得到適用于對(duì)應(yīng)組的兩個(gè)模型.在先驗(yàn)測(cè)試中,這種基于流動(dòng)物理的模型可顯著提升模型的預(yù)測(cè)精度.
近年來(lái),GEP 方法也被應(yīng)用于溫度或標(biāo)量場(chǎng)的建模之中.本文首先針對(duì)豎直平板間自然對(duì)流這一較為簡(jiǎn)單的算例,介紹基于GEP 的標(biāo)量系數(shù)湍流傳熱模型;然后針對(duì)更為復(fù)雜的三維橫向流中的射流算例,介紹張量系數(shù)湍流傳熱建模.
2.2.1 自然對(duì)流的湍流標(biāo)量系數(shù)建模
自然對(duì)流在室內(nèi)通風(fēng)環(huán)境(如雙層玻璃窗)、電子設(shè)備冷卻等場(chǎng)景中十分常見(jiàn).圖9 給出了以浮力驅(qū)動(dòng)豎直充分發(fā)展槽道流中的自然對(duì)流計(jì)算域示意圖,其中兩個(gè)豎直板沿x和z方向采用周期性邊界,y方向左右兩側(cè)為不同溫度的等溫壁面.兩壁面中,左側(cè)高溫?zé)o量綱化后為 Θh=0.5 ,右側(cè)低溫為 Θc=-0.5,有限寬度H=2h.對(duì)于豎直槽道流中空氣流動(dòng),流場(chǎng)特性由Rayleigh 數(shù)Ra=gβΔΘH3/(υκ) 決定.其中ΔΘ=Θh-Θc是兩個(gè)豎直板間無(wú)量綱的溫度差,g為重力加速度,β 為體積膨脹系數(shù),ν 為運(yùn)動(dòng)黏度,κ 為擴(kuò)散系數(shù),空氣普朗特?cái)?shù)Pr=ν/κ=0.709 .
圖9 以豎直槽道流中的自然對(duì)流算例的計(jì)算域Fig.9 Computational domain of natural convection case in a differentially heated vertical channel
Ng 等[49]通過(guò)直接數(shù)值模擬(DNS)中得到了該算例 的平均流場(chǎng)、雷諾應(yīng)力和湍流熱通量.Xu 等[50]以此數(shù)據(jù)集中的算例為訓(xùn)練數(shù)據(jù),提取湍流渦黏系數(shù),并以垂直壁面方向的熱通量為學(xué)習(xí)目標(biāo).基于GEP方法,通過(guò)最小化損失函數(shù)
得到模型方程f1=0.969+2J0,嵌入到SGDH 中為
根據(jù) αt=vt/Prt,這里湍流普朗特系數(shù)等價(jià)于Prt=1/(0.969+2J0).與此對(duì)比,本文以常用的常湍流普朗特?cái)?shù)為基準(zhǔn)模型(baseline model),其中Prt=0.9 .
圖10 湍流普朗特?cái)?shù)的先驗(yàn)預(yù)測(cè)Fig.10 Prediction of turbulent Prandtl number
在后驗(yàn)測(cè)試中,Xu 等[50]在給定DNS 的速度場(chǎng)相關(guān)數(shù)據(jù)的基礎(chǔ)上,求解溫度的對(duì)流-擴(kuò)散輸運(yùn)方程,來(lái)測(cè)試湍流傳熱模型的預(yù)測(cè)效果.該模型被應(yīng)用于3 個(gè)不同Ra算例,圖11 分別展示了后驗(yàn)測(cè)試中基準(zhǔn)模型和GEP 模型預(yù)測(cè)的平均溫度剖面和垂直壁面方向的熱通量.可以看到,基準(zhǔn)模型在平均溫度預(yù)測(cè)上有較大的誤差,且該誤差隨著Ra增大而增大.與此對(duì)比,GEP 模型能夠有效的改善熱通量的預(yù)測(cè),進(jìn)而修正平均溫度剖面,并且計(jì)算結(jié)果與DNS 的真實(shí)值趨于一致.需要強(qiáng)調(diào)的是,上述先驗(yàn)和后驗(yàn)測(cè)試均采用了Ra=2.0×107的數(shù)據(jù)訓(xùn)練得到的GEP 模型,然后在不同下測(cè)試該模型.以上結(jié)果表明,GEP 模型不僅在訓(xùn)練算例中表現(xiàn)優(yōu)秀,而且能夠在其他數(shù)據(jù)集的交叉驗(yàn)證中驗(yàn)證其其有效性.因此,該模型具有較好的泛化能力.
圖11 平均溫度剖面和垂直壁面方向的熱通量的后驗(yàn)預(yù)測(cè)Fig.11 Prediction of mean temperature profiles and wall-normal heat flux
2.2.2 橫向射流的湍流張量擴(kuò)散系數(shù)建模
三維橫向流中的射流 (jet in cross-flow,JICF)在航空發(fā)動(dòng)機(jī)等渦輪葉片中有廣泛應(yīng)用,準(zhǔn)確預(yù)測(cè)其流場(chǎng)和溫度場(chǎng)是工業(yè)界亟待解決的難題.圖12 給出了典型JICF 的算例設(shè)置示意圖.其中,x為流向,y為垂直壁面方向,z為展向,射流孔直徑為d,橫向主流無(wú)量綱溫度 Θin,冷卻射流無(wú)量綱溫度 Θj.流場(chǎng)的主要特征參數(shù)為雷諾數(shù)Rej=Ujd/ν ,入射角 φ,吹風(fēng)比(blowing ratio)BR=Uj/Uin.算例中高溫橫向主流(freestream inflow)與從氣膜孔流出的冷卻氣體射流(coolant inflow)相互作用,形成摻混區(qū)域.由于摻混區(qū)域存在逆梯度輸運(yùn)(counter-gradient transport)、擬序結(jié)構(gòu)形成與破碎等復(fù)雜三維非定常流動(dòng)現(xiàn)象,基于RANS 對(duì)其湍流傳熱的預(yù)測(cè)難度非常大.
圖12 橫流中的射流示意圖Fig.12 The schematic description of jet in crossflow case
Weatheritt 等[33]基于Bodart 等[51]的JICF 大渦模擬數(shù)據(jù)進(jìn)行GEP 湍流傳熱建模,其中展向?qū)挾葹?d,豎向高度為 8.62d,φ=30°,BR=1.0,Rej=5.4×103,這里稱該算例為Bodart 算例.由于流場(chǎng)的復(fù)雜性,標(biāo)準(zhǔn)梯度擴(kuò)散假設(shè)難以準(zhǔn)確預(yù)測(cè)傳熱過(guò)程,因此必須采用湍流擴(kuò)散系數(shù)張量模型對(duì)三維流場(chǎng)中的熱通量進(jìn)行數(shù)據(jù)驅(qū)動(dòng)建模.
可以發(fā)現(xiàn)張量系數(shù)相關(guān)度的排序?yàn)?/p>
在此基礎(chǔ)上,GEP 模型訓(xùn)練中符號(hào)回歸的方程設(shè)置為
此外,損失函數(shù)為
其中,V為三維流場(chǎng)中數(shù)據(jù)集所有離散點(diǎn)數(shù).經(jīng)過(guò)機(jī)器學(xué)習(xí)訓(xùn)練之后,得出如下模型
該模型被應(yīng)用于Bodart 算例.圖13 給出了流向x/d=20處截面(見(jiàn)圖12)上垂直壁面方向的熱通量預(yù)測(cè)結(jié)果,并將GEP 模型、基準(zhǔn)模型(SGDH) 與LES 數(shù)據(jù)作對(duì)比.在不同展向位置(z/d=0.50,0.67,1.10),與SGDH 相比,GEP 模型提高了對(duì)兩個(gè)熱通量分量vθ和wθ 的預(yù)測(cè)精度.此外,沿垂直壁面方向上,GEP模型預(yù)測(cè)的峰值點(diǎn)位置也與LES 結(jié)果較好吻合.
圖13 流向位置x/d=20 的熱通量預(yù)測(cè)Fig.13 The prediction of heat flux vector at x/d=20
此外,該模型還被應(yīng)用于Sakai 等[52]的JICF 大渦模擬算例(Sakai 算例),以測(cè)試基于Bodart 算例的湍流張量擴(kuò)散系數(shù)模型的泛化能力.與Bodart 算例相比,Sakai 算例的幾何尺寸和主要特征參數(shù)均不同,其中,展向?qū)挾葹?3d,豎向高度為 5d,φ =35°,Rej=1.5×104.圖14 比較了兩組不同吹風(fēng)比BR=0.5 和BR=1.0算例中壁面絕熱效率η的流向分布,其中,壁面絕熱效率定義為.當(dāng)BR=1.0(與Bodart 算例一致),GEP 模型對(duì)橫向平均和中線處 η 符合高精度數(shù)據(jù)的總體趨勢(shì),在摻混區(qū)域與參考數(shù)據(jù)非常接近.當(dāng)BR=0.5,GEP 模型對(duì)橫向平均和中線處 η 的預(yù)測(cè)精度均優(yōu)于基準(zhǔn)模型.
圖14 壁面絕熱效率的流向分布Fig.14 Wall adiabatic effectiveness streamwise distribution
綜上所述,無(wú)論是針對(duì)湍流標(biāo)量系數(shù)(如湍流普朗特?cái)?shù))還是對(duì)湍流張量擴(kuò)散系數(shù)進(jìn)行建模,GEP 方法均能提升對(duì)湍流傳熱的預(yù)測(cè)精度,且具有一定程度的泛化能力.
自Weatheritt 和Sandberg[30]首次將GEP 應(yīng)用于湍流建模以來(lái),此方法已在多個(gè)相關(guān)領(lǐng)域得到了成功應(yīng)用.2.1 節(jié)和2.2 節(jié)重點(diǎn)介紹了針對(duì)代數(shù)應(yīng)力模型和湍流傳熱問(wèn)題的典型建模過(guò)程.在此基礎(chǔ)上,在本節(jié)簡(jiǎn)要介紹GEP 方法在湍流建模領(lǐng)域的一些拓展和其他應(yīng)用.
與1.2.1 小節(jié)中的RANS 建模過(guò)程非常類(lèi)似,GEP 方法還可應(yīng)用于大渦模擬亞格子應(yīng)力建模.Schoepplein 等[35]首次將GEP 方法應(yīng)用于湍流預(yù)混火焰的亞格子應(yīng)力建模.與RANS 建模中以平均流場(chǎng)為訓(xùn)練數(shù)據(jù)不同,亞格子應(yīng)力的建?;跒V波后的直接數(shù)值模擬瞬時(shí)流場(chǎng).該訓(xùn)練采用凍結(jié)方法,所得到的亞格子應(yīng)力代數(shù)表達(dá)式在先驗(yàn)測(cè)試中具有較高的準(zhǔn)確度.但是,該模型并未經(jīng)過(guò)后驗(yàn)測(cè)試,且模型的泛化能力也需要進(jìn)一步驗(yàn)證.最近,Reissmann 等[34]嘗試將雙向耦合方法[44]應(yīng)用于泰勒-格林渦的亞格子應(yīng)力建模中.在后驗(yàn)測(cè)試中,所得模型可較為準(zhǔn)確預(yù)測(cè)湍動(dòng)能等統(tǒng)計(jì)量演化過(guò)程.當(dāng)然,這一模型對(duì)充分發(fā)展湍流的適用性還有待一步驗(yàn)證.
在最近的工作中,GEP 方法還被嘗試應(yīng)用于邊界層轉(zhuǎn)捩問(wèn)題的建模[53].基于層流動(dòng)能(laminar kinetic energy,LKE)方法[54],Akolekar 等[53]針對(duì)LKE 方程中的層流渦黏系數(shù)等關(guān)鍵參數(shù)提出修正的代數(shù)模型.在后驗(yàn)測(cè)試中,該模型實(shí)現(xiàn)了對(duì)于低壓渦輪葉柵中分離引起的邊界層轉(zhuǎn)捩的準(zhǔn)確預(yù)測(cè).雖然該模型的魯棒性等問(wèn)題還有待進(jìn)一步研究,但該結(jié)果顯示了GEP 方法在邊界層轉(zhuǎn)捩等復(fù)雜流動(dòng)的建模中具有良好的應(yīng)用前景.
近年來(lái),機(jī)器學(xué)習(xí)輔助的湍流建模研究發(fā)展迅速,然而在訓(xùn)練數(shù)據(jù)和訓(xùn)練方法、模型可解釋性和泛化能力等方面仍然存在較大挑戰(zhàn).本文討論的一系列工作基于GEP 方法在這一領(lǐng)域取得了一定進(jìn)展,實(shí)現(xiàn)了多種類(lèi)型的湍流建模.在訓(xùn)練數(shù)據(jù)方面,這些工作關(guān)注具有較強(qiáng)實(shí)際應(yīng)用背景的復(fù)雜流動(dòng),以發(fā)展實(shí)用工程湍流模型為目標(biāo).在訓(xùn)練方法方面,通過(guò)引入雙向耦合方法等模型訓(xùn)練框架,考慮RANS計(jì)算與高精度訓(xùn)練數(shù)據(jù)的差別,強(qiáng)調(diào)模型的實(shí)際預(yù)測(cè)能力.另外,由于GEP 方法可以提供顯式模型方程,因此模型可解釋性較強(qiáng).
值得注意的是,目前大多數(shù)機(jī)器學(xué)習(xí)方法得到的湍流模型距離真實(shí)工程應(yīng)用還存在明顯的差距.接下來(lái)的工作重點(diǎn)之一是提高模型的泛化能力.一種具有較強(qiáng)應(yīng)用前景的數(shù)據(jù)驅(qū)動(dòng)模型應(yīng)該可以較好地預(yù)測(cè)多種類(lèi)型流動(dòng)(如遠(yuǎn)離壁面的自由剪切流動(dòng)和近壁面湍流)及多種流動(dòng)物理(如強(qiáng)混合、轉(zhuǎn)捩、流動(dòng)分離等).此外,針對(duì)雙向耦合方法,還可以通過(guò)優(yōu)化GEP 算法的收斂速度或種群大小、改善CFD求解器等方式來(lái)提升訓(xùn)練的計(jì)算效率.以上述工作重點(diǎn)為目標(biāo),需要進(jìn)一步完善現(xiàn)有的訓(xùn)練數(shù)據(jù)和訓(xùn)練方法.