任海杰 袁先旭, 陳堅(jiān)強(qiáng), 孫 東 朱林陽 向星皓,

* (空氣動力學(xué)國家重點(diǎn)實(shí)驗(yàn)室,四川綿陽 621000)

? (中國空氣動力研究與發(fā)展中心計(jì)算空氣動力研究所,四川綿陽 621000)

引言

湍流模型是為了封閉Navier-Stokes (N-S)方程中的雷諾應(yīng)力項(xiàng)而額外補(bǔ)充的方程,主要目的是構(gòu)建時(shí)均流動、空間位置與雷諾應(yīng)力張量或湍流渦粘之間的數(shù)學(xué)關(guān)系式[1].雷諾平均(Reynolds averaged Navier-Stokes,RANS)方程是當(dāng)前湍流研究中使用最廣泛的模型之一.RANS 框架的核心是確定流場的平均量,并針對脈動量對于平均場本身的影響進(jìn)行建模;而脈動場則是通過雷諾應(yīng)力張量的散度嵌入到RANS 方程中,RANS 模型正是需要對此張量進(jìn)行建模[2].最早的有效黏度假設(shè)由Boussinesq 提出,該假設(shè)也是使用最為廣泛的有效黏度假設(shè),通常使用的兩方程RANS 模型(例如k-ε,k-ω 模型)通過Boussinesq 渦黏假設(shè)完成雷諾應(yīng)力的封閉.但Tracey 等[3]指出湍流渦黏模型的主要誤差來源是無法解釋雷諾應(yīng)力的各向異性,即Boussinesq 渦黏假設(shè)認(rèn)為湍流黏性系數(shù)是各向同性的標(biāo)量.此外,有諸多學(xué)者根據(jù)推導(dǎo)出的雷諾應(yīng)力輸運(yùn)方程進(jìn)行求解,例如雷諾應(yīng)力方程模型(RSM),閆超等[4]指出RSM 存在雷諾應(yīng)力方程的建模困難、數(shù)值剛性問題較嚴(yán)重、計(jì)算量較大的問題.不同于RSM 的思路,Pope[5]針對各向異性提出了新的有效黏度假設(shè).該假設(shè)的有效性已經(jīng)得到了充分的驗(yàn)證[6],但由于其在N-S 求解器中較差的魯棒性,并未得到廣泛應(yīng)用.這也意味著研究者往往需要從DNS,LES 等高分辨率數(shù)據(jù)中獲得更加準(zhǔn)確的雷諾應(yīng)力各向異性張量.

隨著機(jī)器學(xué)習(xí)方法的不斷發(fā)展,在雷諾應(yīng)力的封閉中結(jié)合機(jī)器學(xué)習(xí)的研究方法日益增加.主要的研究方向分為兩種[7]:采用機(jī)器學(xué)習(xí)方法縮小RANS解與高精度數(shù)據(jù)或?qū)嶒?yàn)值的偏差;或是基于高精度數(shù)據(jù)或?qū)嶒?yàn)值直接構(gòu)建某些湍流變量的代理模型.縮小偏差的方法既可以是改變模型的控制方程形式,也可以構(gòu)建針對RANS 模型偏差函數(shù)進(jìn)行疊加修正.代表性工作包括Tracey 等[8]構(gòu)建了替代SA 模型控制方程中源項(xiàng)的神經(jīng)網(wǎng)絡(luò)模型;Wu 等[9]針對RANS 模型計(jì)算結(jié)果和高分辨率DNS 數(shù)據(jù)之間的雷諾應(yīng)力偏差進(jìn)行建模,提高了原有模型的準(zhǔn)確性.相比于縮小偏差的方法,直接構(gòu)建代理模型的方法同樣頗具亮點(diǎn).Zhu 等[10-11]通過徑向基函數(shù)神經(jīng)網(wǎng)絡(luò)以及深度神經(jīng)網(wǎng)絡(luò)構(gòu)建了渦黏系數(shù)νt 與平均流動變量之間的映射關(guān)系,并將得到的映射關(guān)系與N-S求解器耦合用于封閉湍流模型.而Ling 等[12]針對張量不變性分別在隨機(jī)森林和神經(jīng)網(wǎng)絡(luò)中進(jìn)行了測試.隨機(jī)森林作為決策樹的集合,是將多個(gè)決策樹的預(yù)測組合成一個(gè)模型,其結(jié)構(gòu)并不像神經(jīng)網(wǎng)絡(luò)一樣容易改變,難以將Pope 的本構(gòu)關(guān)系嵌入到隨機(jī)森林中.而神經(jīng)網(wǎng)絡(luò)由于體系結(jié)構(gòu)的靈活性在張量不變性的處理上別具優(yōu)勢.為了在RANS 模型的基礎(chǔ)上,利用機(jī)器學(xué)習(xí)方法提高雷諾應(yīng)力各向異性的計(jì)算精度,Ling 等[7]采用具有伽利略不變性的輸入特征并提出張量基神經(jīng)網(wǎng)絡(luò)(TBNN)架構(gòu).諸多成果表明機(jī)器學(xué)習(xí)方法在湍流模型化工作中具有良好的前景[13].

