任薈穎, 王婧, 王彥棡*

1.中國科學院計算機網絡信息中心, 北京 100190

2.中國科學院大學, 北京 100049

引 言

在過去的十年中，隨著計算機硬件的不斷發(fā)展，CFD(Computational Fluid Dynamics)已經成為用于流體分析的更主要的工具。目前，對復雜流體力學系統(tǒng)進行建模主要有兩種途徑。第一種是基于理論的模型架構，即根據物理問題的控制方程，建立理想的系統(tǒng)描述；第二種則是目前比較流行的數據驅動方法[1-2]，即根據系統(tǒng)仿真或試驗中的樣本數據，直接構造黑箱或灰箱模型。近年來，隨著計算機性能和精細化流動測試手段的發(fā)展，研究者能夠獲得高精度、高時空分辨率的流場信息。如何高效地利用這些大數據，從中提取出關鍵信息，并指導流體力學的發(fā)展，已經成為研究者關注的焦點。作為處理和分析數據的主要手段，數據挖掘、統(tǒng)計學習和機器學習等技術，則為開展此類研究提供了重要基礎[3]。

雷諾時均Navier-Stokes(RANS)方程[4]提供了一種計算時間平均湍流量的有效方法，由于其計算易處理性而被廣泛使用。目前，RANS 模型在實際生產和工程中仍占有重要地位。而在RANS 計算的過程中，能否準確求解雷諾應力至關重要。近些年來，人們開始關注使用機器學習和深度方法來求解RANS 方程。Tracey 等人[5]使用核回歸來模擬雷諾應力各向異性特征值。該方法被證明在大量訓練數據上的應用效果較差。Tracey，Duraisamy 和Alonso[6]后來使用具有單個隱藏層的神經網絡來模擬Spalart Allmaras 模型中的源項。實驗證明這些神經網絡能夠重建這些源項。 Ma，Lu 和Tryggvason[7]使用神經網絡代替流體閉合方程來模擬等溫氣泡流。Zhang 和Duraisamy[8]使用神經網絡來預測湍流產生項的修正因子。該修正因子會影響預測的雷諾應力張量的幅度但不影響其各向異性。 Ling J，Ruiz A 和Lacaze G[9]等人提出使用隨機森林來預測雷諾應力各向異性。該方法被證明不能輕易地對張量數量強制執(zhí)行伽利略不變性，因此在預測完全各向異性張量方面的能力有限。Ling J，Kurzawski A 和Templeton J[10]提出一種新的深度神經網絡結構從高保真仿真數據中學習雷諾應力各向異性張量，并且他們展示了基于深度學習的雷諾應力對不同幾何形狀的流動的預測能力。這種方法被證明對DNS 仿真效果不好，并且只針對于定常流動。Kutz[11]解決了基于DL 的閉合方程應用的幾個開放性挑戰(zhàn)，解決了應該選擇怎樣的訓練數據的問題。Chang 和Dinh[12]分析了五種類型的機器學習框架，討論了如何耦合機器學習框架和流體模擬的問題。Chang 和Dinh[13]后來采用了深度學習的方法來擬合雷諾應力。這種方法被證明對于非定常流動的效果擬合較差，無法擬合出精確的流場效果。

近些年來，自動機器學習技術(Auto Machine Learning，AutoML)有了長足的突破和發(fā)展。在深度學習領域，自動機器學習技術包含自動調節(jié)超參 數[14-15](例如學習率)、自動調整網絡結構[16-17](例如神經元的個數，層數等)，在深度學習領域取得了不錯的成績。由于在湍流問題中，不同的系統(tǒng)初始條件不同，數據的質量千差萬別，難以使用統(tǒng)一的神經網絡進行訓練。因此用深度學習中AutoML 的方法進行擬合，針對不同的數據集自動選取不同的網絡結構和參數。此外本文還通過混合多種初始條件下的數據進行模型的訓練來提高擬合的精度和模型的魯棒性。實驗表明，擬合效果較之前的方法有了較大程度的提升。

1 相關理論及方法

1.1 AutoML 貝葉斯優(yōu)化理論

機器學習中，調參和搜索網絡結構被認為是一項繁瑣且耗費人力的工作。手動調參十分耗時，參數的好壞影響網絡的性能。 網格和隨機搜索需要很長的運行時間，因此誕生了很多自動調參的方法。 長期以來，高斯過程一直被認為是一種針對模型的建模損失函數的好方法[18]。高斯過程(GPs[19])具體描述為建立目標函數的概率模型，用它來選擇效果最好的超參數，在真實的目標函數中進行評價。

