超聲速湍流密度脈動預(yù)測的神經(jīng)網(wǎng)絡(luò)方法

王正魁，靳旭紅，朱志斌，程曉麗

中國航天空氣動力技術(shù)研究院，北京 100074

紅外制導(dǎo)導(dǎo)彈在大氣層內(nèi)飛行時，激波和湍流引起的流場介質(zhì)不均勻分布，使光線產(chǎn)生折射和衍射等現(xiàn)象，從而導(dǎo)致探測圖像發(fā)生畸變，這種現(xiàn)象被稱為氣動光學(xué)效應(yīng)。氣動光學(xué)效應(yīng)會影響光學(xué)導(dǎo)引頭對目標(biāo)的追蹤和識別能力，因此湍流流場及其氣動光學(xué)效應(yīng)研究對精確制導(dǎo)導(dǎo)彈的研究具有重要意義。

國外已有學(xué)者和團(tuán)隊利用此類方法評估或改進(jìn)RANS模型。針對RANS方法中不同湍流模型對特定物理問題的適用性，有文獻(xiàn)利用監(jiān)督學(xué)習(xí)方法評估RANS方法所得結(jié)果的可靠性。其中第1種思路是通過評估所采取假設(shè)的準(zhǔn)確程度來預(yù)測RANS方法所得結(jié)果的可靠度[7-9]；第2種思路是直接分析雷諾應(yīng)力張量的誤差與平均流動特征間的關(guān)系[10]；第3種思路是通過評估輸運方程的準(zhǔn)確性來評估RANS方法的可靠度[11-13]。

Tracey等[14]利用平板邊界層的算例訓(xùn)練神經(jīng)網(wǎng)絡(luò)，使其掌握一方程Spalart-Allmaras模型中輸運方程的規(guī)律。發(fā)現(xiàn)神經(jīng)網(wǎng)絡(luò)模型可預(yù)測平板邊界層和翼型繞流，但對槽道流的預(yù)測效果較差。可見，神經(jīng)網(wǎng)絡(luò)方法在相似流場上也可體現(xiàn)出一定的泛化能力。此外，文獻(xiàn)[15-17]將k-ε模型中的常數(shù)修改為可變系數(shù)，并用神經(jīng)網(wǎng)絡(luò)或貝葉斯方法來減小誤差,發(fā)現(xiàn)系數(shù)具有大幅度的變化，改進(jìn)模型更為準(zhǔn)確。Duraisamy等[18-19]將Spalart-Allmaras模型中耗散項乘以一個函數(shù)。利用貝葉斯方法重構(gòu)出該函數(shù)的空間分布，計算域中較高的函數(shù)值意味著較高耗散率。改進(jìn)模型明顯提高了非平衡湍流邊界層中摩阻系數(shù)的計算精度。采用Sγ代表間歇因子γ的輸運方程中生成項與耗散項之差，從RANS方法得到的參考流場中提取其分布，并利用標(biāo)準(zhǔn)的爬山算法選擇特征量。然后，利用高斯過程和神經(jīng)網(wǎng)絡(luò)方法為Sγ建立模型，得到較好的預(yù)測效果。

考慮湍流的復(fù)雜性、方法的計算量和性能，本文采用神經(jīng)網(wǎng)絡(luò)方法進(jìn)行函數(shù)回歸，并對湍流密度脈動的均方值進(jìn)行建模。首先，利用DNS方法模擬超聲速壁湍流，并統(tǒng)計得到密度脈動均方值，密度和速度等平均場。然后，分別以密度脈動均方值和平均場的特征量作為樣本的輸出和輸入，對神經(jīng)網(wǎng)絡(luò)進(jìn)行訓(xùn)練和測試。因為密度脈動的均方值的影響因素眾多，為提高模型的泛化能力，采用特征選擇方法來篩選出合適的特征量作為輸入。

1 研究方法

1.1 神經(jīng)網(wǎng)絡(luò)模型

為度量模型的性能，先定義性能度量為均方誤差，即損失函數(shù)L(y,d)=(y-d)2，其中y為模型預(yù)測值，d為DNS的統(tǒng)計值。訓(xùn)練的過程即優(yōu)化神經(jīng)網(wǎng)絡(luò)模型，使損失函數(shù)達(dá)到最小。

