針對高階DG數(shù)值格式的非定常流場預(yù)測建模

丁子元，安 慰，劉學(xué)軍，呂宏強

(1.南京航空航天大學(xué)計算機科學(xué)與技術(shù)學(xué)院，模式分析與機器智能工業(yè)和信息化部重點實驗室，南京 211106；2.空氣動力學(xué)國家重點實驗室，綿陽 621000；3.氣動噪聲控制重點實驗室，綿陽 621000；4.軟件新技術(shù)與產(chǎn)業(yè)化協(xié)同創(chuàng)新中心，南京 210023；5.南京航空航天大學(xué)航空學(xué)院，南京 210016)

0 引言

近三十年來，隨著計算能力和數(shù)學(xué)算法的進步，CFD計算技術(shù)取得了長足發(fā)展，基于數(shù)值模擬技術(shù)的設(shè)計和優(yōu)化也日益普及，然而非定常流場的高精度求解計算過程極其耗時，是一個亟待解決的問題。高階間斷伽遼金方法（high-order discontinuous Galerkin,DG）可以在任意幾何形狀的單元獲得高階精度，被認為是解決這一問題最有前途的技術(shù)之一[1]。但是該方法增加了單元內(nèi)部插值函數(shù)的階數(shù)，也提高了單個單元的實際計算量和存儲量。總體而言，雖然得到同等精度數(shù)值結(jié)果的情況下，高階間斷伽遼金所需的未知數(shù)自由度較有限體積法要少得多，但長時間的非定常計算量仍然很大。

許多研究者希望通過降階模型（reduced order models,ROMs）來提高計算效率，比如動態(tài)模態(tài)分解（dynamic mode decomposition,DMD）方法[2]、庫普曼算子（Koopman operator）方法[3]和適用于復(fù)雜的高維度的系統(tǒng)降維提取流場主要特征的本征正交分解（proper orthogonal decomposition,POD）方法[4]。但是邊界條件、雷諾數(shù)等真實流場狀態(tài)對于這類方法的模型精度影響很大，降階模型欠缺魯棒性，并不適合應(yīng)用于所有條件下的流場[5]。

機器學(xué)習(xí)方法可以捕獲流體力學(xué)中主要的動態(tài)特性，近年來相關(guān)研究取得的進展也證明了其能夠作為流場預(yù)測的一個有力工具[6-9]。深度學(xué)習(xí)技術(shù)在近些年被廣泛關(guān)注和應(yīng)用，其不僅能夠挖掘大量數(shù)據(jù)的內(nèi)在規(guī)律，而且能夠構(gòu)建不同變量之間的非線性映射關(guān)系，是一種高效的函數(shù)逼近技術(shù)。近年來，深度學(xué)習(xí)在自然語言處理[10]、醫(yī)療診斷[11]、圖像檢測與識別[12]等相關(guān)領(lǐng)域都取得了突破性進展。更有研究者將深度神經(jīng)網(wǎng)絡(luò)應(yīng)用在基于雷諾平均模擬（Reynolds average Navier-Stokes,RANS）的湍流建模[9,13-14]以及基于大渦模擬（large-eddy simulation,LES）的亞格子應(yīng)力項建模[15-16]。這些嘗試是對湍流模擬中的某些應(yīng)力項進行參數(shù)化建模，能夠通過深度神經(jīng)網(wǎng)絡(luò)提高計算精度，在一定程度上緩解耗時計算的壓力，但往往會忽略非定常流場的時空特征對建模的重要性。

通過提取輸入流場的時間以及空間特征來預(yù)測未來時刻的流場可以看作是視頻幀的預(yù)測。視頻預(yù)測通過對視頻場景的建模幫助機器完成決策，雖然當前的研究還處于早期階段，但是文獻[17-20]利用生成對抗網(wǎng)絡(luò)、卷積神經(jīng)網(wǎng)絡(luò)等方法對這一問題進行了探索，證明了深度學(xué)習(xí)技術(shù)具備學(xué)習(xí)并預(yù)測時空序列信息的潛力。

