基于兩重網(wǎng)格的深度學(xué)習(xí)方法求解定常偏微分方程?

彭湃，馮新龍

(新疆大學(xué) 數(shù)學(xué)與系統(tǒng)科學(xué)學(xué)院，新疆 烏魯木齊 830017)

0 引言

隨著深度學(xué)習(xí)研究的不斷深入,在不同的學(xué)科中產(chǎn)生了深刻變革,如圖像識(shí)別[1]、認(rèn)知科學(xué)[2]等.深度學(xué)習(xí)在處理復(fù)雜的系統(tǒng)模型中的預(yù)測(cè)任務(wù)時(shí),如何賦予深度學(xué)習(xí)這一強(qiáng)大函數(shù)近似器以先驗(yàn)知識(shí)成為學(xué)者的研究重點(diǎn).為了解決這個(gè)問(wèn)題,有些學(xué)者致力于設(shè)計(jì)專(zhuān)門(mén)的神經(jīng)網(wǎng)絡(luò)將給定的先驗(yàn)知識(shí)隱式嵌入到神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)之中,如利用卷積神經(jīng)網(wǎng)絡(luò)處理計(jì)算機(jī)視覺(jué)[3].有些學(xué)者旨在通過(guò)適當(dāng)?shù)膽土P神經(jīng)網(wǎng)絡(luò)的損失函數(shù)賦予神經(jīng)網(wǎng)絡(luò)先驗(yàn)知識(shí),如利用物理信息神經(jīng)網(wǎng)絡(luò)[4](Physics-Informed Neural Networks,PINN)求解偏微分方程[5].

2019年,Raissi等人提出了物理信息神經(jīng)網(wǎng)絡(luò)用于求解偏微分方程的正問(wèn)題及反問(wèn)題.隨后諸多學(xué)者根據(jù)不同的方程及其應(yīng)用背景對(duì)PINN進(jìn)行探索.在應(yīng)用領(lǐng)域,Jin等建立了速度-壓力、速度-渦度形式的NSFnets用于求解層流和槽道湍流問(wèn)題[6].此外,一些學(xué)者將傳統(tǒng)方法求解偏微分方程的技巧融入到PINN之中:如結(jié)合有限體積法在損失函數(shù)中增加界面約束做為正則化因子的DPINN[7],Kharazmi等基于Petrov-Galerkin變分原理提出了VPINN求解偏微分方程[8].通過(guò)分析前人的研究,不難發(fā)現(xiàn)PINN在求解偏微分方程時(shí)難以獲得高精度解.為了改善神經(jīng)網(wǎng)絡(luò)求解偏微分方程的計(jì)算精度,本文提出了修正物理信息神經(jīng)網(wǎng)絡(luò)(Rectified Physics-Informed Neural Networks,RPINN)的解決方案.將兩重網(wǎng)格[9]求解方程的思想融入到PINN求解PDE的框架中,將PINN方法得到數(shù)值解及其梯度傳遞到新的神經(jīng)網(wǎng)絡(luò)之中,并建立相應(yīng)的損失函數(shù),通過(guò)訓(xùn)練新的網(wǎng)絡(luò)從而得到PINN數(shù)值解的修正量.

1 預(yù)備知識(shí)

在這一部分,首先介紹前饋神經(jīng)網(wǎng)絡(luò)以及物理信息神經(jīng)網(wǎng)絡(luò)求解偏微分方程的主要思想,并詳細(xì)闡述了神經(jīng)網(wǎng)絡(luò)求解多目標(biāo)優(yōu)化問(wèn)題時(shí)采取的動(dòng)態(tài)權(quán)重策略.

1.1 物理信息神經(jīng)網(wǎng)絡(luò)

在數(shù)學(xué)意義上,前饋神經(jīng)網(wǎng)絡(luò)可以表示為一個(gè)復(fù)雜的非線(xiàn)性函數(shù)f(x;θ)L,其中L表示神經(jīng)網(wǎng)絡(luò)的層數(shù),x表示神經(jīng)網(wǎng)絡(luò)的輸入向量,θ表示神經(jīng)網(wǎng)絡(luò)的參數(shù),包含神經(jīng)網(wǎng)絡(luò)的權(quán)重和偏置.對(duì)于一個(gè)全連接的前饋神經(jīng)網(wǎng)絡(luò),可將其表示為(1)式:

