林云云，鄭素佩，封建湖，靳 放

（長(zhǎng)安大學(xué) 理學(xué)院，西安 710064）

引言

Burgers 方程是一類常見的一維標(biāo)量雙曲守恒律方程，常被用來模擬均勻密度和均勻壓力下流體的動(dòng)量傳遞，可從Navier-Stokes 方程[1]的速度場(chǎng)中去掉壓力梯度項(xiàng)得到.該方程出現(xiàn)在應(yīng)用數(shù)學(xué)的許多領(lǐng)域，包括流體動(dòng)力學(xué)、氣體動(dòng)力學(xué)和交通流[2]等，因此對(duì)它的數(shù)值算法研究一直是一項(xiàng)非常有意義的工作.近年來，在數(shù)值近似領(lǐng)域出現(xiàn)了各種求解雙曲守恒律方程的高分辨率數(shù)值格式.鄭素佩等[3]基于移動(dòng)網(wǎng)格法提出了一種熵穩(wěn)定格式，該格式有效地提高了二維淺水波方程數(shù)值算法的分辨率.賈豆等[4]提出了一種旋轉(zhuǎn)通量混合格式，該格式提高了二維Euler 方程數(shù)值結(jié)果的分辨率.然而，這些數(shù)值格式仍具有傳統(tǒng)數(shù)值方法適應(yīng)性差、網(wǎng)格剖分、條件限制等問題.隨著深度學(xué)習(xí)在計(jì)算機(jī)視覺和自然語言處理領(lǐng)域的蓬勃發(fā)展，出現(xiàn)了替代傳統(tǒng)數(shù)值格式的無網(wǎng)格算法，使用深度神經(jīng)網(wǎng)絡(luò)求解方程.深度神經(jīng)網(wǎng)絡(luò)作為一種通用的非線性函數(shù)逼近器，通過優(yōu)化神經(jīng)網(wǎng)絡(luò)的參數(shù)來逼近方程的解.在20 多年之前，神經(jīng)網(wǎng)絡(luò)就被用于求解偏微分方程[5-6]，但受限于當(dāng)時(shí)的計(jì)算機(jī)能力，這一技術(shù)未能得到發(fā)展.近年來，隨著機(jī)器學(xué)習(xí)的飛速發(fā)展，這一技術(shù)被重新審視，高普陽等[7]利用卷積神經(jīng)網(wǎng)絡(luò)求解雙曲守恒律方程，通過優(yōu)化傳統(tǒng)有限體積法中的權(quán)重系數(shù)，得到了在粗網(wǎng)格下具有較高精度的新數(shù)值格式，但是該網(wǎng)絡(luò)選取的訓(xùn)練集是由經(jīng)典的數(shù)值格式得到的，因此也會(huì)保留經(jīng)典格式的不足.Raissi 等[8]提出了一種PINN 算法，通過自動(dòng)微分得到帶物理模型約束的神經(jīng)網(wǎng)絡(luò)，將其納入到神經(jīng)網(wǎng)絡(luò)損失函數(shù)中，這種結(jié)構(gòu)化方法引入了一種正則化機(jī)制，從而約束了可容許解的空間.該算法僅需少量的數(shù)據(jù)就可高效地求解偏微分方程解的正問題，以及從給定數(shù)據(jù)中推斷出未知參數(shù)的反問題.然而該算法在求解雙曲守恒律方程時(shí)在間斷附近會(huì)出現(xiàn)抹平現(xiàn)象，很難識(shí)別激波等間斷.Lu 等[9]建立了基于PINN 算法的PYTHON 庫DeepXDE，提出了一種類似于有限元細(xì)化的方法，它是一種新的基于殘差的自適應(yīng)細(xì)化方法.此方法是通過在間斷附近自適應(yīng)的增加訓(xùn)練點(diǎn)，以便很好地捕捉間斷，但需要經(jīng)過“小批量”的梯度下降來判斷間斷的位置，使訓(xùn)練過程變得復(fù)雜.Mao 等[10]利用PINN 算法近似模擬高速氣動(dòng)流動(dòng)的Euler 方程，在間斷附近對(duì)分散的訓(xùn)練點(diǎn)進(jìn)行了隨機(jī)聚類來捕獲解的間斷，以此來提高精度，防止誤差向整個(gè)區(qū)域傳播，但此方法需要事先知道間斷的位置，對(duì)于未知解的問題具有一定的挑戰(zhàn).以上文獻(xiàn)對(duì)于間斷問題的處理方法都是通過在間斷位置增加數(shù)據(jù)集來捕捉間斷，這會(huì)增加訓(xùn)練量，占用大量的存儲(chǔ)空間.LeFloch 等[11]提出了基于擴(kuò)散-色散正則化的數(shù)值方法，并且證明了該方法能很好地逼近非經(jīng)典解，但對(duì)于大振幅的非經(jīng)典激波無法逼近.El 等[12]對(duì)不同正則化方法進(jìn)行了分析，為雙曲守恒律理論和非線性色散波之間提供了明確的聯(lián)系.Clamond 等[13]提出了一種新的淺水(等熵Euler)方程正則化方法，這種正則化的特征是光滑激波的傳播速度與原始的不連續(xù)激波的傳播速度完全相同，但要找到一個(gè)具有所有所需性質(zhì)的正則項(xiàng)并不是一件輕而易舉的事.Minbashian 等[14]利用PINN 算法求解了非凸通量函數(shù)的雙曲守恒律方程，通過自動(dòng)微分將在極限范圍內(nèi)保持平衡的帶有擴(kuò)散和色散正則化的方程納入損失函數(shù)，通過考慮小規(guī)模的擴(kuò)散項(xiàng)和色散項(xiàng)近似求解方程，但是正則化方程中的擴(kuò)散和色散項(xiàng)的系數(shù)很難確定.受以上數(shù)值方法的啟發(fā)，本文提出了一種正則化的PINN反問題方法來逼近無黏Burgers 方程的間斷解.