超參數優(yōu)化的過程可以描述為：

這里(xf)表示在驗證集上評估得出的最低損失函數的值(例如，RMSE 或錯誤率)；x*是產生該最低值的超參數集，并且x可以在超參數域χ中取任意值。簡單來說，我們想找到在驗證集度量標準上產生最佳分數的超參數集合。

高斯過程通過觀察過去的評估結果，計算超參數映射到目標函數得分概率的概率模型：

此函數被稱為目標函數的“代理”。高斯過程方法通過選擇最有潛力的超參數組，通過目標函數進行評估，再帶入到高斯過程模型中，進而選擇下一組超參數，反復迭代。

具體的算法過程可以描述如下：

(1)建立高斯函數模型M，假設損失函數f的輸出值y與超參數x之間的關系服從M；

(2)從M中采集下一組函數值較高的參數x1+t；

(3)將x1+t帶入到損失函數f并計算出分數y1+t；

(4)將x1+t和y1+t更新到M中；

(5)重復步驟(1)到(4)直到滿足停止條件。

貝葉斯優(yōu)化方法在理論上是有效的，因為它以有根據的方式選擇了下一個超參數。它的基本思想是：花更多的時間選擇下一個超參數，以減少對目標函數的調用。實際上，與在目標函數中評估花費的時間相比，選擇下一個超參數所花費的時間是很少甚至無關緊要的。通過評估，從過去的結果選取出看似更有希望的超參數，貝葉斯方法可以在更少的迭代中找到比隨機搜索更好的模型參數設置。

1.2 深度學習在湍流模擬中的應用方法

機器學習與湍流建模相結合的研究工作是流體力學領域新興的熱門研究方向?，F(xiàn)有的研究成果有力地驗證了其有效性和可行性，預示了機器學習在未來湍流模型應用中良好的應用前景[20]。一方面，這些研究工作以數值模擬器或實驗產生的高分辨率的數據作為驅動，一定程度上降低了機器學習模型封閉或湍流模型相關變量模型化的難度，證實了純數據驅動的機器學習黑箱模型在湍流研究應用中的有效性和可行性。另一方面，一些研究者還將機器學習用于描述和量化傳統(tǒng)模型計算結果的不確定度，對未來的建模工作具有很好的指導作用[21]。主要的研究方向和建模流程歸結為圖1所示。

圖1 機器學習應用于湍流研究的主要研究方向及流程Fig.1 Main research direction and process of machine learning applied to turbulence study

2 本文方法

2.1 RANS-DL 模型

我們假設深度學習可以從不同流場區(qū)域中采集的速度場中發(fā)現(xiàn)隱藏的時間導數[13]。因此，我們可以從各種模擬時間步驟中采樣數據，并使用總數據訓練神經網絡。因為深度學習屬于監(jiān)督學習方法，需要輸入和輸出來訓練模型。因此，對于深度學習湍流模型，輸入應該是能表示平均流動特性的流動特征，目標為雷諾應力六個分量。我們基于RANS方程[22]選擇輸入的流動特征。RANS 方程推導輸入輸出的過程如下所示：

圖2 RANS-DL 模型Fig.2 RANS-DL model

這里是平均速度，i，j，k是方向，ρ，，τij，μ和分別為流體密度，平均壓力，雷諾應力分量，粘性系數和Kronecker delta 函數?；谖覀兊募僭O，每個數據集的數據穩(wěn)定，因此不選擇時間導數作為訓練輸入。而壓力是從動量方程中單獨求解的，因此我們進一步發(fā)現(xiàn)壓力項可以從訓練輸入中排除。因此基于深度學習的湍流模型的輸入可以通過速度的空間導數來表示。使用矩陣形式顯示雷諾應力與導數算子、速度的之間的二元乘積，具有的函數關系如公式(6)所示：

由于DNN 能夠較好的解決回歸問題，我們使用DNN 來回歸擬合替雷諾應力。損失函數采用歐幾里德?lián)p失函數，其中，分別為訓練數據的總數，訓練結果和雷諾應力值，如公式(7)所示：

圖2 為模型的具體形式，只采用每個點的9 個速度導數作為輸入流特征，目標是它的六個雷諾應力分量。

2.2 RANS-DL 工作流程