其中:x∈Rn為輸入向量,Wi表示權(quán)重矩陣,bi為偏置向量,σ(·)表示神經(jīng)網(wǎng)絡(luò)的激活函數(shù).在神經(jīng)網(wǎng)絡(luò)中,激活函數(shù)的選取保證了網(wǎng)絡(luò)的非線(xiàn)性特性.在本文中,選取雙曲正切函數(shù)σ(·)=作為神經(jīng)網(wǎng)絡(luò)的激活函數(shù).

物理信息神經(jīng)網(wǎng)絡(luò)旨在求解出一個(gè)連續(xù)的數(shù)值解u(x)使之滿(mǎn)足方程和邊界條件,為此本文定義定常偏微分方程如(2)式所示:

其中:?表示偏微分方程的求解區(qū)域,??為?的邊界,x∈Rn表示空間向量,Nx表示非線(xiàn)性微分算子,u表示方程滿(mǎn)足的數(shù)值解,g(x)表示方程在邊界處滿(mǎn)足的函數(shù).

將神經(jīng)網(wǎng)絡(luò)表示的非線(xiàn)性函數(shù)f(x;θ)L代入到偏微分方程(2)式中,使之成為滿(mǎn)足方程的一個(gè)潛在解.該潛在解滿(mǎn)足殘差rθ(x):

同時(shí)在邊界處滿(mǎn)足:

非線(xiàn)性函數(shù)f(x;θ)L關(guān)于空間向量x的梯度通過(guò)自動(dòng)微分[10]求得.

若非線(xiàn)性函數(shù)f(x;θ)L滿(mǎn)足偏微分方程,則殘差rθ(x)=0,并且有(4)式成立.根據(jù)方程的殘差和邊界條件建立神經(jīng)網(wǎng)絡(luò)的損失函數(shù):

其中:

其中:Lf和Lb分別表示滿(mǎn)足方程和邊界的損失函數(shù).λf和λb分別表示損失函數(shù)Lf和Lb之間的權(quán)重因子,在優(yōu)化過(guò)程中有效避免了梯度病態(tài)現(xiàn)象.{0}表示滿(mǎn)足殘差rθ(x)上的訓(xùn)練集,維數(shù)為為滿(mǎn)足邊界的訓(xùn)練集,維數(shù)為Nb.

當(dāng)為神經(jīng)網(wǎng)絡(luò)建立損失函數(shù)L(x;θ)后,可通過(guò)梯度下降法、Adam算法、L-BFGS-B等一系列優(yōu)化算法優(yōu)化損失函數(shù)L(x;θ),從而得到一組最優(yōu)的神經(jīng)網(wǎng)絡(luò)參數(shù)θ?={W?,b?}.最后將物理信息神經(jīng)網(wǎng)絡(luò)求解得到的數(shù)值解表示為u?(x)=f(x;θ?)L.其具體求解流程如圖1所示.

1.2 動(dòng)態(tài)權(quán)重策略

在神經(jīng)網(wǎng)絡(luò)的訓(xùn)練過(guò)程中,損失函數(shù)Lf(x;θ)和Lb(x;θ)數(shù)值差異較大,在進(jìn)行反向傳播過(guò)程中,往往側(cè)重優(yōu)化Lf(x;θ),忽略在邊界上的損失Lb(x;θ),從而導(dǎo)致梯度病態(tài)現(xiàn)象.為了緩解這一現(xiàn)象,本文對(duì)神經(jīng)網(wǎng)絡(luò)的損失函數(shù)采取動(dòng)態(tài)權(quán)重策略.通過(guò)計(jì)算損失函數(shù)中Lf(x;θ)和Lb(x;θ)關(guān)于神經(jīng)網(wǎng)絡(luò)參數(shù)θ的梯度來(lái)動(dòng)態(tài)調(diào)整權(quán)重λb,從而達(dá)到平衡損失函數(shù)Lf(x;θ)和Lb(x;θ)的作用.