本文基于PINN 算法，并結(jié)合文獻(xiàn)[14]中解決雙曲守恒律問題取決于潛在的正則化，提出了一種求解間斷問題的新方法.在本文中正則項(xiàng)為在極限范圍內(nèi)保持平衡的擴(kuò)散項(xiàng)，利用PINN 算法反向地求解擴(kuò)散正則化方程中的小參數(shù)，繼而用正則化方程的解逼近無黏Burgers 方程的間斷解.最后通過與未正則化方程的數(shù)值結(jié)果進(jìn)行比較，結(jié)果顯示新算法在間斷位置的分辨率更高.

1 擴(kuò)散正則化方程

考慮如下守恒律方程：

對(duì)于方程(1)的初值問題，即使初始條件是光滑的，其解也有可能在某一時(shí)刻發(fā)展成間斷，經(jīng)典解理論不再適用.因此，需要考慮方程的弱解，通常對(duì)于同一初值問題，弱解并不是唯一的，而滿足熵不等式的弱解可使解的唯一性得到保證.受文獻(xiàn)[14]以及Lax 的熵穩(wěn)定條件與“黏性消失”機(jī)制等價(jià)思想的啟發(fā)[15]，為避免非物理解的產(chǎn)生，在求解方程(1)時(shí)會(huì)人為地添加一些擴(kuò)散項(xiàng)來控制熵的變化，當(dāng)擴(kuò)散項(xiàng)系數(shù)趨于零時(shí)的極限解稱為黏性消失解，對(duì)應(yīng)的擴(kuò)散正則化方程可表示為

其中λ表示擴(kuò)散系數(shù)且λ>0，當(dāng)λ →0時(shí)得到方程(1).方程(4)在λ →0下的黏性消失解就是方程(1)的解.

2 間斷問題的PINN反問題求解算法

考慮到雙曲守恒律方程在間斷處會(huì)出現(xiàn)熵增，為了達(dá)到熵穩(wěn)定，需要添加擴(kuò)散項(xiàng)來控制熵的變化，本文通過使用PINN 算法解決擴(kuò)散正則化方程的反問題，訓(xùn)練出滿足物理約束條件的擴(kuò)散項(xiàng)，該方法能確定出擴(kuò)散系數(shù)，且數(shù)值結(jié)果不會(huì)出現(xiàn)抹平和偽振蕩現(xiàn)象.

2.1 網(wǎng)絡(luò)結(jié)構(gòu)