我們假設閉包關系可以從守恒方程實驗數據和數據驅動的機器學習模型中分別導出。如圖3 所示，RANS-DL 工作流具體描述如下：

(1)首先定義了某一問題的物理幾何邊界，如模擬結構和系統(tǒng)特性；

(2)利用該幾何邊界構造了CFD 計算實例，包括初始條件、邊界條件等；

(3) CFD 求解器執(zhí)行模擬；

(4)從模擬中提取和處理相關數據；

(5)使用這些處理后的數據來訓練神經網絡。

2.3 RANS-AutoML 模型

本文采用了基于AutoML 的技術，能夠根據不同的初始條件的數據集自動設計網絡結構以及調整參數。本文采取貝葉斯優(yōu)化器負責自動構建網絡模型結構（例如RANS-DL 模型中隱藏層的節(jié)點數，隱藏層的層數），以及自動搜索模型的參數（例如RANS-DL 模型中的學習率，學習率衰減指數以及L2 正則系數）。 貝葉斯優(yōu)化器需要先用隨機搜索訓練n 個網絡，觀察n 個網絡的性能，從這些過去的經驗中更新優(yōu)化器。主要流程和思想如圖4所示。

具體的流程如下：

(1)利用隨機搜索獲得n 個網絡結構和MSE；

(2)從舊的模型結構，參數以及MSE 值生成下一輪需要搜索的模型結構和參數用于模型訓練；

(3)模型訓練，在驗證集上驗證產生MSE 值；

(4)將本輪的模型結構和參數以及得到的MSE的值再反饋到貝葉斯優(yōu)化器中；

(5) 重復步驟(1)～(4)直到達到搜索時間或者搜索網絡結構得到的MSE 值達到預期，便結束得到最優(yōu)網絡的結構及參數。

圖3 RANS-DL 工作流程Fig.3 RANS-DL workflow

圖4 AutoML 工作流程Fig.4 AutoML workflow

3 實驗設置

3.1 數值實驗

進行數值實驗來產生深度學習模擬雷諾應力需要的數據。采用k-ε模型求解RANS 模擬得到參考解，用于模型的訓練。圖5 描繪了Pitz 和Daily[23]提出的內壁形狀，這是一個2D 的幾何形狀，包括一個向后的臺階和一個噴嘴。系統(tǒng)特性總結在表1 中。由于不穩(wěn)定的流動受到湍流混合層、剪切速率等的影響，因此這種幾何形狀足夠復雜。對于這種幾何形狀，k-ε模型已被驗證是可行的[24]。我們用OpenFOAM[25]中的pimpleFoam 求解器用于生成數據。

使用k-ε模型通過RANS 仿真生成訓練數據。 采用pimpleFoam 求解方程，固定時間步長5×10-10秒。為了觀察不同的系統(tǒng)初始條件對實驗結果的影響，我們設計了8 個數據集，并在表2 中列出。設計實驗如下：

(1)采用A1、A2、A3 分別訓練，并在驗證集B上驗證觀察實驗結果，觀測對訓練數據進行時間和空間的加密對測試結果有無影響；

(2)采用A1、A4、A5/A1、A6、A7 分別/混合訓練，并在驗證集B 上驗證觀察實驗結果，觀測改變訓練數據的入口流速/層流粘性系數對測試結果有無影響；

(3)采用A1、A6、A7 分別/混合訓練，并在驗證集B 上驗證觀察實驗結果，觀測改變訓練數據的層流粘性系數對測試結果有無影響。

表1 系統(tǒng)邊界信息Table 1 System boundary information

3.2 基于AutoML 的雷諾平均湍流模型求解模型

圖5 數值模擬實驗信息Fig.5 Experimental information of numerical simulation

表2 數據集信息Table 2 Dataset information

我們使用Tensorflow[26]來設計DNN 結構，激活函數采用ReLU[27]，并且在每個隱藏層采用BN[28]。 模型采用Adam[29]算法作為深度學習的優(yōu)化算法來更新模型的參數，例如權重w和偏差b等。模型的網絡結果和具體參數，包含：網絡的層數，每層隱藏單元的個數、學習率、L2 正則系數均用前文介紹的貝葉斯優(yōu)化方法進行搜索。

4 實驗結果及分析