對于TBNN,Ling 等[7]使用TBNN 預(yù)測了低速的管流、周期山流動.Fang 等[2]針對槽道流提出了新的神經(jīng)網(wǎng)絡(luò)模型,使其在低速槽道流的預(yù)測效果優(yōu)于TBNN.而該研究也指出,其新的神經(jīng)網(wǎng)絡(luò)未能嵌入旋轉(zhuǎn)不變性,且僅適用于特定流動.張珍等[14]也將TBNN 與一個(gè)預(yù)測渦黏系數(shù)的神經(jīng)網(wǎng)絡(luò)嵌入到了RANS 求解器中,并對低速周期山流動進(jìn)行了預(yù)測,提升了RANS 求解的精度.由于TBNN 基于Pope的有效黏度假設(shè)構(gòu)建,而該假設(shè)是Pope 針對絕大部分不可壓縮流提出的,因此TBNN 的研究工作往往集中在低速的領(lǐng)域.

本文在蒙特利爾大學(xué)和谷歌開發(fā)的Theano 深度學(xué)習(xí)框架(https://pypi.org/project/Theano/)上完成了神經(jīng)網(wǎng)絡(luò)的編譯.以Pope[5]修正的有效黏度模型作為理論基礎(chǔ),并基于Ling 等[7]搭建的TBNN 內(nèi)核構(gòu)建了神經(jīng)網(wǎng)絡(luò)模型.本文通過TBNN 構(gòu)建了從雷諾平均方程(RANS)湍流模型的湍動能、湍流耗散率和速度梯度到高精度數(shù)值解的雷諾應(yīng)力各向異性張量的映射.除了對低速槽道流、低速NACA0012翼型進(jìn)行了預(yù)測以外,還將TBNN 的應(yīng)用范圍從低速拓展到了高超聲速,較精確地預(yù)測了馬赫6 的平板邊界層的雷諾應(yīng)力各向異性張量,并通過一組低速平板進(jìn)行了模型泛化能力的進(jìn)一步驗(yàn)證.驗(yàn)證了TBNN 在寬速域流動下對雷諾應(yīng)力各向異性的預(yù)測能力.對于本文涉及的3 個(gè)TBNN 模型,經(jīng)由槽道流訓(xùn)練得的模型標(biāo)記為TBNN-C,翼型對應(yīng)模型為TBNN-N,高超聲速平板對應(yīng)模型為TBNN-H.

1 計(jì)算方法

1.1 有效黏度假設(shè)

時(shí)均形式的N-S 方程為

主應(yīng)力和切應(yīng)力的不同取決于坐標(biāo)系的選擇,而根據(jù)雷諾應(yīng)力內(nèi)在的差別可以將其區(qū)分為各向同性和各向異性

式中 kδij為各向同性張量,aij為各向異性張量.

為封閉RANS 方程,需要給出雷諾應(yīng)力與平均流場間的關(guān)系,這一封閉的過程可基于有效黏度假設(shè)完成.最早的有效黏度假設(shè)由Boussinesq 提出,該假設(shè)也是使用最為廣泛的有效黏度假設(shè).基于該假設(shè),雷諾應(yīng)力可被定義為

式中k 為湍動能,δij為克羅內(nèi)克符號,μeff為有效黏度,Ui,j和Uj,i為速度梯度.發(fā)現(xiàn)Boussinesq 假設(shè)是基于各向同性假設(shè)的,這也使得該假設(shè)無法準(zhǔn)確地捕捉各向異性.

為了更好的捕獲雷諾應(yīng)力各向異性張量,Pope[5]進(jìn)一步修正了有效黏度模型.將各向異性張量通過湍動能進(jìn)行歸一化

Pope[5]根據(jù)Caley–Hamilton 理論推導(dǎo)了歸一化的各向異性張量b 與基張量之間的本構(gòu)關(guān)系的一般形式

式中 T(1),T(2),···,T(10)為基張量,λ1,λ2···,λ5為張量不變量,二者都是無量綱化后的張量S 和R 相關(guān)的函數(shù).S 為平均應(yīng)變率張量,R 為平均旋轉(zhuǎn)率張量.不同于三維流場中的10 個(gè)基張量,該本構(gòu)關(guān)系在二維流場中僅需要4 個(gè)基張量,具體表達(dá)為

對于張量不變量的構(gòu)造,Johnson[15]枚舉了對稱和反對稱張量的7 個(gè)相關(guān)不變量:Tr(S),Tr(S2),Tr(S3),Tr(R2),Tr(R2S),Tr(R2S2)和Tr(R2SRS2).而Pope 將本構(gòu)關(guān)系中的張量不變量確定為

無量綱化后的張量S 和R 表達(dá)式為

由于張量不變量和基張量的存在,任何滿足該有效黏度模型的本構(gòu)關(guān)系的張量b 都會自動滿足伽利略不變性.

1.2 深度神經(jīng)網(wǎng)絡(luò)

全連接神經(jīng)網(wǎng)絡(luò)的大致結(jié)構(gòu)如圖1 所示.每一個(gè)圓圈代表神經(jīng)元,其中包含兩步計(jì)算:第一步為輸入向量X 和權(quán)值向量W 及偏置b 的線性運(yùn)算a=WX+b,第二步為激活函數(shù)y=h(a)的非線性運(yùn)算.常用的激活函數(shù)h(x) 包括Tanh,Sigmoid,ReLU和階躍函數(shù)等.本文沿用了TBNN 中采用了leaky ReLU 激活函數(shù)[6].

圖1 深度神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)Fig.1 Structure of deep neural network