PINN 算法是一個(gè)簡(jiǎn)單的深度前饋神經(jīng)網(wǎng)絡(luò)，激活函數(shù)運(yùn)用雙曲正切函數(shù).基本結(jié)構(gòu)包括兩個(gè)共享超參數(shù)的神經(jīng)網(wǎng)絡(luò)，一個(gè)是將(x,t)以坐標(biāo)形式輸入的全連接神經(jīng)網(wǎng)絡(luò)，用它來逼近u(x,t)，另一個(gè)是通過對(duì)網(wǎng)絡(luò)函數(shù)應(yīng)用自動(dòng)微分[16]和算術(shù)運(yùn)算所得的帶物理模型的神經(jīng)網(wǎng)絡(luò)f(x,t).這兩個(gè)網(wǎng)絡(luò)對(duì)損失函數(shù)均有貢獻(xiàn)，最后，利用優(yōu)化算法獲得神經(jīng)網(wǎng)絡(luò)的參數(shù)和擴(kuò)散正則化方程的參數(shù).網(wǎng)絡(luò)結(jié)構(gòu)如圖1所示.

圖1 網(wǎng)絡(luò)結(jié)構(gòu)示意圖Fig.1 The diagram of the network structure

為了捕捉間斷，引入擴(kuò)散正則項(xiàng)，定義物理約束條件的神經(jīng)網(wǎng)絡(luò)如下：

在網(wǎng)絡(luò)訓(xùn)練過程中，定義損失函數(shù)如下：

其中Lu指損失函數(shù)中數(shù)據(jù)驅(qū)動(dòng)的部分，Lf指損失函數(shù)中物理模型驅(qū)動(dòng)的部分.損失函數(shù)均使用均方誤差來定義：

2.2 算法步驟

為了使算法更易用程序?qū)崿F(xiàn)，下面總結(jié)出擴(kuò)散正則化PINN 反問題算法的具體步驟：

步驟1選取訓(xùn)練集 T，以坐標(biāo)(x,t)的形式為神經(jīng)網(wǎng)絡(luò)的輸入量.

步驟2構(gòu)造深度神經(jīng)網(wǎng)絡(luò)模型uNN(x,t;θ)來近似方程(1)的解u(x,t)，其中θ={W,b}為參數(shù)，W為權(quán)重，b是偏差.

步驟3通過模型uNN(x,t;θ)計(jì)算網(wǎng)絡(luò)輸出，對(duì)模型uNN(x,t;θ)應(yīng)用自動(dòng)微分技術(shù)并將其代入方程(3)得到殘差網(wǎng)絡(luò)f(x,t).最后判斷網(wǎng)絡(luò)輸出值與訓(xùn)練集的均方誤差、殘差網(wǎng)絡(luò)的均方誤差之和是否小于理想誤差，如果小于理想誤差，停止訓(xùn)練，得到網(wǎng)絡(luò)，否則更新權(quán)重和偏差，重復(fù)步驟3，繼續(xù)訓(xùn)練網(wǎng)絡(luò)，直至小于理想誤差，找到最優(yōu)參數(shù)λ和雙曲守恒律方程初值問題的解.

2.3 模型訓(xùn)練

本文所有算例均使用TensorFlow 1.14 版本來訓(xùn)練.對(duì)網(wǎng)絡(luò)參數(shù)采用Xavier 進(jìn)行初始化，它可以保證輸入和輸出數(shù)據(jù)的方差一致，避免所有輸出值都趨于零.為了使訓(xùn)練加速，訓(xùn)練網(wǎng)絡(luò)時(shí)對(duì)輸入數(shù)據(jù)進(jìn)行歸一化處理，并使用Adam 優(yōu)化器和L-BFGS 優(yōu)化器對(duì)損失函數(shù)進(jìn)行優(yōu)化.其中，Adam 通過自適應(yīng)的學(xué)習(xí)率來加快收斂速度，L-BFGS 是一個(gè)全批量的優(yōu)化算法，具有二階收斂速度和占用內(nèi)存小的優(yōu)勢(shì).

3 數(shù)值算例

下面將以一維和二維無黏Burgers 方程為例驗(yàn)證新方法的性能.本文對(duì)正則化方程和原方程分別采用兩套網(wǎng)絡(luò)進(jìn)行訓(xùn)練和測(cè)試，為使結(jié)果具有可比性，這兩套網(wǎng)絡(luò)將針對(duì)同一問題使用相同的隱藏層數(shù)、神經(jīng)元個(gè)數(shù)和迭代步數(shù).