圖1是雙隱含層神經(jīng)網(wǎng)絡(luò)的示意圖，其中Zm,n代表第m個隱含層的第n個神經(jīng)元，連接不同層神經(jīng)元的箭頭代表突觸權(quán)值。每個神經(jīng)元只能有一個輸出值，但可以有多個輸入值。在第t次訓(xùn)練過程中，先根據(jù)輸入和當(dāng)前權(quán)值計算神經(jīng)元j的誘導(dǎo)局部域，具體表達(dá)式為

(1)

式中：t為當(dāng)前的訓(xùn)練回合數(shù)；wji為神經(jīng)元j與i之間的突觸權(quán)值；bj為神經(jīng)元j的偏置(閾值)；i、j均為神經(jīng)元的編號，且i為j前一隱含層的神經(jīng)元。

計算得到神經(jīng)元的輸出y(t)=φ(V(t))，激活函數(shù)φ代表了神經(jīng)元對輸入信號的響應(yīng)，常用的有Logistic函數(shù)、雙曲正切函數(shù)、線性整流函數(shù)ReLU以及具有自歸一化能力的縮放指數(shù)線性單元SeLU。

文中隱含層選用SeLU，即

(2)

圖1 神經(jīng)網(wǎng)絡(luò)示意圖Fig.1 Schematic diagram of a neural network

輸出層的激活函數(shù)為Logistic函數(shù)，即：φ(V)=(1+exp(-V))-1。

計算隱含層和輸出層神經(jīng)元的輸出時，將前一層的輸出作為其輸入。逐層計算出所有神經(jīng)元的輸出后，根據(jù)實際響應(yīng)y(t)和期望響應(yīng)d(t)調(diào)整突觸權(quán)值w和偏置b。下面給出突觸權(quán)值和偏置的優(yōu)化方法，采用附加動量項的梯度下降法，即：

(3)

式中：κ為梯度下降法的步長；α為動量因子；Δwji(t-1)、Δbj(t-1)為上一次訓(xùn)練中得到的變化量。根據(jù)微分的鏈?zhǔn)揭?guī)則可得

2(yj(t)-dj)φ′(Vj(t))yi(t)

(4)

-2(yj(t)-dj)φ′(Vj(t))

(5)

定義神經(jīng)元j的局域梯度為

(6)

則式(4)～式(5)可簡化為

(7)

將式(7)代入式(3)，并令學(xué)習(xí)率η=2κ可得

Δwji(t)=ηδj(t)yi(t)+αΔwji(t-1)

(8)

Δbj(t)=-ηδj(t)+αΔbj(t-1)

(9)

若神經(jīng)元j為輸出層神經(jīng)元，則可將模型預(yù)測值與訓(xùn)練值之間的誤差代入式(6)，即可直接得到其局域梯度為

δj(t)=(dj(t)-yj(t))φ′(Vj(t))

(10)

若神經(jīng)元j是隱含層神經(jīng)元，不能直接得到其誤差ej(t)，則神經(jīng)元j的局域梯度δj(t)需要由下一層與其直接相連的神經(jīng)元的局域梯度向后遞歸確定，即反向傳播算法。

(11)

為避免梯度下降算法收斂到非最優(yōu)解的局部極小值，文中采用學(xué)習(xí)率退火的方案，即在訓(xùn)練過程中更新學(xué)習(xí)率，即

η=η0/(1+t/τ)

(12)

式中：η0為初始學(xué)習(xí)率；τ為搜索時間常數(shù)。

為降低模型的過擬合風(fēng)險，采用復(fù)雜度正則化約束模型的復(fù)雜度，其思路是修改損失函數(shù)，使其包含標(biāo)準(zhǔn)的性能度量和復(fù)雜度懲罰項，即

L′=L(y,d)+λξ(w)/2

(13)

式中：λ為正則化參數(shù)，代表了性能度量與復(fù)雜度懲罰之間的權(quán)衡，是一個需要使用者調(diào)節(jié)的超參數(shù)；ξ(w)代表了模型的復(fù)雜度，文中取為突觸權(quán)值w的2范數(shù)，則修改后的損失函數(shù)為

(14)