通常采用帶有梯度下降的反向傳播方法來訓(xùn)練神經(jīng)網(wǎng)絡(luò),本文選擇了Adam 方法.該方法通過迭代調(diào)整神經(jīng)網(wǎng)絡(luò)中各節(jié)點(diǎn)的權(quán)重,從而在各向異性張量的預(yù)測值和真實(shí)值之間提供最低的均方誤差,以使模型更適合訓(xùn)練數(shù)據(jù).此外隱藏層的層數(shù)、每層的寬度也會顯著影響神經(jīng)網(wǎng)絡(luò)的性能.神經(jīng)網(wǎng)絡(luò)的層數(shù)和節(jié)點(diǎn)數(shù)較多時(shí)能容納更復(fù)雜的數(shù)據(jù),但也容易造成過擬合.

1.3 張量基神經(jīng)網(wǎng)絡(luò)

TBNN 目標(biāo)是確定式(6)中的g(n)(λ1,λ2,···,λ5)函數(shù).一旦確定了函數(shù),則可以使用公式(6)求解張量b.對于機(jī)器學(xué)習(xí)而言,將相同的流動在不同坐標(biāo)系下的流場用于機(jī)器學(xué)習(xí)時(shí),機(jī)器學(xué)習(xí)模型可能會產(chǎn)生不同的預(yù)測,因此在使用TBNN 時(shí)保證神經(jīng)網(wǎng)絡(luò)的伽利略不變性是必要的.

基于Pope 的有效黏度模型,將5 個(gè)張量不變量λ1,λ2,···,λ5作為TBNN 的輸入變量,標(biāo)量函數(shù)g(n)作為輸出變量,從而避免在神經(jīng)網(wǎng)絡(luò)中引入張量運(yùn)算,保證了其伽利略不變性.在完成了TBNN 模型的構(gòu)建后,結(jié)合基張量進(jìn)行式(7)的運(yùn)算,即可得到歸一化的雷諾應(yīng)力各向異性張量b.TBNN 的結(jié)構(gòu)圖如圖2 所示.

圖2 張量基神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)Fig.2 Structure of tensor basis neural network (TBNN)

此外為了衡量訓(xùn)練效果,引入均方根誤差作為損失函數(shù)

學(xué)習(xí)率采用學(xué)習(xí)率衰減(learning rate decay),初始學(xué)習(xí)率為0.01,最小學(xué)習(xí)率為1 × 10-6.

并根據(jù)Banerjee 等[16]提出的湍流不變量與特征值方法的各向異性特性的表述,對b 和其特征值ξ1≥ξ2≥ξ3添加約束,從而對預(yù)測好的各向異性進(jìn)行后處理

2 數(shù)據(jù)集

Ling 等[7]的工作對多種基礎(chǔ)流動進(jìn)行了研究,諸如管流(Reb=3500)、周期山流動等.張珍等[14]使用EVNN 和TBNN 的組合神經(jīng)網(wǎng)絡(luò)進(jìn)一步研究了周期山流動.為了驗(yàn)證基于Pope 有效黏度模型的TBNN在寬速域和小樣本的前提下的適用范圍,本文選取了低速槽道流、低速翼型以及高超聲速平板三類算例展開研究.由于獲得的高精度數(shù)據(jù)均對應(yīng)二維,本文采用翼型、槽道流以及高超平板的二維算例,以精確匹配高精度DNS 和LES 數(shù)據(jù)的二維工況.