考慮動(dòng)態(tài)權(quán)重策略算法[11]如下:首先確定損失函數(shù)Lf的權(quán)重λf=c,其中c為常數(shù).當(dāng)采用梯度下降算法優(yōu)化損失函數(shù)

時(shí),神經(jīng)網(wǎng)絡(luò)參數(shù)θ的更新公式為:

其中:θk表示在第k次迭代過(guò)程中神經(jīng)網(wǎng)絡(luò)的參數(shù),超參數(shù)η表示神經(jīng)網(wǎng)絡(luò)的學(xué)習(xí)率,?(·)表示梯度算子,?θLf表示損失函數(shù)Lf(x;θ)關(guān)于參數(shù)θk的梯度值.然后計(jì)算中間變量

其中:μ是一個(gè)取值范圍為(0,1)的超參數(shù).同時(shí)定義未采取動(dòng)態(tài)權(quán)重策略為(D0)策略.

2 基于兩重網(wǎng)格的深度學(xué)習(xí)算法

兩重網(wǎng)格算法最初來(lái)源于對(duì)方程的迭代求解,它將細(xì)網(wǎng)格方程限制到粗網(wǎng)格上進(jìn)行求解,然后將所得的解延拓到細(xì)網(wǎng)格之上,最后與原來(lái)的解組合形成細(xì)網(wǎng)格上的精確解.本文提出的基于兩重網(wǎng)格的深度學(xué)習(xí)算法來(lái)源于兩重網(wǎng)格求解方程的思想,將PINN得到的數(shù)值解及其梯度傳遞到新的神經(jīng)網(wǎng)絡(luò)之中,通過(guò)為新的神經(jīng)網(wǎng)絡(luò)建立損失函數(shù)以此獲得數(shù)值解的修正量.

首先,將PINN求解偏微分方程得到的數(shù)值解u?(x)=f(x;θ?)L代入到偏微分方程(2)式中,可得(12)式:

(12)式表明由PINN計(jì)算得到的數(shù)值解u?(x)與真解u存在誤差ε=u?u?(x),接下來(lái)為誤差ε建立新的前饋神經(jīng)網(wǎng)絡(luò)ε(x;θ)L,使之滿(mǎn)足(13)式:

根據(jù)(13)式,構(gòu)建神經(jīng)網(wǎng)絡(luò)ε(x;θ)L與PINN數(shù)值解之間的殘差方程(x;t):

并根據(jù)(13)和(14)式,為神經(jīng)網(wǎng)絡(luò)ε(x;θ)L建立損失函數(shù)

其中:

和為權(quán)重因子.{0}表示神經(jīng)網(wǎng)絡(luò)ε(x;θ)L滿(mǎn)足殘差(x)的訓(xùn)練集,為提高計(jì)算精度,訓(xùn)練集的選取與神經(jīng)網(wǎng)絡(luò)f(x;θ)L的訓(xùn)練集相互獨(dú)立,{()}表示神經(jīng)網(wǎng)絡(luò)ε(x;θ)L在邊界處的訓(xùn)練集.在神經(jīng)網(wǎng)絡(luò)ε(x;θ)L的損失函數(shù)Lε(x;θ)中,u?()和u?()為神經(jīng)網(wǎng)絡(luò)f(x;θ)L傳遞到神經(jīng)網(wǎng)絡(luò)ε(x;θ)L的數(shù)值解,數(shù)值解u?關(guān)于向量x的梯度和二階導(dǎo)數(shù)通過(guò)自動(dòng)微分求解.

在本文中,采用Adam算法優(yōu)化神經(jīng)網(wǎng)絡(luò)的損失函數(shù)ε(x;θ)L,同時(shí)為了有效獲得神經(jīng)網(wǎng)絡(luò)的最優(yōu)參數(shù),采取學(xué)習(xí)率衰減策略,即經(jīng)過(guò)迭代步長(zhǎng)N次迭代后,學(xué)習(xí)率η衰減為原來(lái)的ρ倍.通過(guò)對(duì)神經(jīng)網(wǎng)絡(luò)不斷優(yōu)化,得到滿(mǎn)足數(shù)值解u?(x)的修正量ε?=ε(x;θ)L.從而將RPINN求解偏微分方程的數(shù)值解表示為:uR=u?(x)+ε(x;θ)L,具體求解流程如圖2所示.