根據(jù)修改后的損失函數(shù)式(14)，其對突觸權(quán)值的偏導(dǎo)數(shù)應(yīng)由式(7)修改為

(15)

突觸權(quán)值的變化量由式(8)修改為

Δwji(t)=η(δj(t)yi(t)-λwji(t-1))+

αΔwji(t-1)

(16)

同時，為提高模型學(xué)習(xí)效率、降低收斂到局部極小值的風(fēng)險，文中又采用批梯度下降法優(yōu)化模型參數(shù)，在每次訓(xùn)練過程中使用的訓(xùn)練樣本數(shù)為批尺寸。本文在訓(xùn)練初期采用較小的批尺寸，訓(xùn)練速度快；隨著模型的收斂，逐漸增加批尺寸，在訓(xùn)練后期得到穩(wěn)定的收斂結(jié)果。

1.2 樣本數(shù)據(jù)獲取方法

考慮到神經(jīng)網(wǎng)絡(luò)的訓(xùn)練需要大量的可靠數(shù)據(jù)，本文采用DNS方法計算得到超聲速平板邊界層的密度脈動、平均場以及湍流統(tǒng)計量，并提取特征量。為獲取較高的空間分辨率，采用加權(quán)基本無振蕩(Weighted Essentially Non-Oscillatory, WENO)的七階格式[20-21]，時間方向采用顯式的二階時間差分格式求解。

固壁邊界條件為速度無滑移、等溫壁；上游采用可壓縮邊界層的層流剖面作為入口邊界，并施加周期性的速度擾動[22]；遠(yuǎn)場和下游為出口邊界條件，并預(yù)置緩沖層來抑制數(shù)值擾動反射的影響；展向采用周期性邊界條件。

計算區(qū)域的流向長度LX=203.7 mm、法向高度LY=7.62 mm、展向?qū)挾萀Z=5.08 mm。網(wǎng)格在X和Z方向是均勻分布的，Y方向網(wǎng)格高度隨著壁面法向距離的增加而增加。網(wǎng)格點數(shù)為2 000×200×200，邊界層首層網(wǎng)格高度為ΔYmin=1.9×10-3mm，湍流區(qū)中0.33

圖2給出了DNS模擬的速度與理論值的比較，可見DNS統(tǒng)計平均得到的流向速度在黏性底層、對數(shù)律層均與壁面律相符得較好，從而驗證了流向速度的可靠性。進(jìn)一步在相同來流條件下，將DNS模擬的密度時均值與文獻(xiàn)[24]的結(jié)果進(jìn)行比較，如圖3所示(圖中ρ∞為來流密度)。可見，本文模擬的密度與文獻(xiàn)結(jié)果相符得較好，進(jìn)一步驗證了訓(xùn)練數(shù)據(jù)所得流場計算方法和結(jié)果的可靠性。

表1 不同算例的來流條件Table 1 Free-stream conditions of different cases

圖2 Van Driest平均速度型的驗證(Re=2.50×107)Fig.2 Validation of Van Driest mean velocity distribution (Re=2.50×107)

流場參數(shù)的時均值通過統(tǒng)計平均獲取，即

(17)

式中：M為統(tǒng)計的時間樣本數(shù)；θi為流場中速度、溫度、密度等瞬時量。

(18)

(19)

(20)

不同的特征量之間可能相差多個數(shù)量級，湍動能和湍流耗散率等數(shù)量級大的特征對神經(jīng)網(wǎng)絡(luò)的影響過大，不利于訓(xùn)練神經(jīng)網(wǎng)絡(luò)，因此文中采用歸一化方法消除特征量之間的量級差異。本文選用MapAverage函數(shù)將輸入特征線性地變化到均值為0、標(biāo)準(zhǔn)差為1的狀態(tài)，即

(21)

圖4給出X=139.7 mm截面處DNS統(tǒng)計的密度脈動隨Y+的變化。可見密度脈動在線性底層(Y+10)較小，在過渡層(10500)迅速減弱。密度脈動在邊界層內(nèi)的分布，與Wang K和Wang M[25]分析不同區(qū)域的光程差作用所得結(jié)論有著十分相似的規(guī)律。

圖4 密度脈動隨Y+的變化Fig.4 Density fluctuation as a function of Y+

2 特征選擇及泛化能力評估