針對低速工況,選擇NACA0012 翼型以及槽道流作為研究算例.槽道流的DNS 數(shù)據(jù)來自Moser等[17-19],用于TBNN 訓(xùn)練的為Reτ=395 工況下湍流充分發(fā)展區(qū)域的一條截線上136 個(gè)點(diǎn)的流場數(shù)據(jù),Reτ=590 工況對應(yīng)的136 個(gè)點(diǎn)將作為外推的測試集驗(yàn)證TBNN 預(yù)測能力;針對Reτ=590,Krogstad 等[20]在風(fēng)洞實(shí)驗(yàn)中通過熱線風(fēng)速儀測量了該工況的雷諾應(yīng)力各向異性張量,相關(guān)實(shí)驗(yàn)中的18 個(gè)流場點(diǎn)可作為測試集的驗(yàn)證,與TBNN 預(yù)測結(jié)果和DNS 結(jié)果形成充分的對比,進(jìn)一步衡量TBNN 預(yù)測能力.翼型來自Vinuesa 等[21-22]通過LES 求解的NACA0012 翼型,工況為Rec=4×105,AoA=0°,不可壓縮.數(shù)據(jù)包含了流場上翼面的1800 個(gè)數(shù)據(jù)點(diǎn),其中80%的點(diǎn)作為訓(xùn)練數(shù)據(jù),余下20%的點(diǎn)將通過散點(diǎn)圖和云圖的形式對TBNN 的預(yù)測效果進(jìn)行驗(yàn)證,對比測試中尚不包含實(shí)驗(yàn)結(jié)果.

高超聲速流動涉及了強(qiáng)激波、強(qiáng)逆壓梯度、強(qiáng)壓縮效應(yīng)等諸多因素[23],因此Pope[5]針對大部分不可壓縮流動提出的本構(gòu)關(guān)系在高超是否適用需要重新判斷.本文選擇一組Ma=6 的平板作為研究算例.該算例來自Sun 等[24]最新計(jì)算的高超聲速平板,以邊界層內(nèi)的流場為研究對象,并于該平板上選取了8 條截線共計(jì)1540 個(gè)流場點(diǎn)進(jìn)行研究,其中的6 條截線用于訓(xùn)練,相關(guān)DNS 計(jì)算方法[25]以及平板計(jì)算結(jié)果[24]見以下文章.而相較于低速槽道流的實(shí)驗(yàn),高超聲速工況下的雷諾應(yīng)力難以使用熱線風(fēng)速儀進(jìn)行測量.其它測量手段如納米示蹤平面激光散射技術(shù)(NPLS)[26]可以獲得超聲速、高超聲速條件下的雷諾應(yīng)力分布,但由于本文需要進(jìn)行邊界層內(nèi)不同站位法向方向雷諾應(yīng)力各向異性分量的定量對比,故僅采用DNS 計(jì)算結(jié)果與本文計(jì)算結(jié)果進(jìn)行對照.

本文中用于TBNN 的訓(xùn)練集包含兩部分.第一部分為RANS 求解所得結(jié)果,其中湍流模型采用k-ε,SST 模型進(jìn)行求解.根據(jù)流場的湍動能、湍流耗散率和速度梯度構(gòu)造無量綱化后的張量S 和R,進(jìn)而計(jì)算出5 個(gè)張量不變量作為神經(jīng)網(wǎng)絡(luò)輸入.第二部分為高精度的雷諾應(yīng)力數(shù)據(jù),經(jīng)求解得雷諾應(yīng)力各向異性張量,作為TBNN 的輸出.在經(jīng)過訓(xùn)練集完成訓(xùn)練后,便可運(yùn)用TBNN 模型對篩選出的測試集進(jìn)行預(yù)測.

在完成RANS 計(jì)算后,可基于RANS 求解結(jié)果并根據(jù)公式(12)對雷諾應(yīng)力各向異性張量進(jìn)行初步的反推[7]

通過此式亦可對Boussinesq 渦黏假設(shè)在各向異性的預(yù)測能力進(jìn)行初步的判斷.

3 數(shù)值模擬結(jié)果

本文主要研究雷諾應(yīng)力各向異性張量的預(yù)測,需要完成RANS 求解、訓(xùn)練集構(gòu)建以及TBNN 的訓(xùn)練和預(yù)測.

3.1 低速槽道流

低速槽道流的RANS 計(jì)算由k-ε 模型完成.選擇Reτ=395 的工況作為訓(xùn)練集,并將訓(xùn)練好的TBNN 模型命名為TBNN-C.將TBNN-C 用于Reτ=590 的工況進(jìn)行外推算例的預(yù)測.開源的DNS 數(shù)據(jù)中[17-19]給出了完全發(fā)展后的一個(gè)截線,共有136 個(gè)流場點(diǎn);兩種工況在RANS 計(jì)算的網(wǎng)格細(xì)節(jié)如表1 所示.

表1 中Lx為計(jì)算域在x 方向的長度,δ 為槽道的半高,Ly=2δ 為y 方向的長度.Δx+為x 方向兩個(gè)網(wǎng)格點(diǎn)的距離差,Δyc+為y 方向中心位置兩個(gè)網(wǎng)格點(diǎn)的距離差.此處針對RANS 所用網(wǎng)格的進(jìn)行了網(wǎng)格無關(guān)性驗(yàn)證.在確保了網(wǎng)格不會對求解精度造成影響后,選取其中湍流充分發(fā)展的位置作截線,與DNS 數(shù)據(jù)中的流場點(diǎn)進(jìn)行匹配.由于槽道流的對稱性,只取下半部分作為研究對象,網(wǎng)格尺寸為1024 ×257.在保證RANS 求解的速度型與DNS 相近之后,將RANS 求得的速度梯度、湍動能和湍流耗散率以及DNS 中的雷諾應(yīng)力張量輸入到TBNN-C 中.隱藏層的層數(shù)分別設(shè)置為2,4 和6 層并對訓(xùn)練結(jié)果的均方根誤差進(jìn)行對比,最終選擇6 層,每層20 個(gè)節(jié)點(diǎn).