有研究者將上面提到的相關(guān)技術(shù)成功應(yīng)用于流場預(yù)測。Lee等[21]將基于質(zhì)量和動量守恒的損失函數(shù)加入到神經(jīng)網(wǎng)絡(luò)中，運用生成對抗網(wǎng)絡(luò)（generative adversarial networks,GAN）等四種不同的方法預(yù)測了不同雷諾數(shù)下的非定常流場。Han等[22]構(gòu)建了一個由卷積神經(jīng)網(wǎng)絡(luò)（convolutional neural networks,CNN）、卷 積 長 短 期 記 憶 網(wǎng) 絡(luò)（convolutional long short term memory neural network,ConvLSTM）和反卷積神經(jīng)網(wǎng)絡(luò)（deconvolutional neural network,DeCNN）組成的混合神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)，降低了非定常流場維度并且捕獲到動態(tài)時空特征，成功預(yù)測了不同工況下未來時刻的流場。文獻[23]提出一種改進的生成對抗網(wǎng)絡(luò)，可以生成給定時間狀態(tài)下的非定常周期性流場預(yù)測結(jié)果。然而上述工作存在局限性，即訓(xùn)練樣本是由計算網(wǎng)格上有限個采樣點上的物理量構(gòu)成，若要得到流場中所有位置的物理量，一般采取的是插值方法，這會造成潛在誤差。

以上工作均將流場作為圖像來處理，區(qū)別于之前的工作，本文基于高階間斷伽遼金方法的隱式大渦模擬（implict large eddy simulation,ILES）求解過程，提出一個全新的建模思路，將求解中產(chǎn)生的高階多項式系數(shù)（即高階間斷伽遼金方法的原始變量）依據(jù)網(wǎng)格單元的空間相鄰關(guān)系排列，并按照時序關(guān)系組成序列數(shù)據(jù)。具體來說，基于高階間斷伽遼金方法的原始數(shù)據(jù)結(jié)構(gòu)構(gòu)建了一種可行的新型深度神經(jīng)網(wǎng)絡(luò)，通過直接預(yù)測非定常的高階多項式系數(shù)（即原始未知數(shù)）得到全息的非定常流場。本文首先給出了混合深度神經(jīng)網(wǎng)絡(luò)的結(jié)構(gòu)設(shè)計，然后介紹了針對高階間斷伽遼金數(shù)據(jù)結(jié)構(gòu)的樣本集構(gòu)造、訓(xùn)練和測試結(jié)果。最后給出了結(jié)論和展望。

1 混合深度神經(jīng)網(wǎng)絡(luò)

1.1 整體結(jié)構(gòu)設(shè)計

為了利用深度神經(jīng)網(wǎng)絡(luò)捕捉到高雷諾數(shù)下非定常流場內(nèi)隱含的時空特征，并以此來預(yù)測未來時刻的流場狀態(tài)，我們設(shè)計了一個由三維卷積、殘差網(wǎng)絡(luò)和通道級別注意力機制三大模塊組成的新型深度神經(jīng)網(wǎng)絡(luò)，如圖1所示。

其中，三維卷積模塊由三個三維卷積層組成，通過三維濾波器與時間上連續(xù)的多個幀組成的立方體進行卷積操作，從而獲取到相鄰幀之間的時間特征。在進行二維殘差模塊之前，需要先將上一個模塊的輸出矩陣重塑為殘差模塊的輸入矩陣所對應(yīng)的維度，然后通過由六層殘差網(wǎng)絡(luò)堆疊而成的殘差網(wǎng)絡(luò)模塊進一步提取特征表示之間的空間相關(guān)性。需要注意的是，在提取空間特征時沒有涉及任何下采樣或者池化操作，因為這些操作會導(dǎo)致張量的尺寸變小，從而使網(wǎng)絡(luò)對空間特征不敏感[24]。最后，為了識別和量化特征在局部空間中作用的不同程度，我們引入一個通道級別的注意力機制模塊，自適應(yīng)地重新調(diào)整每個通道的重要程度，進一步提升性能，且保證最終的預(yù)測值和網(wǎng)絡(luò)的輸入有相同的維度。下面將對構(gòu)成深度神經(jīng)網(wǎng)絡(luò)整體結(jié)構(gòu)的三大模塊及其作用進行詳細闡述。

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

1.2 三維卷積模塊

三維卷積最早是被提出并應(yīng)用于識別人體運動行為[18]，因為它可以從多個連續(xù)的相鄰幀之間捕獲隱含的動態(tài)信息，所以在視頻理解領(lǐng)域，也能夠幫助機器更好地學(xué)習(xí)到數(shù)據(jù)中包含的時序信息，發(fā)現(xiàn)視頻幀前后之間的相關(guān)性。