算例1滿足連續(xù)初始條件的一維無黏Burgers 方程

計(jì)算區(qū)間[?1,1]，方程(2)滿足如下初始條件：

u(x,0)=?sin(πx).

其精確解為

u(x,t)=?sin(πx?u(x,t)t),

滿足Neumann 邊界條件，時(shí)間計(jì)算到t=0.96.訓(xùn)練集的創(chuàng)建是從精確解的整個(gè)時(shí)間和空間中隨機(jī)選取N=2 000個(gè)點(diǎn)，使用此訓(xùn)練集訓(xùn)練一個(gè)8 層的神經(jīng)網(wǎng)絡(luò)，每個(gè)隱藏層都含有60 個(gè)神經(jīng)元.該網(wǎng)絡(luò)使用LBFGS 迭代50000 步，經(jīng)訓(xùn)練知：添加的擴(kuò)散項(xiàng)的系數(shù)λ=0.000 74，并預(yù)測(cè)了t=0.3,0.6,0.9時(shí)刻的解.如圖2所示，該算例的解在有限的時(shí)間內(nèi)會(huì)產(chǎn)生激波，可以直觀地看出，對(duì)于初值連續(xù)的方程，PINN 算法在未正則化方程的間斷位置雖能捕捉激波，但會(huì)出現(xiàn)振蕩，而新方法在間斷處無偽振蕩的產(chǎn)生且能很好地捕捉激波.

圖2 例1 的數(shù)值結(jié)果對(duì)比圖：(a)未正則化方程的精確解與預(yù)測(cè)解；(b)正則化方程的精確解與預(yù)測(cè)解Fig.2 Comparison of numerical results for example 1:(a)exact and predicted solutions to unregularized equations;(b)exact and predicted solutions to regularized equations

算例2滿足間斷初始條件的一維無黏Burgers 方程

計(jì)算區(qū)間[?1,1]，方程(2)滿足如下初始條件：

其精確解為

圖3 例2 的數(shù)值結(jié)果對(duì)比圖：(a)未正則化方程的精確解與預(yù)測(cè)解；(b)正則化方程的精確解與預(yù)測(cè)解Fig.3 Comparison of numerical results for example 2:(a)exact and predicted solutions to unregularized equations;(b)exact and predicted solutions to regularized equations

算例3滿足間斷初始條件的一維無黏Burgers 方程

其精確解為

滿足Neumann 邊界條件，時(shí)間計(jì)算到t=1.5.訓(xùn)練集的創(chuàng)建是從精確解的整個(gè)時(shí)間和空間中隨機(jī)選取N=3 000個(gè)點(diǎn)，使用該訓(xùn)練集訓(xùn)練一個(gè)8 層的神經(jīng)網(wǎng)絡(luò)，每個(gè)隱藏層都含有60 個(gè)神經(jīng)元.先使用Adam 迭代10000步，再使用L-BFGS 迭代50000 步，經(jīng)訓(xùn)練知：添加的擴(kuò)散項(xiàng)系數(shù)λ=0.000 34，并預(yù)測(cè)了t=0.45,0.90,1.35時(shí)刻的解.如圖4所示，該算例的精確解是由兩個(gè)激波組成的，隨著時(shí)間的演變，在某一時(shí)刻會(huì)合成一個(gè)激波.可以直觀地看出，PINN 算法對(duì)未正則化方程雙激波情形會(huì)出現(xiàn)嚴(yán)重的抹平現(xiàn)象且伴有振蕩，而新方法在多個(gè)間斷位置都能很好地捕捉激波.

圖4 例3 的數(shù)值結(jié)果對(duì)比圖：(a)未正則化方程的精確解與預(yù)測(cè)解；(b)正則化方程的精確解與預(yù)測(cè)解Fig.4 Comparison of numerical results for example 3:(a)exact and predicted solutions to unregularized equations;(b)exact and predicted solutions to regularized equations

算例4滿足連續(xù)初始條件的二維無黏Burgers 方程