表1 槽道流RANS 計(jì)算的網(wǎng)格信息Table 1 Mesh information for channel flow in RANS

該工況的實(shí)驗(yàn)由Krogstad 等[20]在回流式風(fēng)洞中完成,通過兩個(gè)平行的光滑平板形成了槽道,試驗(yàn)段長5 m,入口面積1.35 m × 0.10 m,通過熱線風(fēng)速儀測量速度脈動,從而獲得雷諾應(yīng)力.圖3 為Reτ=590 的工況下,雷諾主應(yīng)力的各分量隨著y+的變化.y+的計(jì)算公式如式(14) 所示,其中 uτ為摩擦速度(即湍流壁面附近的黏性速度量級),v 為運(yùn)動黏度

圖3 Reτ=590 雷諾應(yīng)力各向異性分量隨y+值的分布Fig.3 Reynolds normal stress anisotropy components vs y + at Reτ=590

區(qū)別于雷諾應(yīng)力張量,從圖3 結(jié)果我們可以發(fā)現(xiàn),TBNN-C 的預(yù)測結(jié)果在y+較小的區(qū)域(即y+＜5的黏性子層區(qū)域)有較大偏離,但在過渡子層以及完全湍流區(qū)域TBNN-C 的預(yù)測效果與DNS 以及風(fēng)洞實(shí)驗(yàn)的結(jié)果相近,預(yù)測效果良好,能夠揭示雷諾正應(yīng)力的各向異性部分在槽道流中的發(fā)展規(guī)律.在y+＞5的范圍內(nèi),TBNN-C 預(yù)測的bii與DNS 結(jié)果的整體誤差僅有10%左右,與實(shí)驗(yàn)相同點(diǎn)位上TBNN-C 的預(yù)測誤差均不超過10%;而與TBNN 相比,RANS 的預(yù)測結(jié)果存在量級上的差異,也說明了基于Boussinesq渦黏假設(shè)的湍流k-ε 湍流模型無法準(zhǔn)確捕捉雷諾應(yīng)力各向異性張量.

此外,針對槽道流的TBNN-C 模型僅根據(jù)Reτ=395 的136 個(gè)流場點(diǎn)進(jìn)行訓(xùn)練,在樣本極小的情況下,訓(xùn)練完成的TBNN-C 模型針對Reτ=590 的不同流場外推預(yù)測依舊展現(xiàn)出極佳的泛化能力.這也說明對于槽道流這種外形簡單的研究對象,Pope[5]提出的本構(gòu)關(guān)系能夠充分詮釋其流動機(jī)理,從而極大的提升了TBNN 模型的泛化能力.

3.2 低速NACA0012 翼型

本文所研究的NACA 翼型為NACA0012.由于NACA0012 在零度攻角下的對稱性,僅研究其上翼面.大渦模擬方法(LES)作為研究復(fù)雜湍流問題的重要工具,其求解結(jié)果擁有較高的精度[27-28],可作為神經(jīng)網(wǎng)絡(luò)的訓(xùn)練標(biāo)簽.在開源的LES 數(shù)據(jù)中,共有1800 個(gè)流場點(diǎn)[21-22],皆位于機(jī)翼上翼面近壁區(qū).翼型RANS 計(jì)算選擇的湍流模型為SST 模型.在完成RANS 計(jì)算后,將其中的流場點(diǎn)提取并與LES 流場點(diǎn)一一匹配.為驗(yàn)證RANS 計(jì)算精度,圖4 中對比了RANS 解得的速度場云圖以及LES 的1800 個(gè)數(shù)據(jù)點(diǎn)插值所得速度場云圖,RANS 的計(jì)算結(jié)果能夠初步的反映LES 結(jié)果中的速度變化趨勢.

圖4 RANS 與LES 的速度場對比Fig.4 Comparison of velocity field for RANS and LES

圖4 RANS 與LES 的速度場對比(續(xù))Fig.4 Comparison of velocity field for RANS and LES (continued)

在LES 數(shù)據(jù)中選擇1440 個(gè)流場點(diǎn)作為訓(xùn)練集對TBNN 進(jìn)行訓(xùn)練,訓(xùn)練完成后的TBNN 模型命名為TBNN-N.將額外的360 個(gè)流場點(diǎn)作為測試集輸入到TBNN-N 模型中,檢驗(yàn)神經(jīng)網(wǎng)絡(luò)對于翼型的預(yù)測效果.測試集測點(diǎn)分布如圖5 所示.

圖5 NACA0012 測試集流場點(diǎn)分布Fig.5 Distribution of NACA0012 test data set