三維卷積作用于由多個相鄰幀堆疊在一起形成的立方體，然后在這個立方體中利用三維卷積核進行計算。如圖2所示，左側(cè)一列在時間上是連續(xù)的，右側(cè)的每個特征圖與上一層中兩個相鄰幀相連，且特征圖中每一位置處的值，都是通過與上一層三個連續(xù)幀的同一處的感受野經(jīng)過卷積操作之后所得，以此捕捉動態(tài)信息。在本文的深度神經(jīng)網(wǎng)絡(luò)中，位于第i層第j 個特征圖上的點（x,y,z）處的值，其數(shù)學(xué)表達式為：

圖2 三維卷積Fig.2 Three-dimensional convolution

上式中的W、H、D分別表示三維卷積核中的寬度、長度、高度三個維度的大小，tanh是非線性激活函數(shù)，是連接到上一層中第 k個特征圖的三維卷積核中的(w,h,d)位 置處的值，bij代表偏置參數(shù)。

假設(shè)圖1的三維卷積模塊中，總共堆疊了Lc個三維卷積層，下面列出兩層之間的數(shù)學(xué)表達式：

其中，x(l-1)為第l 層的輸入，?是對應(yīng)的三維卷積操作，Wl和Bl分別表示三維卷積核以及偏置參數(shù)。

1.3 殘差網(wǎng)絡(luò)模塊

殘差網(wǎng)絡(luò)可以很好地解決網(wǎng)絡(luò)退化問題。在三維卷積模塊充分提取時間維度上的信息之后，我們采取了文獻[25]中提出的一種不會給網(wǎng)絡(luò)增加額外計算量、更加容易訓(xùn)練的殘差網(wǎng)絡(luò)結(jié)構(gòu)。

如圖3所示，每一個殘差模塊由兩個激活函數(shù)和兩個二維卷積操作構(gòu)成，將殘差單元中的映射表示成f ，并且在三維卷積模塊之后堆疊Lr個殘差單元，形式上為：

其中，θl表示第l個殘差單元中所有可供學(xué)習(xí)的參數(shù)的集合， x(Lc+l-1)是 第l個殘差單元的輸入。

圖3 二維殘差模塊結(jié)構(gòu)Fig.3 Structureof the two-dimensional residual module

1.4 注意力機制模塊

注意力機制的靈感來源于對人類視覺的研究，也就是將有限的視覺信息處理資源優(yōu)先分配給關(guān)鍵的部分，并忽略掉不相關(guān)的信息。

在本文中，為了進一步提升模型的能力，一個叫做SENET（squeeze-and-excitation networks）的軟注意力機制模塊[26]（如圖4）被我們嵌入到整體網(wǎng)絡(luò)結(jié)構(gòu)的最后一部分，用于探索并自動量化通道級別上的特征對于每個區(qū)域的貢獻程度。

圖4 注意力機制模塊結(jié)構(gòu)Fig.4 Structure of the attention mechanism module

2 實驗設(shè)計

2.1 隱式大渦模擬的高階間斷伽遼金求解

非定常可壓縮N-S方程形式如下：

二維情況下守恒變量U為：

其中，ρ表示密度，u和 v表示速度的正交分量，E代表單位總能量。

本文采用高階精度間斷伽遼金法[27]對兩種不同工況下的二維圓柱繞流進行ILES數(shù)值模擬。算例1工 況 為Re=3 900,Ma=0.2，α=0°；算 例2為Re=2×104,Ma=0.2，α=0°。計算域的劃分采用了結(jié)構(gòu)化網(wǎng)格，計算網(wǎng)格如圖5所示。在網(wǎng)格單元內(nèi)部用高階多項式對守恒變量進行表達：

其中， uj(t)為 高階系數(shù)形式的原始變量，?j(x)為對應(yīng)的基函數(shù)，N是對應(yīng)階數(shù)的基函數(shù)個數(shù)。

首先采用間斷伽遼金法對N-S方程進行空間離散，本文采用BR2格式[28]，得到如下離散形式：

其中，Θ=?U為 引入的輔助變量，F(xiàn)包括對流項和擴散項。Ωe表示單元， ?Ωe表示單元Ωe的邊界，n表示單元邊界 ?Ωe的外法矢。

圖5 計算網(wǎng)格示意圖Fig.5 Schematic of the computational grid

整個離散系統(tǒng)可簡寫為如下形式：

