趙 暾,周 宇,程艷青,錢煒祺,*
(1. 空氣動力學(xué)國家重點實驗室,綿陽 621000;2. 中國空氣動力研究與發(fā)展中心,綿陽 621000)
自從Hinton[1]在2006年提出深度學(xué)習(xí)的概念之后,最近幾年深度學(xué)習(xí)在科學(xué)領(lǐng)域和商業(yè)領(lǐng)域已經(jīng)取得了巨大的成功。一般來說,深度學(xué)習(xí)需要大量的訓(xùn)練樣本,通常應(yīng)用在大數(shù)據(jù)領(lǐng)域,而且由于過度依賴訓(xùn)練數(shù)據(jù),導(dǎo)致泛化能力差,同時還具有以下缺點:缺乏推理和因果關(guān)系的表達能力而導(dǎo)致表征能力弱;無法解釋形成決策的推理過程;很難在計算過程中導(dǎo)入選擇性注意力機制。
在實際工程領(lǐng)域,受成本限制,所獲取到的樣本數(shù)據(jù)量是有限的。因此,很難將深度學(xué)習(xí)技術(shù)直接應(yīng)用于小樣本問題領(lǐng)域。但是這些數(shù)據(jù)通常都包含由物理定律、經(jīng)驗規(guī)則或者領(lǐng)域知識構(gòu)成的先驗信息。如何將深度學(xué)習(xí)方法有效運用到包含先驗信息的小樣本數(shù)據(jù)領(lǐng)域,相關(guān)學(xué)者做了大量研究。其中,美國布朗大學(xué)Karniadakis教授團隊做出了較為出色的工作。他們通過將蘊藏在小樣本數(shù)據(jù)中的先驗信息編碼到深度網(wǎng)絡(luò)中,提出了一種新的深度學(xué)習(xí)算法—內(nèi)嵌物理機理神經(jīng)網(wǎng)絡(luò)[2](Physics-informed neural networks,PINNs),該網(wǎng)絡(luò)可以由神經(jīng)網(wǎng)絡(luò)對輸入及模型參數(shù)進行自動微分所獲取。PINNs通常受到數(shù)據(jù)中包含的物理定律的約束,而這些物理定律一般可以建模成時間相關(guān)的非線性偏微分方程。這樣就可以利用PINNs發(fā)展出一類新的偏微分方程數(shù)值求解方法以及一類新的系統(tǒng)辨識和模型反演的數(shù)據(jù)驅(qū)動方法。PINNs方法的核心思想是將上述原問題(即求解偏微分方程的正問題或者逆問題)轉(zhuǎn)換為優(yōu)化問題求解。與傳統(tǒng)的偏微分方程數(shù)值求解方法以及系統(tǒng)辨識方法相比,PINNs具有以下幾個特點:
1)PINNs可以作為訓(xùn)練辨識一體化的方法,即在求解偏微分方程的同時,根據(jù)附加信息直接辨識出方程中的未知參數(shù)。
2)與傳統(tǒng)的基于網(wǎng)格的偏微分方程數(shù)值求解方法相比,PINNs是無網(wǎng)格的,這樣可以節(jié)省大量的網(wǎng)格繪制與生成時間[3],而且還可以避免求解高維PDEs時會遇到的維度災(zāi)難問題。
3)與傳統(tǒng)的基于變分原理求解偏微分方程方法相比,PINNs可以求解不滿足變分原理的偏微分方程[4]。
4)由于PINNs采用了自動微分,因此可以避免傳統(tǒng)數(shù)值方法中由差分格式帶來的誤差影響。而且自動微分是直接對神經(jīng)網(wǎng)絡(luò)的輸出而不是對數(shù)據(jù)進行微分,因此,PINNs對包含噪聲的輸入數(shù)據(jù)具有較強的魯棒性。
5)PINNs可以有效地解決不適定的偏微分方程逆問題。
6)由于PINNs中的損失函數(shù)是非凸的,會存在多個局部極小值,因此基于梯度下降的優(yōu)化算法幾乎一定會陷入某個極小值,從而導(dǎo)致PINNs的求解結(jié)果不唯一。
7)PINNs對于每個新的實例(例如,方程參數(shù)變化)就需要重新訓(xùn)練網(wǎng)絡(luò)。后續(xù),Lu等[5]提出了DeepONet去解決這個問題。該網(wǎng)絡(luò)可以直接對數(shù)據(jù)中的非線性算子進行近似,例如,如果數(shù)據(jù)是從偏微分方程獲得,則可以直接學(xué)習(xí)出積分算子,從而近似出偏微分方程的解。
Jagtap等[6]將PINNs中的激活函數(shù)修改為自適應(yīng)激活函數(shù),大大加速了訓(xùn)練收斂速度。具體做法是在激活函數(shù)中引入一個超參數(shù),并將該參數(shù)作為神經(jīng)網(wǎng)絡(luò)參數(shù)進行優(yōu)化。在PINNs的理論基礎(chǔ)層面,Shin等[4]針對橢圓型和拋物線型偏微分方程,證明了PINNs訓(xùn)練得到的網(wǎng)絡(luò)可以一致收斂到偏微分方程的解。
經(jīng)過最近兩年的發(fā)展,PINNs發(fā)展出了一大類家族—cPINNs[7]、fPINNs[8]、nPINNs[9]、sPINNs[10-11]、vPINNs[12]以及XPINNs[13]等。前面四種可以分別用來求解離散域非線性守恒律方程、分?jǐn)?shù)階偏微分方程、非局部橢圓型偏微分方程、隨機偏微分方程,vPINNs則是以Petrov-Galerkin方法為基礎(chǔ)求解偏微分方程,而XPINNs可以求解任意復(fù)雜幾何域上的非線性偏微分方程。
Lu等[14]基于Tensorflow開發(fā)出了一個Python工具包—DeepXDE,該工具包對上述算法進行了整合,并且還可以用來進行多保真度數(shù)據(jù)建模,以及從數(shù)據(jù)中直接學(xué)習(xí)非線性算子和函數(shù)關(guān)系。更進一步地,PINNs由于其原理簡單,靈活度高,效果良好,已被英偉達、ANSYS、西門子等公司采用。英偉達已經(jīng)開發(fā)出了SimNet[15]工具包,用來求解偏微分方程的正問題、逆問題以及數(shù)據(jù)融合問題。
雖然PINNs方法目前發(fā)展迅速,但是對求解正問題時的精度以及求解逆問題時的魯棒性還需要進一步分析。本文則針對一維熱傳導(dǎo)方程對上述問題進行了研究。
下面以求解非線性偏微分方程的正問題和逆問題為例,闡釋PINNs的基本網(wǎng)絡(luò)結(jié)構(gòu)及其訓(xùn)練方式。
非線性偏微分方程的一般形式可以表示為:
其中,u是上述偏微分方程的解, N[u,λ]為帶有參數(shù)λ的非線性算子, Ω 是RD的 子集,tT為終端時刻,h(x)為方程的初始條件,g(x,t)為方程的邊界條件。
首先,建立一個深度網(wǎng)絡(luò)n(x,t;θ)去近似u,并定義偏微分方程的殘差為PINN的輸出:
上式中的偏導(dǎo)數(shù)都可以通過自動微分得到。相應(yīng)地,損失函數(shù)可以給定為:
損失函數(shù)的訓(xùn)練集分為三部分,第一部分為初始條件和邊界條件的采樣點集,第二部分為對偏微分方程殘差進行采樣的配置點集。配置點集一般有兩種選取方式,第一種為類網(wǎng)格點選取方式,第二種為偽隨機離散取點方式,如圖1所示。當(dāng)針對高維非常規(guī)區(qū)域上偏微分求解問題,偽隨機離散取點方式具有更強的靈活性。第三部分為已知的數(shù)據(jù)樣本集,該部分標(biāo)簽數(shù)據(jù)用來對方程中的參數(shù)進行辨識。
圖1 配置點離散方式Fig. 1 The distribution of collocation points
從損失函數(shù)式(3)的組成可以看出,偏微分殘差項起到了正則化作用,這樣可以增強網(wǎng)絡(luò)n(x,t;θ)的泛化能力。而從另一個角度看,通過自動微分構(gòu)建偏微分方程殘差并將其引入到損失函數(shù)的過程,可以看作是將該偏微分方程表示的物理規(guī)律嵌入到深度網(wǎng)絡(luò)的過程。
最后,基于該損失函數(shù),通過梯度下降方法對深度網(wǎng)絡(luò)參數(shù)θ和方程參數(shù)λ進行訓(xùn)練。
最終,通過訓(xùn)練就可以求得使損失函數(shù)(3)最小的 參 數(shù) θ?和 λ?,而 方 程(1)的解可以用神經(jīng)網(wǎng)絡(luò)=n(x,t;θ?)來表示。PINNs的網(wǎng)絡(luò)結(jié)構(gòu)及訓(xùn)練過程可以用圖2表示。
圖2 PINNs基本原理圖Fig. 2 A sketch of PINNs
在物體邊界給定為第二類邊界條件后,可得到如下的一維熱傳導(dǎo)控制方程[16]:
邊界條件:
初始條件:
本節(jié)的目的是在已知邊界條件和初始條件的情況下,利用PINN求解導(dǎo)熱系數(shù)k= 1時方程(5)的解。
首先設(shè)計一個深度網(wǎng)絡(luò)n(x,t;θ)去近似方程(5)的解,即:
該網(wǎng)絡(luò)包含2個輸入,1個輸出,3個隱層且每個隱層包含40個節(jié)點,隱層節(jié)點的激活函數(shù)設(shè)置為sigmoid。隨后對n(x,t;θ)的輸出進行自動微分,定義PINN為:
接著,選擇如下均方誤差的和作為損失函數(shù):
其中,Nb和N0分別為邊界點和初始點訓(xùn)練集大小,Nf為方程配置點訓(xùn)練集大小。在本文中,首先給定Nb和N0分 別為80和100,以及在時間空間域 Ω×R1內(nèi)采用類網(wǎng)格點方式選擇Nf= 8000個配置點。接著,采用Adam優(yōu)化器并設(shè)置學(xué)習(xí)率為0.001對所建的網(wǎng)絡(luò)n(x,t;θ)進 行訓(xùn)練。即在訓(xùn)練過程中,通過給定邊界條件數(shù)據(jù)集、初始條件數(shù)據(jù)集{ti,xi;Ti}以及方程殘差點集 {ti,xi;fi},通過優(yōu)化尋找最接近觀測數(shù)據(jù)集以及偏微分方程殘差的參數(shù) θ?所確定的深度網(wǎng)絡(luò)。最后,利用訓(xùn)練完的網(wǎng)絡(luò)預(yù)測區(qū)域(x,t)∈[0,1]×(0,1]內(nèi)的溫度值。
還需要特別注意的是,由于PINN會從不同的網(wǎng)絡(luò)初值收斂到不同的解,因此,在本文中,隨機初始化20次訓(xùn)練PINN,損失函數(shù)隨迭代步數(shù)變化情況如圖3所示。
圖3 20次不同初始權(quán)重網(wǎng)絡(luò)的損失函數(shù)隨迭代步數(shù)變化Fig. 3 Loss functions with 20 different initial weights versus the iteration step
從圖3可以看出,從不同的初始權(quán)重出發(fā),網(wǎng)絡(luò)經(jīng)過大約10000次左右的迭代訓(xùn)練后,損失函數(shù)都可以收斂。
然后利用訓(xùn)練完的網(wǎng)絡(luò)分別預(yù)測溫度場,最后將平均值作為最終結(jié)果,并將其與溫度場的解析解進行對比,如圖4所示。與解析解相比,利用PINN預(yù)測整體溫度場的均方誤差約為2.2865× 10?7。
圖4 基于PINN的一維熱傳導(dǎo)方程求解結(jié)果與解析解對比Fig. 4 Comparison of the solution obtained by PINN and the analytical solution of the one-dimensional heat conduction equation
通過對PINN方法分析可知,其求解誤差主要包括有限訓(xùn)練點集導(dǎo)致的采樣誤差、神經(jīng)網(wǎng)絡(luò)的逼近誤差以及陷入局部最優(yōu)解導(dǎo)致的優(yōu)化誤差等三部分組成。下面在神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)固定的情況下,即逼近誤差一定,分別從空間和時間兩個尺度分析采樣誤差以及優(yōu)化誤差對精度的影響,并且與控制體積方法的求解結(jié)果進行對比分析。為分析空間尺度影響,采用以下配置點選取方式:固定100個時間點,分別選取20、40、80、160、320、640和1280個空間點;而為分析時間尺度影響,則采用以下配置點選取方式:固定80個空間點,分別選取25、50、100、200、400、800和1600個時間點,具體如表1和表2所示。
表1 分析空間尺度離散影響的采樣點選取方式Table 1 The selection method of sampling points to analyze the influence of spatial discretization
表2 分析時間尺度離散影響的采樣點選取方式Table 2 The selection method of sampling points to analyze the influence of temporal discretization
本文訓(xùn)練環(huán)境為英偉達RTX2060顯卡,顯存6 GB,運行頻率1800 MHz。在該環(huán)境中隨機初始化20次訓(xùn)練網(wǎng)絡(luò),訓(xùn)練時間隨采樣點數(shù)的變化如圖5所示??梢钥闯?,隨著采樣點數(shù)的增加,訓(xùn)練時間呈線性增長的趨勢。
圖5 訓(xùn)練時間隨采樣點數(shù)的變化Fig. 5 The training time as a function of the sampling number
利 用 訓(xùn) 練 完 的 網(wǎng) 絡(luò) 計 算 其 在x= 0 m、0.1 m、0.25 m、0.6 m及1 m處的平均預(yù)測值及其對應(yīng)的均方誤差,并與控制體積的求解結(jié)果相對比,如圖6所示,圖中的黑色虛線表示的是O(N?2)收斂速率,其中N為采樣點數(shù)。
圖6 PINN方法與控制體積方法求解精度比較Fig. 6 CA comparison of the solution accuracy of between the PINN and the control volume method
通過對圖6的分析,可以得出以下幾個結(jié)論:
1)當(dāng)采樣點數(shù)較少時,時間和空間尺度采樣點數(shù)對PINN方法的求解精度影響與控制體積方法類似,即在空間采樣點固定情況下,求解誤差按照O(N?2)速率收斂,而且時間尺度的采樣誤差相較空間尺度占主要地位。
這是由于在損失函數(shù)(10)中,時間尺度采樣點對應(yīng)的有左邊界殘差和右邊界殘差兩項,而空間尺度采樣點只對應(yīng)一項初始條件殘差。
2)PINN方法在采樣點數(shù)較少時,采樣誤差為主要誤差來源,隨著采樣點數(shù)增加,其求解精度隨之提高,但是當(dāng)總采樣點數(shù)約大于8000時,兩種配置點采樣方式的求解精度幾乎收斂到同一水平。
這是由于過多的采樣點數(shù)會導(dǎo)致過于復(fù)雜的損失函數(shù),從而使優(yōu)化誤差增大,這讓優(yōu)化誤差成為主要誤差來源。
下面,固定8000個配置點,分別分析初始和邊界條件的采樣點對求解精度的影響。左右邊界分別選取100個采樣點,初始條件分別選取20、40、80、160、320、640和1280個采樣點;初始條件選取80采樣點,左右邊界分別選取25、50、100、200、400、800、1600個采樣點時,方程(5)求解精度如圖7所示。
從圖7可以得出與圖6一致的結(jié)論,即采樣點數(shù)少時,采樣誤差主導(dǎo)了收斂誤差,當(dāng)采樣點數(shù)增加到一定程度后,優(yōu)化誤差成為了收斂誤差的主要組成部分。
圖7 配置點固定時PINN方法求解精度分析Fig. 7 The solution accuracy of the PINN with fixed collocation points
本節(jié)針對熱傳導(dǎo)方程(5),在已知初始條件、混合邊界條件:
以及Nt個測點溫度的情況下,對方程中的導(dǎo)熱系數(shù)k進行辨識,其中,T0為基于2.1節(jié)方程的解析解。通過對上節(jié)方法稍加改動,就可以對參數(shù)k進行辨識。具體實現(xiàn)步驟如下:
首先,在構(gòu)建PINN時,將參數(shù)k作為待優(yōu)化變量,即PINN的輸出為:
接著,將Nt個測點溫度作為標(biāo)簽數(shù)據(jù),定義如下?lián)p失函數(shù):
最后,通過Adam算法對所構(gòu)建的網(wǎng)絡(luò)進行優(yōu)化,得到使損失函數(shù)最小的k?。
本文選取網(wǎng)絡(luò)結(jié)構(gòu)如下:2個輸入及1個輸出,6個隱層且每個隱層包含20個節(jié)點,隱層節(jié)點的激活函數(shù)設(shè)置為tanh。在固定網(wǎng)絡(luò)初始參數(shù)的情況下,同時給定無噪聲、10%噪聲以及50%噪聲測點數(shù)據(jù),經(jīng)過10次蒙特卡洛仿真,PINN訓(xùn)練過程中的損失函數(shù)變化如圖8所示??梢钥闯觯S著數(shù)據(jù)中噪聲的增大,損失函數(shù)的收斂值也變大。
圖8 給定不同噪聲數(shù)據(jù)時損失函數(shù)隨迭代步數(shù)變化Fig. 8 Loss functions versus the iteration step for training data with different noisy levels
隨后,將PINN方法與共軛梯度法(CGM)的辨識結(jié)果進行對比,如表3所示。從表3可以看出,基于PINN的辨識結(jié)果與CGM結(jié)果精度相當(dāng)。對于50%的噪聲數(shù)據(jù),基于PINN的辨識結(jié)果有1.253%的相對誤差,這說明利用PINN的參數(shù)辨識方法對噪聲數(shù)據(jù)具有較強的魯棒性。這是因為,第一,損失函數(shù)中包含初始、邊界條件以及偏微分方程殘差構(gòu)成的正則化項,可以防止網(wǎng)絡(luò)對噪聲過擬合;第二,由于采用了自動微分,一方面在計算溫度的偏導(dǎo)數(shù)時,不會引入標(biāo)簽溫度的噪聲,另一方面,也可以避免傳統(tǒng)數(shù)值方法中低精度差分格式對辨識結(jié)果的影響[17]。在計算過程中注意到,網(wǎng)絡(luò)的初始權(quán)重對結(jié)果影響較大,當(dāng)選取不恰當(dāng)時,會導(dǎo)致訓(xùn)練發(fā)散。
表3 導(dǎo)熱系數(shù)k辨識結(jié)果Table 3 Identification results of the thermal conductivity k
下面研究小導(dǎo)熱系數(shù)以及邊界條件對基于PINN的辨識方法的影響。
已知導(dǎo)熱系數(shù)k= 0.01時x= 0.1 m處的溫度場測點數(shù)據(jù),以及邊界條件(11),利用PINN方法的導(dǎo)熱系數(shù)辨識結(jié)果如表4所示。從表4可以看出,在給定較小導(dǎo)熱系數(shù)的測點數(shù)據(jù)時,基于PINN方法可以對參數(shù)k進行高精度的辨識。
表4 導(dǎo)熱系數(shù)k辨識結(jié)果Table 4 Identification results of the thermal conductivity k
當(dāng)給定第二類邊界條件(6)時,將導(dǎo)熱系數(shù)k= 1時x= 1 m處的Nt個溫度場測點數(shù)據(jù)作為樣本值,定義如下?lián)p失函數(shù):
基于PINN方法的導(dǎo)熱系數(shù)辨識結(jié)果如表5所示。從表5可以看出,當(dāng)熱傳導(dǎo)方程采用第二類邊界條件時,基于PINN方法的導(dǎo)熱系數(shù)的辨識精度與采用混合邊界條件熱傳導(dǎo)方程的辨識精度相當(dāng)。
表5 導(dǎo)熱系數(shù)k辨識結(jié)果Table 5 Identification results of the thermal conductivity k
針對一維熱傳導(dǎo)方程的正問題和逆問題,本文介紹了一種內(nèi)嵌物理機理神經(jīng)網(wǎng)絡(luò)的求解方法。該方法將神經(jīng)網(wǎng)絡(luò)作為解的近似,并在給定數(shù)據(jù)和偏微分方程的約束下,對網(wǎng)絡(luò)進行訓(xùn)練。通過將偏微分方程作為正則化部分加入到損失函數(shù)中,使得上述深度網(wǎng)絡(luò)具有較強的泛化性能以及一定程度上的可解釋性。
通過分析熱傳導(dǎo)方程求解的收斂誤差得到,采樣點較少時,采樣誤差是收斂誤差的主要組成部分,而當(dāng)采樣點較多時,主導(dǎo)誤差變?yōu)閮?yōu)化誤差。再通過對方程中的未知參數(shù)辨識的魯棒性進行分析得知,由于損失函數(shù)中包含了初始條件、邊界條件以及方程殘差構(gòu)成的正則項,因此,可以防止網(wǎng)絡(luò)對噪聲出現(xiàn)過擬合,從而增強辨識方法的魯棒性。但是,基于PINN的導(dǎo)熱系數(shù)辨識方法在給定不恰當(dāng)?shù)某瑓?shù)時,可能會導(dǎo)致訓(xùn)練發(fā)散。