3 數(shù)值算例

為了驗(yàn)證結(jié)合動(dòng)態(tài)權(quán)重策略的深度學(xué)習(xí)方法在提高計(jì)算精度上的有效性,本文給出了若干數(shù)值實(shí)驗(yàn),分別為對(duì)流擴(kuò)散反應(yīng)方程和三維Navier-Stokes方程的求解,并且定義相對(duì)L2誤差為

例1對(duì)流擴(kuò)散反應(yīng)方程屬于流體力學(xué)中的基本方程.大氣污染、河流污染等諸多變化過(guò)程都可以用對(duì)流擴(kuò)散反應(yīng)方程來(lái)表示.對(duì)于求解此類(lèi)方程,常用的求解方法有有限元法、多重網(wǎng)格法.下面本文通過(guò)求解一個(gè)二維的對(duì)流擴(kuò)散反應(yīng)方程來(lái)驗(yàn)證基于兩重網(wǎng)格的深度學(xué)習(xí)算法的有效性.定義對(duì)流擴(kuò)散反應(yīng)方程如下所示[12]

其中:α=0.1是擴(kuò)散項(xiàng)參數(shù),β=(2,3),f是滿(mǎn)足方程的源項(xiàng),求解區(qū)域?yàn)?=[0,1]2.定義方程滿(mǎn)足的真解為:

在PINN求解對(duì)流擴(kuò)散反應(yīng)方程時(shí),可根據(jù)方程滿(mǎn)足的殘差及其邊界建立神經(jīng)網(wǎng)絡(luò)的損失函數(shù):

其中:

Nf,Nb分別為在求解區(qū)域?以及邊界上選取訓(xùn)練點(diǎn)的個(gè)數(shù).

在本算例中,PINN方法中神經(jīng)網(wǎng)絡(luò)f(x;θ)L的參數(shù)設(shè)置如下:神經(jīng)網(wǎng)絡(luò)的網(wǎng)絡(luò)結(jié)構(gòu)為2+[50]×4+1,即輸入層含有2個(gè)輸入神經(jīng)元,中間隱藏層的個(gè)數(shù)為4,每層包含有50個(gè)神經(jīng)元,輸出層含有1個(gè)輸出神經(jīng)元.訓(xùn)練集的維數(shù)為:Nf=10 500,Nb=2 800.初始學(xué)習(xí)率η設(shè)置為η=10?3,衰減步長(zhǎng)為N=1 000,學(xué)習(xí)率的衰減率為ρ=0.9.動(dòng)態(tài)權(quán)重策略中采取的超參數(shù)為μ=0.5,λf=1,采取Adam優(yōu)化算法的迭代次數(shù)為8 000，并采用L-BFGS-B算法加速收斂.

在RPINN方法中神經(jīng)網(wǎng)絡(luò)f(x;θ)L的參數(shù)設(shè)置如下:神經(jīng)網(wǎng)絡(luò)的網(wǎng)絡(luò)結(jié)構(gòu)為2+[50]×4+1,訓(xùn)練集的維數(shù)為:Nf=8 000,Nb=2 000.初始學(xué)習(xí)率η設(shè)置為η=10?3,衰減步長(zhǎng)為N=1 000,學(xué)習(xí)率的衰減率為ρ=0.9.動(dòng)態(tài)權(quán)重策略中采取的超參數(shù)為μ=0.5,λf=1,采取Adam優(yōu)化算法的迭代次數(shù)為3 000.并將PINN方法訓(xùn)練得到滿(mǎn)足偏微分方程的數(shù)值解記作.