神經(jīng)網(wǎng)絡(luò)模型可能對“見過的”訓(xùn)練樣本預(yù)測得非常好，而對“沒見過”的新樣本預(yù)測誤差很大，即過擬合現(xiàn)象。因此，為合理地評估模型的泛化能力，采用留出法對模型進(jìn)行訓(xùn)練和測試，即將整個樣本空間分為訓(xùn)練集、驗證集和測試集。在訓(xùn)練階段，利用訓(xùn)練集中的樣本(q,d)調(diào)整模型參數(shù)，在每訓(xùn)練一定次數(shù)后利用驗證集的樣本評估模型；在訓(xùn)練階段結(jié)束后，選取在驗證集上性能最佳的模型，并將測試樣本的q輸入到所得模型中，將其預(yù)測的密度脈動y與d進(jìn)行比較，并評估模型的預(yù)測能力。文中不考慮轉(zhuǎn)捩問題，從完全湍流區(qū)(152.4 mm

2.1 平均影響值算法

特征選擇可降低問題和模型的復(fù)雜度，緩解過擬合問題。參考已有模型和先驗知識，選取了19個待選特征量。先以全部待選特征作為輸入量，對模型進(jìn)行訓(xùn)練和測試；接著利用平均影響值算法，得到所有待選特征對模型輸出的影響值，并進(jìn)行排序；再假設(shè)影響值越大的特征越重要，逐一添加影響值大的待選特征，并對模型進(jìn)行訓(xùn)練和測試，直到只剩下一個特征輸入量時為止；確定在篩選特征過程中驗證誤差最小的模型，并將對應(yīng)特征作為選擇出的特征輸入量。本文采用的特征選擇過程如圖5所示。

平均影響值算法可估計特征輸入量對輸出量的影響，即影響值，其絕對值表示影響程度，符號代表相關(guān)的方向。先以全部特征量作為輸入，并訓(xùn)練得到一個模型；然后對其原有輸入值逐個加減10%，并利用訓(xùn)練得到的模型預(yù)測對應(yīng)的輸出；其預(yù)測差值的絕對值在訓(xùn)練集的170 000個樣本上的平均值即為該特征量對輸出的影響值，對其按絕對值由大到小排序，即得到表2。表中v為法向速度，μt為湍流黏性系數(shù)，T為靜溫，ω為比耗散率。為降低特征量出現(xiàn)同時增大或減少的風(fēng)險，隨機(jī)地選擇部分特征量取其相反數(shù)，如表2中-ρ、-T等。

圖5 基于平均影響值算法的特征選擇過程Fig.5 Process of feature selection based on mean impact value algorithm

表2 不同特征的影響值Table 2 Impact values of different features

平均影響值算法得到的影響值的絕對值代表了各個特征的重要程度，在每添加一個特征后對模型進(jìn)行訓(xùn)練和測試，并記錄對應(yīng)的模型測試誤差。下面定義兩個誤差。

1) 均方根誤差，用以衡量模型的準(zhǔn)確性，其表達(dá)式為

(22)

2) 最大誤差，用以衡量模型的穩(wěn)定性，其表達(dá)式為

(23)

實驗過程中，選取不同輸入特征量時，模型的均方根誤差變化如圖6所示。由圖中可見，輸入特征數(shù)為9、11、15和19時均方根誤差比較小，其中誤差最小的模型為15個輸入量，即：-?ρ/?Y,?ρ/?X,-?2ρ/?Y2,u,v,-?u/?X,-?u/?Y,?v/?Y,Y,-T,-ρ,k,ε,μt,μ。對應(yīng)的均方根誤差為1.81×10-4。通過在實驗中調(diào)節(jié)超參數(shù)，確定5個隱含層，從低到高分別有19、15、12、9和5個神經(jīng)元，2范數(shù)的正則化參數(shù)取λ=3×10-8，動量因子α=0.2，初始批尺寸為50，初始學(xué)習(xí)率為η0=0.2。