計(jì)算區(qū)域[0,4]×[0,4]，方程(3)滿足初始條件u0(x)=0.5+sin(π(x+y)/2).采用周期邊界條件，時(shí)間計(jì)算到t=1.5/π.訓(xùn)練集的創(chuàng)建是從參考解(熵穩(wěn)定格式[19]網(wǎng)格加密)的整個(gè)時(shí)間和空間中隨機(jī)選取N=4 000個(gè)點(diǎn)，使用該訓(xùn)練集訓(xùn)練一個(gè)6 層的神經(jīng)網(wǎng)絡(luò)，每個(gè)隱藏層都含有40 個(gè)神經(jīng)元.先使用Adam 優(yōu)化器迭代50000 步，接著使用L-BFGS 迭代50000 步，經(jīng)訓(xùn)練知：添加的擴(kuò)散項(xiàng)系數(shù)λ=0.000 328，并預(yù)測(cè)了t=0.16,0.32,1.5/π時(shí)刻的解.如圖5所示，為了直觀比較，該算例畫出了等值線圖，可以看出對(duì)于連續(xù)初值問題，PINN 算法在未正則化方程和正則化方程的間斷位置捕捉激波的能力相當(dāng).

圖5 例4 的數(shù)值結(jié)果對(duì)比圖：(a)未正則化方程的精確解與預(yù)測(cè)解；(b)正則化方程的精確解與預(yù)測(cè)解Fig.5 Comparison of numerical results for example 4:(a)exact and predicted solutions to unregularized equations;(b)exact and predicted solutions to regularized equations

算例5滿足間斷初始條件的二維無黏Burgers 方程

計(jì)算區(qū)域[0,1]×[0,1]，方程(3)滿足如下初邊值條件：

其精確解為

滿足Neumann 邊界條件，時(shí)間計(jì)算到t=0.5.訓(xùn)練集的創(chuàng)建是從精確解的整個(gè)時(shí)空域中隨機(jī)選取N=5 000個(gè)點(diǎn)，使用該訓(xùn)練集訓(xùn)練一個(gè)9 層的神經(jīng)網(wǎng)絡(luò)，每個(gè)隱藏層含有60 個(gè)神經(jīng)元.先使用Adam 優(yōu)化器迭代50000步，接著使用L-BFGS 迭代50000 步，經(jīng)訓(xùn)練知：添加的擴(kuò)散項(xiàng)系數(shù)λ=0.000 08，并預(yù)測(cè)了t=0.25,0.35,0.45時(shí)刻的解.如圖6所示，該算例的精確解含有激波，可以直觀地看出，隨著時(shí)間的推進(jìn)，PINN 算法在未正則化方程的間斷位置出現(xiàn)嚴(yán)重抹平現(xiàn)象，而新方法能很好地捕捉激波.

圖6 例5 的數(shù)值結(jié)果對(duì)比圖：(a)未正則化方程的精確解與預(yù)測(cè)解；(b)正則化方程的精確解與預(yù)測(cè)解Fig.6 Comparison of numerical results for example 5:(a)exact and predicted solutions to unregularized equations;(b)exact and predicted solutions to regularized equations

4 結(jié)論

本文提出了一種擴(kuò)散正則化的PINN反問題求解算法，通過優(yōu)化損失函數(shù)來減少真解與數(shù)值解之間的誤差，訓(xùn)練出滿足物理約束條件的參數(shù)，從而得到了一種有效求解無黏Burgers 方程的新方法.數(shù)值試驗(yàn)表明，該方法對(duì)連續(xù)初值和間斷初值問題都能很好地捕捉激波.在具體的數(shù)值模擬中發(fā)現(xiàn)，在訓(xùn)練數(shù)據(jù)受到噪聲干擾的情況下，精度會(huì)有所降低，擴(kuò)散系數(shù)也隨之發(fā)生微小變化.如何在存在噪聲數(shù)據(jù)的情況下仍能得到高精度的數(shù)值結(jié)果，還需進(jìn)一步研究.在下一步的研究中，可將該神經(jīng)網(wǎng)絡(luò)模型得到的擴(kuò)散項(xiàng)作為構(gòu)造傳統(tǒng)的高分辨率數(shù)值格式的參考依據(jù)，另外，也可將這種正則化處理的方法應(yīng)用到一般的雙曲守恒律方程中.