神經(jīng)網(wǎng)絡(luò)ε(x;θ)L的參數(shù)設(shè)置如下:神經(jīng)網(wǎng)絡(luò)的網(wǎng)絡(luò)結(jié)構(gòu)為2+[50]×4+1,訓(xùn)練集的維數(shù)為:=2 500,=800.初始學(xué)習(xí)率η設(shè)置為η=10?3,學(xué)習(xí)率的衰減步長(zhǎng)為N=1 000,學(xué)習(xí)率的衰減率為ρ=0.9.動(dòng)態(tài)權(quán)重策略中采取的超參數(shù)為μ=0.5,=1,采取Adam優(yōu)化算法進(jìn)行15 000次迭代,并采取L-BFGS-B算法加速網(wǎng)絡(luò)的收斂.最后將神經(jīng)網(wǎng)絡(luò)ε(x;θ)L得到的修正量記作ε?=ε(x;θ)L,從而RPINN方法得到數(shù)值解為u=+ε?.

下面本文將多重網(wǎng)格、PINN、RPINN方法得到的數(shù)值解與真解進(jìn)行對(duì)比分析.值得說(shuō)明的是:在多重網(wǎng)格求解方程時(shí),迭代方法采用SOR迭代,PINN方法選取訓(xùn)練集的維數(shù)為13 300,具體實(shí)驗(yàn)數(shù)據(jù)如表1所示.通過(guò)分析可知:與多重網(wǎng)格相比，RPINN方法利用更少的數(shù)據(jù)獲得高精度的數(shù)值解.同時(shí),采取動(dòng)態(tài)權(quán)重策略(D1)能夠有效緩解梯度病態(tài)現(xiàn)象,得到高精度的數(shù)值解.RPINN與PINN方法相比,計(jì)算精度提高了3至8倍.關(guān)于RPINN(D1)、PINN(D1)、多重網(wǎng)格求解得到的數(shù)值解與真解的圖像如圖3所示.

表1 多重網(wǎng)格、PINN和RPINN求解對(duì)流擴(kuò)散反應(yīng)方程之間的數(shù)據(jù)對(duì)比

然后,本文給出了RPINN求解對(duì)流擴(kuò)散反應(yīng)方程時(shí),采取動(dòng)態(tài)權(quán)重策略(D0)、(D1)的RPINN關(guān)于神經(jīng)網(wǎng)絡(luò)參數(shù)θ的梯度的直方分布圖,分別如圖4和圖5所示.結(jié)果表明結(jié)合動(dòng)態(tài)權(quán)重策略的RPINN中的值均趨向于0，并且數(shù)值之間差異小，驗(yàn)證了動(dòng)態(tài)權(quán)重策略能有效緩解梯度病態(tài)現(xiàn)象.

最后,本文分析了RPINN的神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)對(duì)修正量ε?訓(xùn)練結(jié)果的影響.首先以神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)為2+[50]×4+1的PINN訓(xùn)練得到的數(shù)值解作為基準(zhǔn),通過(guò)改變第二個(gè)神經(jīng)網(wǎng)絡(luò)的結(jié)構(gòu)對(duì)同一數(shù)值解進(jìn)行修正.具體的神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)設(shè)計(jì)和誤差如表2所示.結(jié)果表明:隨著神經(jīng)網(wǎng)絡(luò)的網(wǎng)絡(luò)層數(shù)和每層神經(jīng)元個(gè)數(shù)的增加,對(duì)修正量ε?訓(xùn)練效果越好.

表2 不同網(wǎng)絡(luò)結(jié)構(gòu)對(duì)數(shù)值精度影響的比較

例2Navier-Stokes方程是一類(lèi)描述流體流動(dòng)的非線(xiàn)性方程,在流體力學(xué)中占有重要地位.本算例中,定義穩(wěn)態(tài)不可壓縮的Navier-Stokes方程如(26)式所示:

其中:u=[u,v,w]T是無(wú)量綱的速度向量,p是無(wú)量綱壓力,雷諾數(shù)Re=UL/v,其中U和L分別是特征速度和長(zhǎng)度,v是流體的運(yùn)動(dòng)粘度,f=[f1,f2,f3]為滿(mǎn)足方程的源項(xiàng),Γ表示狄利克雷邊界條件.求解區(qū)域?=[?1,1]3,定義方程滿(mǎn)足的真解為：

其中:a=d=1.