在對隱藏層的層數(shù)進(jìn)行了測試后,選擇三層隱藏層,每層20 個(gè)節(jié)點(diǎn)進(jìn)行訓(xùn)練.訓(xùn)練步數(shù)為8000 步,最終訓(xùn)練完成時(shí)訓(xùn)練集的均方根誤差為0.052,訓(xùn)練的樣本數(shù)為1440,樣本數(shù)量較小.圖6 和圖7 給出了測試集的360 個(gè)流場點(diǎn)的各向異性分量b11,b22的分布,并選取關(guān)鍵區(qū)域的流場點(diǎn)插值繪制云圖.插值方法為Matlab 中的V4 雙調(diào)和樣條插值(biharmonic spline interpolation).

圖6 雷諾正應(yīng)力各向異性分量b11 結(jié)果Fig.6 The result of Reynolds normal stress anisotropy component b11

圖7 歸一化雷諾正應(yīng)力各向異性分量b22 云圖分布Fig.7 The result of normalized Reynolds normal stress anisotropy component b22

對于NACA0012 翼型,基于RANS 求解的平均速度場幾乎完全無法捕捉雷諾應(yīng)力各向異性分量,而基于RANS 結(jié)果預(yù)測的TBNN-N 各向異性分布在個(gè)別流場點(diǎn)的預(yù)測較差,但在重點(diǎn)關(guān)注區(qū)域TBNN-N 能夠較好地預(yù)測雷諾應(yīng)力各向異性.同樣在低速工況下,TBNN-N 針對翼型的預(yù)測相比于槽道流精度有所下降,但在小樣本的前提下TBNNN 取得的預(yù)測效果尚可接受,較RANS 結(jié)果顯著提升.要對翼型進(jìn)行精準(zhǔn)預(yù)測,可考慮進(jìn)一步擴(kuò)大樣本范圍,以提升TBNN 的預(yù)測精度.

3.3 高超聲速平板

針對低速工況的模擬方法面對高超聲速往往會產(chǎn)生不適用性[29-30],故通過高超聲速平板來衡量TBNN 在高超聲速下預(yù)測的精確度.當(dāng)前所擁有的平板DNS 數(shù)據(jù)為二維算例[24],來流參數(shù)中馬赫數(shù)為6,單位雷諾數(shù)為12 000,在湍流區(qū)域共有8 條截線.選取截線位置邊界層內(nèi)的數(shù)據(jù)用于本算例,其中B11 距平板前緣點(diǎn)51.68%處,各截線位置以及邊界層厚度如圖8 所示.

圖8 平板DNS 數(shù)據(jù)各截線位置以及邊界層厚度Fig.8 Position of transversals and boundary layer thickness of DNS data

在8 條截線中,選擇B13 和B16 作為測試集,其余6 組截線數(shù)據(jù)作為訓(xùn)練集對TBNN 進(jìn)行訓(xùn)練,訓(xùn)練完成后的模型命名為TBNN-H.訓(xùn)練集中共有1155 個(gè)流場點(diǎn),用于預(yù)測的B13 和B16 截線分別有187 和199 個(gè)流場點(diǎn).

高超聲速平板的RANS 計(jì)算由Fluent 完成,使用可壓縮k-ε 模型作為湍流模型.根據(jù)邊界層厚度在RANS 結(jié)果中完成B11 到B18 的各截線位置的匹配.RANS 計(jì)算的網(wǎng)格尺寸為3869 × 320,在完成網(wǎng)格無關(guān)性驗(yàn)證后,保證相同站位RANS 與DNS 邊界層厚度一致,對比結(jié)果如圖9 所示.此外,如圖9還給出了邊界層內(nèi)近壁區(qū)域的網(wǎng)格加密.

圖9 Ma6 平板RANS 邊界層網(wǎng)格以及邊界層厚度Fig.9 The mesh for RANS and boundary layer thickness

根據(jù)DNS 計(jì)算結(jié)果,截線B13 的邊界層厚度δ=6.77 mm,B16 的邊界層厚度δ=7.91 mm.在邊界層厚度相同的情況下,相比于DNS 計(jì)算結(jié)果,RANS的速度發(fā)展較慢.以此RANS 結(jié)果為基礎(chǔ)的TBNN-H依舊取得了極佳的預(yù)測效果,如圖10 和圖11 所示.(由于流場點(diǎn)過多,進(jìn)行數(shù)據(jù)點(diǎn)顯示控制以美觀曲線,每兩個(gè)點(diǎn)顯示一個(gè)點(diǎn))

圖10 截線B13 的雷諾主應(yīng)力各向異性分量的分布Fig.10 Reynolds normal stress anisotropy components on transversals B13

圖11 截線B16 的雷諾主應(yīng)力各向異性分量的分布Fig.11 Reynolds normal stress anisotropy components on transversals B16