為驗證神經(jīng)網(wǎng)絡(luò)算法性能、程序的可靠性，先利用單一流場對模型進(jìn)行訓(xùn)練和測試，評估結(jié)果如圖7所示。可以看出，在密度脈動較大的區(qū)域內(nèi)，神經(jīng)網(wǎng)絡(luò)方法較傳統(tǒng)的AMS模型具有更高的預(yù)測精度。即使訓(xùn)練集和測試集不同，神經(jīng)網(wǎng)絡(luò)模型仍然能給出較好的預(yù)測結(jié)果，顯示了一定的泛化能力。且在單一流場中特征選擇對神經(jīng)網(wǎng)絡(luò)模型的預(yù)測效果影響不大，因為過擬合現(xiàn)象不嚴(yán)重。

圖6 不同輸入特征量的模型誤差Fig.6 Errors of models with different input features

圖7 特征選擇對單一流場預(yù)測的影響Fig.7 Effect of feature selection on prediction of a single flow field

模型MSEMEAMS7.87×10-42.27×10-3神經(jīng)網(wǎng)絡(luò),19個輸入量3.61×10-51.10×10-4神經(jīng)網(wǎng)絡(luò),15個輸入量3.85×10-51.17×10-4

表3給出不同模型的測試誤差，在特征選擇前后，神經(jīng)網(wǎng)絡(luò)模型采用19個和15個輸入量時在測試集上的均方根誤差分別為AMS模型的4.59%和4.89%，在剖面上的最大誤差分別為AMS模型的4.85%和5.15%。因此，相對于AMS模型，神經(jīng)網(wǎng)絡(luò)模型的預(yù)測精度顯著提高。

神經(jīng)網(wǎng)絡(luò)模型在單一流場中對“沒見過”的測試樣本給出了很好的預(yù)測結(jié)果。然而，單一流場中的測試樣本與訓(xùn)練樣本之間具有很大的相似性，無法充分評估模型的泛化能力。而且，僅利用單一流場的數(shù)據(jù)訓(xùn)練神經(jīng)網(wǎng)絡(luò)模型，也無法獲得較高的通用性，工程應(yīng)用價值較低。

圖8給出了多流場訓(xùn)練模型的預(yù)測結(jié)果，從圖中可見，神經(jīng)網(wǎng)絡(luò)模型能夠發(fā)現(xiàn)密度脈動的規(guī)律，且預(yù)測效果明顯優(yōu)于AMS模型。但在0.6 mm

圖8 特征選擇對多流場預(yù)測的影響Fig.8 Effect of feature selection on prediction of multiple flow fields

模型MSEME神經(jīng)網(wǎng)絡(luò),19個輸入量1.84×10-46.13×10-4神經(jīng)網(wǎng)絡(luò),15個輸入量1.81×10-45.93×10-4

表4定量地比較特征選擇前后，神經(jīng)網(wǎng)絡(luò)模型的預(yù)測誤差。在特征選擇前后，神經(jīng)網(wǎng)絡(luò)模型采用19個和15個輸入量在測試集上的均方根誤差分別為AMS模型的23.4%和15.9%，在剖面上的最大誤差分別為AMS模型的27.0%和26.1%。因此，相比于AMS模型，神經(jīng)網(wǎng)絡(luò)模型的預(yù)測結(jié)果更為準(zhǔn)確和穩(wěn)定，且特征選擇進(jìn)一步提高了模型的性能和精度。

2.2 加入先驗信息

上述模型中，輸入量Y是一個與流場空間位置有關(guān)的特征量，為提高神經(jīng)網(wǎng)絡(luò)模型對復(fù)雜外形的適應(yīng)性，對Y是否作為特征輸入?yún)?shù)進(jìn)行了進(jìn)一步的評估。

在實驗中重新調(diào)節(jié)超參數(shù)，確定5個隱含層，從低到高分別有13、12、10、10和8個神經(jīng)元，2范數(shù)的正則化參數(shù)取λ=2.7×10-8，動量因子α=0.2，初始批尺寸為50，初始學(xué)習(xí)率為η0=0.2。下面對比是否考慮特征量Y時模型的測試結(jié)果。

圖9分別給出了?ρ/?X、?ρ/?Y、μ、u、v、Y、S中是否考慮參數(shù)Y時，神經(jīng)網(wǎng)絡(luò)模型對所有測試樣本的預(yù)測值隨著真實值的變化情況，可見圖中所有預(yù)測點都在斜率為1的直線附近變化，顯示了神經(jīng)網(wǎng)絡(luò)模型很好的預(yù)測能力。比較是否考慮特征量Y情況下的測試結(jié)果可知，忽略特征Y會使得模型對中等強(qiáng)度密度脈動的預(yù)測誤差增加。