為了更好的展示實驗的結果，如圖5 所示，我們選取了驗證集上流場的不同區(qū)域來觀察模型訓練的效果。這些區(qū)域包含入口，臺階流，流場中較為平緩的區(qū)域以及出口處的區(qū)域，比較能代表流場的情況，從而反應出模型訓練的效果。

4.1 AutoML 實驗結果

網絡結構搜索以及參數調整的小時數為48 小時，本文做實驗的機器配置：CPU 為雙核4 線程， GPU 為Tesla V100。最終AutoML 的訓練效果如圖6 所示。經過訓練，本文最終采用的模型參數如表3所示。

圖6 AutoML 的效果Fig.6 Effect of AutoML model

表3 模型參數Table 3 model parameter

如圖7 所示，與Chang 和Dinh[10]采用的神經網絡結構來擬合雷諾應力的效果對比，我們用AutoML 獲得的網絡結構擬合效果更為優(yōu)異。

圖7 AutoML 的網絡結構與之前方法的對比Fig.7 Comparison between AutoML network structure and previous methods

4.2 時間信息對模型的影響

我們假設深度學習可以從不同流場區(qū)域中采集的速度場中發(fā)現(xiàn)隱藏的時間導數，因此需要驗證訓練好的模型能否重建驗證集上的不同時刻的雷諾應力的參考解。 如圖8 所示，在模擬過程中，每個雷諾應力張量分量的MSE 隨時間保持穩(wěn)定。因此可以得出結論，DNN 模型可以學習到流場數據中隱含的時間導數信息。

圖8 雷諾應力張量分量隨時間變化擬合效果Fig.8 MSE along with the time during simulation

4.3 網格疏密對模型的影響

為了檢驗數值實驗中網格的疏密對深度學習網絡模型的訓練效果有無影響，我們將模型分別以四倍的網格差在不同疏密程度的數據上進行訓練，同時為了降低數據規(guī)模變化對模型訓練造成的影響，將稀疏網格的數據在時間采樣上加密了四倍，在驗證集上進行驗證。如圖9 所示，從粗網格生成的訓練數據的訓練結果與從細網格生成的數據訓練結果幾乎相同。

圖9 不同疏密網格數據訓練模型的效果Fig.9 MSE with different grid quantity

4.4 擾動入口流速對模型的影響

為了檢驗入口流速的改變是否會對模型的訓練結果造成影響，將模型分別在入口流速為5m/s、20m/s 以及混合多種流速的情況下進行訓練，并在驗證集上進行驗證。如圖10 所示，模型在較小的入口流速下訓練模型，會傾向于預測較小的結果；模型在較大的入口流速下訓練模型，預測的結果也會偏大；而在混合入口流速下，訓練的效果更好，模型的魯棒性也更好。

4.5 擾動層流粘性系數對模型的影響

為了檢驗層流粘性系數的改變是否會對模型的訓練結果造成影響，將模型分別在層流粘性系數為 0.5e-5m/s、1e-5m/s、2e-5m/s 以及混合多種層流粘性系數的情況下進行訓練，并在驗證集上進行驗證。如圖11 所示，模型在較小的層流粘性系數下訓練模型，會傾向于預測較大的結果；模型在較大的層流粘性系數下訓練模型，預測的結果會偏??；而在混合層流粘性系數下，訓練的效果更好，模型的魯棒性也更好。

5 結論與展望

本文提出了一種用AutoML 方法建模RANS 模型 的方法，該方法使用AutoML 方法來再現(xiàn)RANS(k-ε) 數值模擬的解。通過實驗結果的分析, 可以看出基于 AutoML 的方法能夠自動構建網絡結構和調整參數, 在保證準確率的情況下減少人力的工作。另外，本文還發(fā)現(xiàn)了一些有意義的結論：速度場的一階空間導數的流動特征對于重建RANS 結果是必要且充分的，模型能學到流場數據中隱含的時間導數信息；另外，在混合入口流速、層流粘性系數的數據上訓練，能使得模型的魯棒性更好. 目前，系統(tǒng)的針對湍流機器學習的特征選取的研究還相對較少。在現(xiàn)有工作積累的經驗基礎上，機器學習在未來的湍流模型化中必將扮演著更加重要的角色。

圖10 不同入口流速下數據訓練模型的效果Fig.10 The effect of data training model under different inlet velocity

圖11 不同層流粘性系數下數據訓練模型的效果Fig.11 The effect of data training model under different laminar viscosity coefficient

利益沖突聲明

所有作者聲明不存在利益沖突關系。