在本算例中,PINN方法中神經(jīng)網(wǎng)絡(luò)參數(shù)設(shè)置如下:神經(jīng)網(wǎng)絡(luò)的網(wǎng)絡(luò)結(jié)構(gòu)為3+[50]×4+4,即輸入層包含3個(gè)神經(jīng)元,中間隱藏層的個(gè)數(shù)為4,每層包含50個(gè)神經(jīng)元,輸出層包含4個(gè)神經(jīng)元.訓(xùn)練集的維數(shù)為:Nf=11 375,Nb=4 200.初始學(xué)習(xí)率η設(shè)置為η=10?3,衰減步長(zhǎng)為N=1 000,學(xué)習(xí)率的衰減率為ρ=0.9.在動(dòng)態(tài)權(quán)重策略中采取的超參數(shù)為μ=0.5,λf=1,Adam算法的迭代次數(shù)為8 000,采用L-BFGS-B算法加速收斂.

在RPINN方法中神經(jīng)網(wǎng)絡(luò)f(x;θ?)L的參數(shù)設(shè)置如下:神經(jīng)網(wǎng)絡(luò)的網(wǎng)絡(luò)結(jié)構(gòu)為3+[50]×4+4，訓(xùn)練集的維數(shù)為:Nf=8 000,Nb=3 000.初始學(xué)習(xí)率η設(shè)置為η=10?3，衰減步長(zhǎng)為N=1 000,學(xué)習(xí)率衰減率為ρ=0.9.在動(dòng)態(tài)權(quán)重策略中采取的超參數(shù)為μ=0.5,λf=1,Adam 算法的迭代次數(shù)為4 000.

神經(jīng)網(wǎng)絡(luò)ε(x;θ)L的參數(shù)設(shè)置如下:神經(jīng)網(wǎng)絡(luò)的網(wǎng)絡(luò)結(jié)構(gòu)為3+[40]×4+4,訓(xùn)練集的維數(shù)為:=3 375,=1 200.初始學(xué)習(xí)率η設(shè)置為η=10?3,衰減步長(zhǎng)為N=1 000,學(xué)習(xí)率的衰減率為ρ=0.9.動(dòng)態(tài)權(quán)重策略中采取的超參數(shù)為μ=0.5,=1,Adam算法的迭代次數(shù)為12 000,采用L-BFGS-B算法加速收斂.經(jīng)過(guò)對(duì)神經(jīng)網(wǎng)絡(luò)ε(x;θ)L不斷訓(xùn)練,最終得到數(shù)值解的修正量為ε?=ε(x;θ)L,從而將RPINN方法得到數(shù)值解記為u=+ε?.

最后比較了PINN、RPINN得到數(shù)值解在z=0平面上與真解的誤差,如圖6所示.實(shí)驗(yàn)結(jié)果表明,PINN方法和RPINN方法在求解偏微分方程時(shí),與求解區(qū)域內(nèi)部相比,在邊界上的誤差偏大.有關(guān)速度場(chǎng)的數(shù)值解與真解的誤差如表3所示.通過(guò)比較分析可知,RPINN方法與PINN方法相比,能夠有效提高數(shù)值精度.

表3 PINN、RPINN方法求解Navier-Stokes方程的數(shù)值比較

4 結(jié)論

本文根據(jù)兩重網(wǎng)格求解偏微分方程的思想,對(duì)物理信息神經(jīng)網(wǎng)絡(luò)進(jìn)行改進(jìn),提出了基于兩重網(wǎng)格的深度學(xué)習(xí)方法求解定常偏微分方程.然后針對(duì)優(yōu)化過(guò)程中出現(xiàn)的梯度病態(tài)現(xiàn)象,采取了動(dòng)態(tài)權(quán)重策略.最后通過(guò)若干數(shù)值實(shí)驗(yàn)驗(yàn)證,本文提出的基于兩重網(wǎng)格的深度學(xué)習(xí)方法與原始的物理信息神經(jīng)網(wǎng)絡(luò)相比,在計(jì)算精度上提高了3至8倍.關(guān)于如何設(shè)計(jì)出快速穩(wěn)定的深度學(xué)習(xí)方法以減少計(jì)算時(shí)間以及進(jìn)一步提高數(shù)值精度將是以后研究工作的重點(diǎn).