其 中，u=[u1,u2,···,uk,···,uNele]T是 全 局 自 由 度矢量，uk代表單元k 的自由度，具體表示為uk=其中下標Nele代表全局單元的個數(shù)，下標Ne、Nd分別表示單元內(nèi)方程個數(shù)和每個變量的自由度。R(u)=[R1,R2,...,Rk,···,RNele]T為全局殘值矢量，Rk代表單元k的殘值矢量，M為全局的質(zhì)量矩陣。

然后進行隱式時間離散，殘值項取tn+1時間步變量進行計算，并對式（9）采用向后歐拉差分法，得到如下形式：

定義tn+1步 的非定常殘值Re為：

在每個時間步tn+1采用牛頓迭代方法進行求解，直到殘值 Re(uk)→0。具體步驟如下：

得到上式中的原始變量uj(t)之后，代入公式（6）得到流場中每個單元內(nèi)的物理量，即可獲得全息流場信息。

2.2 樣本集的構(gòu)造

我們的工作目標是利用計算機視覺領(lǐng)域的先進方法結(jié)合間斷伽遼金法數(shù)據(jù)結(jié)構(gòu)的特點，通過先前時間段內(nèi)的流場信息來預(yù)測未來時刻的流場。此處需要特別指出的是本文預(yù)測的對象為高階間斷伽遼金數(shù)值格式的原始變量，而非前人工作中廣泛應(yīng)用的圖像形式。

首先以Re=3 900的圓柱繞流算例為例，詳細闡述深度學(xué)習(xí)數(shù)據(jù)集的構(gòu)造，且在所有的實驗中，流場數(shù)據(jù)的生成方式、數(shù)據(jù)集的構(gòu)造方式都相同。

由于深度神經(jīng)網(wǎng)絡(luò)在圖像領(lǐng)域的應(yīng)用已非常廣泛和成熟，因此本文在用于訓(xùn)練以及測試的樣本集的構(gòu)造方面參考了圖像數(shù)據(jù)集的形式。圖5是圓柱繞流算例的計算網(wǎng)格，網(wǎng)格總量為476，我們將每一個時刻的流場信息均勻標定在大小為22×26的矩形空間內(nèi)。在數(shù)值模擬實驗中，所有的網(wǎng)格都有各自唯一確定的編號，且單元網(wǎng)格為四邊形，與之相鄰的四個網(wǎng)格（若該網(wǎng)格處于邊界，則相鄰網(wǎng)格個數(shù)小于4）是唯一確定的，所以最終矩形空間的大小由網(wǎng)格單元的數(shù)量以及相鄰關(guān)系唯一確定。如圖6所示，每一個小方塊類似于像素點，相鄰像素點代表相鄰的網(wǎng)格單元，其中的數(shù)字代表對應(yīng)的網(wǎng)格單元序號（經(jīng)過網(wǎng)格單元的重構(gòu)，圖5中的圓柱和網(wǎng)格外圍部分退化成了圖6右上區(qū)域的0號單元區(qū)域，不代表任何一個計算網(wǎng)格）。單個時間步內(nèi)，按如下方法取牛頓法迭代計算所得系數(shù)：在每個單元中提取計算守恒變量U（由四個分量組成）所需的插值多項式系數(shù)uj(t)(j=1,···,10)，采用三階格式的情況下每個單元包含40個高階表達系數(shù)。因此采用三維形式的數(shù)據(jù)22×26×40表示每個時刻的流場信息，最終將數(shù)據(jù)按照時間順序排列獲得總體數(shù)據(jù)集。

2.3 模型訓(xùn)練

本文中的兩個圓柱繞流算例皆具有準周期性流動的特征，一個周期T內(nèi)大約有250個時間步求解，其中每個時間步內(nèi)取7個式（12）～式（15）牛頓迭代步中的收斂結(jié)果，即單個時間步內(nèi)有大小為7×22×26×40的數(shù)據(jù)作為樣本，每一個算例的總體數(shù)據(jù)集都包含四個周期，由7000個樣本組成。

數(shù)據(jù)歸一化是模型訓(xùn)練前的一項基礎(chǔ)工作，可以消除維數(shù)對預(yù)測結(jié)果的影響，加快模型的收斂速度[28]。本實驗中由于數(shù)據(jù)本身具有正負屬性，為了防止數(shù)據(jù)特性的丟失，沒有采取Min-Max歸一化或Z-score歸一化，而是將所有用作神經(jīng)網(wǎng)絡(luò)輸入的樣本都在[-1,1]范圍內(nèi)進行縮放。