針對高超平板的TBNN-H 使用了4 個(gè)隱藏層,每層18 個(gè)神經(jīng)元.在訓(xùn)練完成后,訓(xùn)練集的均方根誤差(root mean square error,RMSE)維持在0.01 左右,而訓(xùn)練好的TBNN-H 對于測試集的預(yù)測表現(xiàn)良好,RMSE 約為0.015 左右.從結(jié)果中可以發(fā)現(xiàn),在y+較小的區(qū)域(即y+＜5 的黏性子層區(qū)域)TBNN-H 的預(yù)測結(jié)果有一定偏離;與槽道流相同,TBNN-H 在過渡子層以及完全湍流區(qū)域TBNN-H 的預(yù)測效果與DNS 的結(jié)果相近,預(yù)測效果良好.這也說明了TBNN 的預(yù)測在一定程度上可能會受到k-ε 模型的固有限制,因此在湍流邊界層的黏性子層的預(yù)測效果較差.

此外,從圖10 和圖11 的結(jié)果中我們可以發(fā)現(xiàn)TBNN-H 對于高超聲速平板邊界層依舊具有良好的預(yù)測能力,盡管Pope 的本構(gòu)關(guān)系針對不可壓縮流動提出,但基于該本構(gòu)關(guān)系構(gòu)造的TBNN-H 依舊能夠?qū)?qiáng)壓縮性的高超聲速平板流動雷諾應(yīng)力各向異性張量進(jìn)行預(yù)測.而這也意味著通過高超聲速樣本構(gòu)建的TBNN 模型可以用來提高可壓縮模型的針對性,為高超聲速湍流模型的定制化提供方法基礎(chǔ)[1].

3.4 模型泛化能力驗(yàn)證

神經(jīng)網(wǎng)絡(luò)模型的泛化能力至關(guān)重要,本節(jié)將對訓(xùn)練好的TBNN 模型能否應(yīng)用于與訓(xùn)練算例不同的算例展開討論.在3.1 節(jié)中,TBNN-C 模型由Reτ=395 的低速槽道流完成訓(xùn)練,并對Reτ=590 的槽道流工況進(jìn)行了較好的預(yù)測,初步驗(yàn)證了TBNN 的外推能力.為了進(jìn)一步驗(yàn)證其泛化能力,本節(jié)選取了一組Reθ=1100 的低速平板作為驗(yàn)證算例,該平板的DNS 結(jié)果由Jimenez 等[31]計(jì)算.將Reτ=395 的低速槽道流作為訓(xùn)練數(shù)據(jù),通過訓(xùn)練好的TBNN-C 模型對該低速平板進(jìn)行預(yù)測.

驗(yàn)證算例的開源DNS 數(shù)據(jù)[31]給出了邊界層厚度為δ99=2.756 8 的對應(yīng)截線上的數(shù)據(jù)點(diǎn),其中邊界層內(nèi)的流場點(diǎn)共有130 個(gè).該算例的RANS 計(jì)算由standard k-ε 模型完成,網(wǎng)格參數(shù)見表2.

表2 驗(yàn)證算例RANS 計(jì)算的網(wǎng)格信息Table 2 Mesh information for validation case in RANS

表中Lx為計(jì)算域在x 方向的長度,Ly為y 方向的長度,Δx+為x 方向兩個(gè)網(wǎng)格點(diǎn)的距離差.RANS計(jì)算的y 方向的網(wǎng)格分布與DNS 計(jì)算保持一致.此處針對RANS 所用網(wǎng)格的進(jìn)行了網(wǎng)格無關(guān)性驗(yàn)證.在確保了網(wǎng)格不會對求解精度造成影響后,選取其中湍流充分發(fā)展的位置作截線,與DNS 數(shù)據(jù)中的流場點(diǎn)進(jìn)行匹配以便于最終進(jìn)行對比.

將RANS 算得的平板邊界層的湍動能、湍流耗散率以及速度梯度輸入到TBNN-C 模型中,得到預(yù)測的結(jié)果.同時(shí)將預(yù)測結(jié)果與DNS 結(jié)果以及通過公式(13)解得的RANS 結(jié)果進(jìn)行對比,如圖12 所示.

圖12 Reθ=1100 平板邊界層的雷諾主應(yīng)力各向異性分量的分布Fig.12 Reynolds normal stress anisotropy components on Reθ=1100 boundary layer

相比于3.1 節(jié)的結(jié)果,基于槽道流的136 個(gè)流場點(diǎn)訓(xùn)練的TBNN-C 模型對于不同工況的平板邊界層的預(yù)測精度盡管稍有下降,但與DNS 結(jié)果相比誤差仍不超過20%,依舊較為準(zhǔn)確的預(yù)測了雷諾應(yīng)力各向異性張量,驗(yàn)證了TBNN 模型的泛化能力.而使用上述模型預(yù)測翼型以及高超聲速平板邊界層,結(jié)果顯示模型的預(yù)測的精度大幅下滑,在一定程度上說明了模型性能對于訓(xùn)練數(shù)據(jù)的依賴.盡管TBNN 是基于Pope 的本構(gòu)關(guān)系構(gòu)建,但依舊會受到神經(jīng)網(wǎng)絡(luò)特性的固有限制.