圖10對比了兩個模型對同一個剖面的預(yù)測結(jié)果。從圖9和圖10中可見，加入相關(guān)模型的先驗知識，特征選擇后7個輸入特征的模型可以得到比圖8中15個輸入特征更低的預(yù)測誤差。而忽略特征量Y則降低了模型的預(yù)測精度。

模型定量的模擬誤差見表5，從表中可見，考慮特征量Y(共7個輸入量)時，神經(jīng)網(wǎng)絡(luò)模型在測試集上的均方根誤差和在剖面上的最大誤差分別為1.25×10-4和4.31×10-4，是19個輸入量時模型的67.9%和70.3%。所以，在特征選擇過程中，加入相關(guān)的先驗信息可提高模型性能，在一定程度上避免了模型過擬合現(xiàn)象的發(fā)生。在忽略特征量Y后，神經(jīng)網(wǎng)絡(luò)模型在測試集上的均方根誤差和剖面上的最大誤差分別增加了24.0%和44.8%，特征量Y對模型泛化能力的影響需要在更廣泛的算例中進(jìn)一步評估。

圖9 特征量Y對多流場測試的影響(回歸線)Fig.9 Effect of feature Y on test of multiple flow fields (regression line)

圖10 特征量Y對多流場預(yù)測的影響(剖面圖)Fig.10 Effect of feature Y on predicting multiple flow fields (profile)

模型MSEME神經(jīng)網(wǎng)絡(luò),考慮特征量Y1.25×10-44.31×10-4神經(jīng)網(wǎng)絡(luò),不考慮特征量Y1.55×10-46.24×10-4

以考慮特征量Y時的模型為例，評估輸入量誤差對預(yù)測結(jié)果的影響。研究表明，模型對流向速度u的誤差最敏感，其10%的誤差對預(yù)測結(jié)果的影響也在10%左右。因此，總體來說該模型的魯棒性較高。

僅利用CPU、單線程地預(yù)測完全湍流中240 000個樣本的流場密度脈動共用時124 s，顯著低于DNS和LES方法對計算機(jī)存儲空間和計算量的要求，滿足工程應(yīng)用對計算效率的需求。

3 結(jié) 論

針對氣動光學(xué)問題中對湍流密度脈動預(yù)測能力不足的問題，本文提出超聲速湍流密度脈動預(yù)測的神經(jīng)網(wǎng)絡(luò)方法。先通過直接數(shù)值模擬方法得到不同來流單位雷諾數(shù)下超聲速平板邊界層的平均場、湍流參數(shù)，以及密度脈動，再利用神經(jīng)網(wǎng)絡(luò)方法從流場數(shù)據(jù)中挖掘規(guī)律，建立密度脈動模型，并評估其泛化能力。數(shù)值實驗的結(jié)果表明：

1) 預(yù)測密度脈動的神經(jīng)網(wǎng)絡(luò)模型包含5個隱含層，從低到高分別有13、12、10、10和8個神經(jīng)元，采用的初始批尺寸為50，2范數(shù)正則化參數(shù)λ=2.7×10-8，動量因子α=0.2，初始學(xué)習(xí)率為η0=0.2， 輸入?yún)?shù)為?ρ/?X、?ρ/?Y、μ、u、v、Y和S這7個特征量，訓(xùn)練模型并取得了較好的預(yù)測精度。

2) 相對于傳統(tǒng)模型，該方法在極大地降低對湍流理論等先驗信息依賴的同時，也顯著提高了預(yù)測能力。

3) 對不同流場中密度脈動的評估結(jié)果表明，本文所發(fā)展的方法具有一定的泛化能力。

4) 利用模型先驗信息對輸入?yún)?shù)進(jìn)行特征選擇，模型的精確度和穩(wěn)定性得到了進(jìn)一步提升。

本文探索了神經(jīng)網(wǎng)絡(luò)方法在超聲速湍流問題中的應(yīng)用，研究是初步的，模型的泛化能力和對其他外形的適應(yīng)能力還有待于進(jìn)一步評估。