本文通過最小化損失函數(shù)訓(xùn)練深度網(wǎng)絡(luò)模型，損失函數(shù)采用了L1損失和L2損失結(jié)合的形式，數(shù)學(xué)表達式為：

其中，N為訓(xùn)練數(shù)據(jù)個數(shù)，∧y為預(yù)測的值，y為通過流場求解器計算所得的對應(yīng)真實值。利用自適應(yīng)矩估計[29]和反向傳播算法[30]對網(wǎng)絡(luò)進行訓(xùn)練。前者可以根據(jù)訓(xùn)練數(shù)據(jù)迭代更新神經(jīng)網(wǎng)絡(luò)權(quán)值；后者是一種廣泛使用的參數(shù)學(xué)習(xí)方法，它通過誤差反向傳播校正每個神經(jīng)元的權(quán)值，最終通過最小化損失函數(shù)來得到最優(yōu)的參數(shù)。

深度神經(jīng)網(wǎng)絡(luò)的結(jié)構(gòu)參數(shù)見表1，其中三維卷積部分的填充方式選擇的是“SAME”。在時間維度上，三維卷積層1、2、3的卷積步長分別為1、2、1，根據(jù)式（17）特征圖大小的計算方法，得到了如表1所示信息。

圖6 網(wǎng)格相鄰關(guān)系Fig.6 Neighboring relationship of the grids

需要特別指出的是，隱式時間離散過程中式（9）中的全局自由度矢量 u在 時間步tn+1可以表示為：

因為計算網(wǎng)格總量為476，本文采取三階格式，所以上式中 Nele=476 ， N=40。從式（12）～式（15）可以發(fā)現(xiàn)時間步tn+1時 刻的uk+1只與u0和uk這兩個變量相關(guān)。為了和數(shù)值模擬過程中流場原始變量的產(chǎn)生保持一致性，混合深度神經(jīng)網(wǎng)絡(luò)將u0和uk疊加作為輸入，即經(jīng)過樣本集構(gòu)造輸入數(shù)據(jù)尺寸為22×26×40×2， 模 型 的 輸 出 為 時 間 步tn+1時 刻 的uk+1， 數(shù)據(jù)尺寸為 22×26×40。

表1 深度神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)以及參數(shù)設(shè)置Table 1 Structureand parameter setting of the deep neural network

優(yōu)化后的模型能夠代替隱式大渦模擬的高階間斷伽遼金求解過程中復(fù)雜的空間離散和隱式時間離散計算過程。本文利用訓(xùn)練好的網(wǎng)絡(luò)來實現(xiàn)動態(tài)預(yù)測，將預(yù)測結(jié)果轉(zhuǎn)化為高階系數(shù)，然后代入公式（6）可計算出流場中任意單元內(nèi)的物理量，從而得到全息流場信息。隨著時間的推進，循環(huán)地將預(yù)測的結(jié)果作為下一次的輸入，因此網(wǎng)絡(luò)可以在未來時刻流場未知的情況下，通過前面時間段內(nèi)已知的流場來完成預(yù)測。

2.4 整體工作流程

圖7展示了所發(fā)展方法的整體工作流程圖。對2.1節(jié)中所要求解的控制方程進行隱式大渦模擬。首先利用間斷伽遼金法對N-S方程進行空間離散，然后采取隱式時間離散方式完成非定常時間推進。其中單個時間步內(nèi)采用牛頓法進行迭代計算，由此產(chǎn)生的大型線性系統(tǒng)利用預(yù)處理的GMERS方法計算，從而得到訓(xùn)練所需的數(shù)據(jù)，即高階表達系數(shù)。值得注意的是，相應(yīng)單元上構(gòu)造出來的高階函數(shù)可以表征網(wǎng)格內(nèi)任一點的流場信息。按照2.2節(jié)中的描述完成樣本集的構(gòu)造，之后訓(xùn)練神經(jīng)網(wǎng)絡(luò)直至滿足精度要求為止。

下面的章節(jié)中，我們將利用兩個不同來流工況的算例對所提出方法的準確性及泛化性進行討論。

3 實驗結(jié)果及分析