4 結(jié)論

本文寬速域下TBNN 對雷諾應(yīng)力各向異性張量的預(yù)測結(jié)論如下.

(1) 無論是低速還是高超聲速,基于Boussinesq有效黏度假設(shè)的RANS 模型均難以準(zhǔn)確地捕捉雷諾應(yīng)力各向異性張量.

(2) 基于Pope[5]提出的有效黏度假設(shè)構(gòu)造出的TBNN-C 對于槽道流的預(yù)測結(jié)果較好,在訓(xùn)練樣本僅有136 個(gè)流場點(diǎn)的極小樣本的情況下,依舊可以準(zhǔn)確地預(yù)測不同雷諾數(shù)下槽道流的雷諾應(yīng)力各向異性張量.TBNN 預(yù)測結(jié)果與DNS 誤差在10%左右,部分點(diǎn)位與風(fēng)洞實(shí)驗(yàn)結(jié)果誤差亦在10%以內(nèi),泛化能力較好.

(3) TBNN-N 對于NACA0012 翼型的預(yù)測在個(gè)別流場點(diǎn)存在偏差.在選取了流場中1440 個(gè)點(diǎn)作為訓(xùn)練集之后,TBNN 對于同流場額外的360 個(gè)點(diǎn)的預(yù)測結(jié)果尚可.在小樣本的訓(xùn)練前提下,TBNN-N 能夠?qū)ACA0012 的關(guān)鍵區(qū)域進(jìn)行較為準(zhǔn)確的預(yù)測.與槽道流相比,TBNN 對于復(fù)雜流場需要增大樣本量以提升準(zhǔn)確性.

(4) 對于高超聲速平板,在以湍流區(qū)域部分位置進(jìn)行訓(xùn)練后,TBNN-H 能夠較好地給出流場中其他位置的雷諾應(yīng)力各向異性張量分布.盡管Pope[5]有效黏度假設(shè)是針對不可壓縮流提出的,TBNN-H 依舊可以在小樣本的前提下對高超聲速平板邊界層內(nèi)的流場進(jìn)行精準(zhǔn)的預(yù)測.

(5) 在低速槽道流和低/高超聲速平板中,TBNN在黏性子層(即y+＜5 區(qū)域)表現(xiàn)較差,猜測原因可能是訓(xùn)練集的RANS 部分來自于k-ε 模型求解,TBNN 的預(yù)測能力會受到RANS 結(jié)果的固有限制.

(6) 低速槽道流訓(xùn)練的TBNN-C 模型能夠較為精確的預(yù)測低速平板算例,模型的泛化能力得到了驗(yàn)證.而相同模型對于翼型和高超聲速平板的預(yù)測出現(xiàn)的預(yù)測精度下滑則在一定程度上說明了模型性能對于訓(xùn)練數(shù)據(jù)的依賴.

本文在寬速域下,對TBNN 預(yù)測能力進(jìn)行了充分的驗(yàn)證.相比于Pope 針對不可壓縮工況提出的本構(gòu)關(guān)系,神經(jīng)網(wǎng)絡(luò)憑借出色的從數(shù)據(jù)中提取信息的能力在映射關(guān)系中學(xué)習(xí)到了可壓縮相關(guān)的信息,使得基于Pope 本構(gòu)關(guān)系的TBNN 能夠?qū)Ω叱曀俟r進(jìn)行較好的預(yù)測.盡管神經(jīng)網(wǎng)絡(luò)的“黑匣子”性質(zhì)使得TBNN 難以像本構(gòu)關(guān)系一樣具有可解釋性,但本工作對于高超聲速的雷諾應(yīng)力各向異性的求解仍具有參考價(jià)值.

在下一步的工作中,將針對更多的高超聲速復(fù)雜工況展開研究,同時(shí)增加三維的相關(guān)預(yù)測.在進(jìn)一步完成TBNN 在復(fù)雜高超聲速工況的適用性驗(yàn)證后,將通過單向耦合的方式將訓(xùn)練好的TBNN 模型與可壓縮湍流模型以及N-S 求解器耦合,嘗試改善對于流動的預(yù)測精度.而耦合對應(yīng)的穩(wěn)定性與收斂性的問題將是研究的重點(diǎn)與難點(diǎn).此外,還將通過在神經(jīng)網(wǎng)絡(luò)中添加物理約束,構(gòu)建物理驅(qū)動的張量基神經(jīng)網(wǎng)絡(luò)(physics informed tensor based neural network),進(jìn)一步增強(qiáng)TBNN 在小樣本的泛化能力.

致 謝

感謝Julia Ling 于https://github.com/tbnn/tbnn開源的TBNN 內(nèi)核代碼;感謝Ricardo Vinuesa 提供的Naca0012 翼型LES 數(shù)據(jù)和幫助;感謝Robert D.Moser 提供的槽道流DNS 數(shù)據(jù)以及Jimenez 提供的低速平板DNS 數(shù)據(jù).