為驗證本文所提出的深度神經(jīng)網(wǎng)絡(luò)的預(yù)測準確性，我們首先利用Re=3 900 、Ma=0.2、α=0°工況下隱式大渦模擬的高階間斷伽遼金方法求解所得的數(shù)據(jù)對網(wǎng)絡(luò)進行離線訓(xùn)練。然后用優(yōu)化完之后的模型進行動態(tài)預(yù)測，將先前時間步下預(yù)測出的流場數(shù)據(jù)放到輸入中，實現(xiàn)遞歸預(yù)測接下來時間步的流場。最后將其預(yù)測結(jié)果與CFD計算結(jié)果做對比。圖8～圖10分別是緊接著訓(xùn)練所用數(shù)據(jù)50個時間步之后、100個時間步之后、200個時間步之后的流場情況對比。圖11和圖12分別是三個時間節(jié)點下Ma和 p的預(yù)測誤差展示圖，對比之后發(fā)現(xiàn)混合深度神經(jīng)網(wǎng)絡(luò)預(yù)測結(jié)果與CFD計算結(jié)果吻合較好。

圖8 模型預(yù)測與CFD計算在50個時間步之后的流場對比（Re = 3 900）Fig.8 Flow field comparison between the modeling prediction and the CFD simulation after 50 time steps (Re = 3 900)

選取如圖13所示的兩個點P1（1.5，0）和P2（5，0.5）來展示時間序列的預(yù)測精度。將深度神經(jīng)網(wǎng)絡(luò)在選定位置上預(yù)測所得的三個物理量u、v、p和真實樣本進行比較，得到結(jié)果如圖14。其中相對誤差曲線圖中的縱坐標表示平均百分比誤差（mean absolute percentage error,MAPE），表達式如下：

圖9 模型預(yù)測與CFD計算在100個時間步之后的流場對比（Re = 3 900）Fig.9 Flow field comparison between the model prediction and the CFD simulation after 100 time steps(Re = 3 900)

圖10 模型預(yù)測與CFD計算在200個時間步之后的流場對比（Re = 3 900）Fig.10 Flow field comparison between the model prediction and the CFD simulation after 200 time steps(Re = 3 900)

圖11 Ma 的絕對預(yù)測誤差（Re = 3 900）Fig.11 Absolute prediction error of Ma between the predicted and CFD results(Re = 3 900)

圖12 p 的絕對預(yù)測誤差（Re = 3 900）Fig.12 Absolute prediction error of p between the predicted and CFD results(Re = 3 900)

圖13 計算域中兩個指定觀測點所處位置Fig.13 Locationsof two specified observation pointsin the computational domain

圖14 預(yù)測過程中兩個指定位置處流場變量的演化（Re = 3 900）Fig.14 Evolution of flow field variables at the two specified observation pointsduring the prediction process(Re = 3 900)

為了進一步探索本文所提出模型的泛化能力，我們將神經(jīng)網(wǎng)絡(luò)用于預(yù)測流動特征更復(fù)雜的非定常流場。雖然同樣是圓柱繞流，但是該算例具有更高的雷諾數(shù)： Re=2×104，Ma=0.2，α=0°?；旌仙疃壬窠?jīng)網(wǎng)絡(luò)預(yù)測結(jié)果與CFD計算結(jié)果的對比如圖15～圖17所示。圖18和圖19分別表示三個時間節(jié)點下Ma和p的絕對預(yù)測誤差。可以發(fā)現(xiàn)，Ma、p這兩個流場變量的預(yù)測與CFD模擬流場高度吻合，特別是圖15所示較少時間步之后的預(yù)測，能夠清晰地觀測到大尺度渦從物面開始脫落以及尾部區(qū)域小尺度結(jié)構(gòu)的耗散。但是隨著遞歸預(yù)測過程的推進，來自前面時間步的細微誤差將被累積，使得預(yù)測結(jié)果與真實CFD計算結(jié)果相比出現(xiàn)不可避免的差異。例如圖16（d）中邊界層信息相比圖16（c）的真實情況產(chǎn)生了局部細微變化，圖17（b）中尾部的渦與圖17（a）中同一位置存在很小色差。

圖15 模型預(yù)測與CFD計算在50個時間步之后的流場對比（Re = 2×104）Fig.15 Flow field comparison between the model prediction and the CFD simulation after 50 time steps (Re = 2×104)

圖16 模型預(yù)測與CFD計算在100個時間步之后的流場對比（Re = 2×104）Fig.16 Flow field comparison between the model prediction and the CFD simulation after 100 time steps(Re = 2×104)

圖17 模型預(yù)測與CFD計算在200個時間步之后的流場對比（Re = 2×104）Fig.17 Flow field comparison between the model prediction and the CFD simulation after 200 timesteps(Re = 2×104)

將深度神經(jīng)網(wǎng)絡(luò)在選定位置上預(yù)測所得的三個物理量u、v、p和真實樣本進行比較，具體結(jié)果展示如圖20，所有預(yù)測結(jié)果與CFD計算結(jié)果吻合較好，表明神經(jīng)網(wǎng)絡(luò)可以在整體上預(yù)測流場特征，而且可以在時空維度上精確地預(yù)測流場的變化。

在計算耗時方面，CFD計算以及深度網(wǎng)絡(luò)訓(xùn)練與測試的具體信息如表2所示（CFD計算的硬件是Intel Xeon E5649 2.53 GHz 4核并行，采用的是并行計算方式；神經(jīng)網(wǎng)絡(luò)預(yù)測的硬件是Intel Xeon Silver 4210R 2.40 GHz 40核并行，未采用并行計算方式），可以發(fā)現(xiàn)在第二個算例中CFD數(shù)值模擬一個周期需要7500 s，基于深度學(xué)習(xí)策略的非定常流場預(yù)測方法訓(xùn)練完成之后預(yù)測一個周期所需時間為60 s。所有算例中模型預(yù)測所需時間降低了兩個量級以上，可以顯著節(jié)約計算成本。

圖18 Ma 的絕對預(yù)測誤差（Re = 2×104）Fig.18 Absolute prediction error of Ma between the predicted and CFD results(Re = 2×104)

圖19 p 的絕對預(yù)測誤差（Re = 2×104）Fig.19 Absolute prediction error of p between the predicted and CFD results(Re = 2×104)

圖20 預(yù)測過程中兩個指定位置處流場變量的演化（Re = 2×104）Fig.20 Evolution of flow field variablesat the two specified observation pointsduring the prediction (Re = 2×104)

表2 CFD計算和深度神經(jīng)網(wǎng)絡(luò)訓(xùn)練與測試時長Table 2 Training and testing time of the CFD calculation and the deep neural network modeling

4 結(jié)論

目前采用機器學(xué)習(xí)法預(yù)測非定常流場的工作廣泛采用圖像形式表征流場信息，該路線的優(yōu)勢在于可以充分利用圖像處理領(lǐng)域被廣泛應(yīng)用的機器學(xué)習(xí)方法，缺點在于預(yù)測的信息只能用圖像表達。區(qū)別于已有工作，本文采用機器學(xué)習(xí)方法直接預(yù)測高階間斷伽遼金形式的流場原始變量，保持機器學(xué)習(xí)模型預(yù)測的流場信息形式與CFD方法完全一致。

本文首先針對高階間斷伽遼金的數(shù)據(jù)特征設(shè)計了適用于深度神經(jīng)網(wǎng)絡(luò)輸入輸出的數(shù)據(jù)結(jié)構(gòu)，然后將三維卷積網(wǎng)絡(luò)、殘差網(wǎng)絡(luò)、注意力機制相結(jié)合，設(shè)計了一個新型的深度神經(jīng)網(wǎng)絡(luò)。通過對雷諾數(shù)3900和雷諾數(shù)20000（該工況下的流動具有典型的較復(fù)雜多尺度流動特征）下的兩個圓柱繞流算例進行數(shù)值模擬得到用于訓(xùn)練網(wǎng)絡(luò)的數(shù)據(jù)集，將訓(xùn)練完成的網(wǎng)絡(luò)用于預(yù)測未來時刻的流場。將預(yù)測結(jié)果與直接用CFD求解器計算的結(jié)果進行對比，非定常流場預(yù)測結(jié)果與本文采用的高精度間斷伽遼金法流場求解器計算的云圖結(jié)果高度一致，典型觀測點的預(yù)測結(jié)果與真實結(jié)果的平均百分比誤差在0.2%之下，反映出本文提出的新型深度神經(jīng)網(wǎng)絡(luò)對預(yù)測復(fù)雜非線性問題的潛力。本文發(fā)展的方法優(yōu)點在于直接預(yù)測原始變量而非流場圖像，因此獲得的是全息的流場信息，但也存在目前普遍都有的長周期預(yù)測精度的問題。因此進一步提升小樣本量情況下的長周期預(yù)測能力是未來研究的一個重要方向。