應(yīng)用自適應(yīng)節(jié)點(diǎn)生成物理信息網(wǎng)絡(luò)計(jì)算地震波旅行時(shí)

摘要：求解程函方程能夠獲得震源定位、層析成像等地球物理反演所需的地震波旅行時(shí)，常用算法包括快速推進(jìn)法（FMM）和快速掃描法（FSM）等。物理信息神經(jīng)網(wǎng)絡(luò)（PINN）是一種新穎的無(wú)網(wǎng)格方法，可將偏微分方程中的微分形式約束條件融入到神經(jīng)網(wǎng)絡(luò)的損失函數(shù)中，從而獲得帶物理信息約束的神經(jīng)網(wǎng)絡(luò)。文中聚焦訓(xùn)練過(guò)程中的節(jié)點(diǎn)優(yōu)化配置，采用基于殘差分布的自適應(yīng)采樣方法改善PINN的訓(xùn)練效果，提出了基于自適應(yīng)節(jié)點(diǎn)生成的物理信息網(wǎng)絡(luò)旅行時(shí)計(jì)算方法。Marmousi模型和起伏地表模型的測(cè)試結(jié)果均表明，該方法相較于固定節(jié)點(diǎn)生成方法具有更穩(wěn)定的訓(xùn)練過(guò)程并且旅行時(shí)計(jì)算結(jié)果能保持較高的精度。

關(guān)鍵詞：旅行時(shí)計(jì)算，程函方程，物理信息神經(jīng)網(wǎng)絡(luò)（PINN），殘差，自適應(yīng)節(jié)點(diǎn)生成，無(wú)網(wǎng)格方法中圖分類號(hào)：P631 文獻(xiàn)標(biāo)識(shí)碼：A DOI：10.13810/j.cnki.issn.1000-7210.20240086

Traveltime calculation method with physics-informed neural network byresidual-based adaptive node generation

TANG Jie12， WANG Haicheng12，F(xiàn)AN Zhonghao12，PAN Deng12，REN Limin3，ZHANG Jingdong （1.SchoolofGeosciencesChina UniversityofPetroleu（EastChia）Qingdao，hndong2658O，Chna;2.Laboratoryforie MineralResourcesQingdaoNationalLaboratoryforMarineScienceand TechnologyQingdao，Shandong 266071，China ;3.BGP Equipment，CNC，Zhuozhou，HebeiO751，China;4.GeologicalResearchCenterGPIc.，ZhuozhouCNPC，HebeiOCina）

Abstract： Seismic wave traveltime information required for geophysical inversion such as source localization and tomographic imaging can be obtained by solving the eikonal：equation. Common algorithms for this purpose include the fast marching method （FMM） and the fast sweeping method （FSM）.Physics-informed neural network （PINN） is a novel mesh-free method that incorporates the constraints of partial differential equations into the loss function of the neural network，thus embedding physical information into the network. Focusing on optimizing node distribution during training，this study adopts an adaptive sampling strategy based on residual distributio to improve the training performance of PINN and proposes a travel-time calculation method using PINN with adaptive node generation. Application tests on the Marmousi model and an iregular topography model show that， compared with the fixed node-generation method，the proposed approach yields a more stable training process and maintains high accuracy in traveltime calculations.

Keywords：traveltime calculation，eikonal equation，physics-informed neural network（PINN），residual，adaptive nodegeneration，mesh-free method

唐杰，王海成，范忠豪，等，應(yīng)用自適應(yīng)節(jié)點(diǎn)生成物理信息網(wǎng)絡(luò)計(jì)算地震波旅行時(shí)[J].石油地球物理勘探， 2025，60（4） ：840-851. TANG Jie， WANG Haicheng，F(xiàn)AN Zhonghao，et al.Traveltime calculation method with physics-informed neural network by residual-based adaptive node generation[J]. Oil Geophysical Prospecting，2025，60（4）：840-851.

0 引言

波初至旅行時(shí)，是地震資料處理和成像所需的重要信息[1-3]。傳統(tǒng)的程函方程求解主要有兩種常用算法，分別是快速推進(jìn)法（FastMarchingMethod，F(xiàn)MM）和快速掃描法（FastSweepingMethod，程函方程（EikonalEquation）常被用來(lái)計(jì)算地震

FSM）。Sethian4采用迎風(fēng)差分格式離散程函方程，再通過(guò)FMM計(jì)算旅行時(shí)分布。FSM將旅行時(shí)場(chǎng)傳播方向分成多組，對(duì)每組分別利用Gauss-Seidel迭代方法求解離散化后的方程[5-6]。這兩種方法都能夠計(jì)算初至旅行時(shí)，但是對(duì)于結(jié)構(gòu)復(fù)雜模型，在精度和計(jì)算效率方面仍有提升空間。

近年來(lái)，深度學(xué)習(xí)方法得到了廣泛的研究及應(yīng)用。神經(jīng)網(wǎng)絡(luò)可以看作一個(gè)通用的非線性函數(shù)逼近器[]，使用神經(jīng)網(wǎng)絡(luò)求解偏微分方程（Partial Dif-ferentialEquation，PDE）亦成為研究熱點(diǎn)[8-9]。Raissi等[10開(kāi)發(fā)了基于物理信息神經(jīng)網(wǎng)絡(luò)的深度學(xué)習(xí)框架，可用來(lái)求解PDE。Smith等提出了一種求解程函方程的深度學(xué)習(xí)網(wǎng)絡(luò)EikoNet，能夠獲得三維速度結(jié)構(gòu)中的初至旅行時(shí)場(chǎng)，但是在高維情況下面臨計(jì)算效率和精度的挑戰(zhàn)。Waheed等12將因式分解后的程函方程結(jié)合物理信息網(wǎng)絡(luò)計(jì)算旅行時(shí)，需利用局部自適應(yīng)激活函數(shù)提高網(wǎng)絡(luò)收斂效果和旅行時(shí)求解精度，如果激活函數(shù)設(shè)計(jì)不合適，會(huì)影響在復(fù)雜地形條件下收斂效果。Grubas等[3提出一種求解程函方程的全連接網(wǎng)絡(luò)架構(gòu)，可以在存在焦散情況下保證旅行時(shí)求解的穩(wěn)定性，但計(jì)算復(fù)雜度較高。

物理信息神經(jīng)網(wǎng)絡(luò)（Physics-InformedNeuralNetworks，PINN）將非線性微分方程編碼為神經(jīng)網(wǎng)絡(luò)框架的一個(gè)組成部分，將PDE中的微分形式作為約束條件融入神經(jīng)網(wǎng)絡(luò)的損失函數(shù)，從而獲得帶物理機(jī)制約束的神經(jīng)網(wǎng)絡(luò)[14]。這樣訓(xùn)練出來(lái)的神經(jīng)網(wǎng)絡(luò)不僅能夠逼近觀測(cè)數(shù)據(jù)，而且能夠滿足PDE所遵循的物理規(guī)律[15]。PINN求解程函方程是一種無(wú)網(wǎng)格方法，需要選取節(jié)點(diǎn)進(jìn)行訓(xùn)練，節(jié)點(diǎn)的選取方法會(huì)影響網(wǎng)絡(luò)的應(yīng)用效果，以往的研究對(duì)于節(jié)點(diǎn)選取方法分析較少。近年來(lái)國(guó)內(nèi)外學(xué)者開(kāi)發(fā)了許多節(jié)點(diǎn)生成方法，如基于重要性的采樣和加權(quán)采樣方法等[16-17]。Das 等[18]比較了PINN的幾種常用節(jié)點(diǎn)生成方法，結(jié)果表明Hammersley采樣方法與其他方法相比具有明顯優(yōu)勢(shì)。Peng等[19]引入了具有局部規(guī)則性的變密度節(jié)點(diǎn)生成方法，能夠自適應(yīng)地生成具有可變密度的局部準(zhǔn)均勻節(jié)點(diǎn)。但這兩種節(jié)點(diǎn)方法并未用于物理信息網(wǎng)絡(luò)的旅行時(shí)計(jì)算。此外，物理信息網(wǎng)絡(luò)旅行時(shí)計(jì)算在訓(xùn)練階段通常保持訓(xùn)練節(jié)點(diǎn)位置不變，沒(méi)有考慮訓(xùn)練過(guò)程中殘差的空間分布特征，導(dǎo)致部分區(qū)域的計(jì)算結(jié)果存在較大偏差。

針對(duì)PINN計(jì)算旅行時(shí)固定節(jié)點(diǎn)分布導(dǎo)致的局部誤差過(guò)大的問(wèn)題，本文提出了一種基于殘差分布的自適應(yīng)節(jié)點(diǎn)細(xì)化調(diào)整方案。對(duì)于復(fù)雜速度模型，基于該方案使用PINN求解程函方程，能夠有效獲得高精度的旅行時(shí)結(jié)果，模型測(cè)試驗(yàn)證了該方法的有效性。

1理論

1.1 程函方程因式分解

程函方程可用于計(jì)算高頻近似下的地震波旅行時(shí)。對(duì)于各向同性速度模型，程函方程可表示為

式中： τ 為地震波旅行時(shí)； ?=（?_x，?_y） 為空間坐標(biāo) x 處的梯度算子； v（x） 為速度模型。該方程定義了受震源位置 x_s=（x_s，y_s） 處旅行時(shí) τ（x_s）=0 約束的固定震源在速度模型 v（x） 中的初至波旅行時(shí)場(chǎng)。

具有單個(gè)隱層和有限個(gè)神經(jīng)元的神經(jīng)網(wǎng)絡(luò)可以用來(lái)表示任意精度的有界連續(xù)函數(shù)。對(duì)于程函方程而言，物理信息網(wǎng)絡(luò)將模型的空間位置信息以及速度信息轉(zhuǎn)化為旅行時(shí)信息。

程函方程的解在源點(diǎn)處的奇異性會(huì)降低該區(qū)域解的精度[20-21]。使用因式分解方法可以避免這一問(wèn)題[22-23]。將旅行時(shí)分解為兩個(gè)函數(shù)，其中一個(gè)函數(shù)是解析指定的，另一個(gè)函數(shù)在源鄰域中是光滑的，即

式中： R=|x-x_s| 為距離變量，即模型速度為1時(shí)的解 ;f（?） 是神經(jīng)網(wǎng)絡(luò)輸出 ，f（x） 表示位置 x 處的神經(jīng)網(wǎng)絡(luò)預(yù)測(cè)的結(jié)果； γ（?） 是邊界函數(shù)，即非均勻性導(dǎo)致的小偏差。由于做了因式分解，所以將 和 R 相乘后得到 x 處的旅行時(shí)。

考慮到旅行時(shí)必須為正，即 τ（x）gt;0，?x≠x_sc （204號(hào)由于 R?0 ，則邊界函數(shù) γ 必須始終是正的。假設(shè)單位速度為最小值，那么 R 限制了 γ 值域?yàn)?（0，1] 任何速度模型都有速度范圍 v_minmax ，即允許用范圍 [1/v_max，1/v_min] 約束 γ ，則限制了旅行時(shí)解的范圍 τ∈[R/v_max，R/v_min]， ，即限制在最快和最慢的解之間，這種方法能夠加快物理信息網(wǎng)絡(luò)的訓(xùn)練效率。

邊界函數(shù)滿足

（3）式中 σ 可以是在[0，1]范圍內(nèi)定義的任意函數(shù)，本文選取

1.2基于殘差分布的自適應(yīng)節(jié)點(diǎn)調(diào)整

PINN將PDE作為神經(jīng)網(wǎng)絡(luò)的物理約束。與經(jīng)典計(jì)算方法相比，需要平衡神經(jīng)網(wǎng)絡(luò)的逼近誤差和泛化誤差。其中逼近誤差主要來(lái)自神經(jīng)網(wǎng)絡(luò)的逼近能力，泛化誤差是模型對(duì)未知數(shù)據(jù)的預(yù)測(cè)能力，主要與訓(xùn)練中參與計(jì)算損失函數(shù)的節(jié)點(diǎn)有關(guān)[24-25]。PINN在選擇殘差計(jì)算節(jié)點(diǎn)方面有很大的靈活性，以往的做法是在訓(xùn)練開(kāi)始時(shí)選擇固定的節(jié)點(diǎn)位置，可以是規(guī)則網(wǎng)格點(diǎn)或隨機(jī)點(diǎn)，在訓(xùn)練過(guò)程中不再更改它們的位置。殘差通常在計(jì)算域中隨機(jī)分布，但對(duì)于某些區(qū)域呈現(xiàn)陡峭梯度的PDE，殘差會(huì)出現(xiàn)分布異常，需要在陡峭梯度區(qū)域及附近布設(shè)更多的節(jié)點(diǎn)，即對(duì)殘差較高的區(qū)域進(jìn)行細(xì)化，通過(guò)引入具有局部規(guī)則性的變密度節(jié)點(diǎn)生成方法改進(jìn)訓(xùn)練節(jié)點(diǎn)的分布，能夠在使用相同數(shù)量訓(xùn)練節(jié)點(diǎn)的情況下獲得更精確的結(jié)果。

固定節(jié)點(diǎn)的生成方法包括隨機(jī)采樣和Hammersley采樣。隨機(jī)采樣方法在指定范圍內(nèi)隨機(jī)生成點(diǎn)，由于點(diǎn)的位置隨機(jī)，可能會(huì)出現(xiàn)點(diǎn)聚集或空白區(qū)域。隨機(jī)采樣方法適用于對(duì)點(diǎn)分布均勻性要求不高的情況。Hammersley采樣是一種準(zhǔn)隨機(jī)采樣方法，生成的點(diǎn)能更均勻地覆蓋整個(gè)采樣區(qū)域，并且在空間分布上更有規(guī)律，能夠有效減少采樣點(diǎn)之間的間隙和重疊。

給定樣本的總數(shù)為 N ，可以生成在 [0，1]× [0，1]分布的二維點(diǎn)集。Hammersley法生成的二維點(diǎn)集滿足

式中： 為樣點(diǎn)序號(hào)， a_m 是 i 用二進(jìn)制表示時(shí)的各個(gè)系數(shù)，即 i= a₀+a₁×2+a₁×2²+…+a_n×2^m 。根據(jù)具體的模型尺寸可以對(duì)節(jié)點(diǎn)與節(jié)點(diǎn)之間的距離進(jìn)行縮放。

基于殘差分布的自適應(yīng)節(jié)點(diǎn)調(diào)整方案將殘差分布與計(jì)算域中動(dòng)態(tài)改變采樣半徑的方法相結(jié)合，能夠在訓(xùn)練過(guò)程中基于殘差分布動(dòng)態(tài)改變節(jié)點(diǎn)的位置和密度，進(jìn)而能夠局部細(xì)化高殘差區(qū)域上的節(jié)點(diǎn)密度，改進(jìn)PINN訓(xùn)練效果。本文采用Fornbergamp;Flyer法生成基于殘差分布的自適應(yīng)節(jié)點(diǎn)。該方法是一種變密度的網(wǎng)格節(jié)點(diǎn)生成方法，在殘差大處節(jié)點(diǎn)密集，在殘差小處局部節(jié)點(diǎn)稀疏2，即節(jié)點(diǎn)間距隨殘差的變化調(diào)整，同時(shí)能避免節(jié)點(diǎn)分布過(guò)稀或過(guò)密。變密度網(wǎng)格節(jié)點(diǎn)生成的步驟如下。

（1）如圖1a所示，隨機(jī)將矩形中的一些節(jié)點(diǎn)初始化為潛在點(diǎn)集 A ，并將所選節(jié)點(diǎn)集 B 初始化為空集。隨后在潛在節(jié)點(diǎn)集 A 中找出 y 軸方向最低的節(jié)點(diǎn)標(biāo)紅，并將該點(diǎn)作為一個(gè)新的點(diǎn)添加到節(jié)點(diǎn)集 B 中。

（2）如圖1b所示，以在步驟（1）中所選中的節(jié)點(diǎn)為中心，通過(guò)其左右兩側(cè)的兩個(gè)相鄰點(diǎn)確定一個(gè)圓，圓的半徑可控制生成節(jié)點(diǎn)的密度。

（3）如圖1c所示，保留作為圓心的節(jié)點(diǎn)，刪除該圓內(nèi)其他所有節(jié)點(diǎn)。在該圓的上半圓周上生成5個(gè)新的等間距節(jié)點(diǎn)，并將它們添加到點(diǎn)集 A 。

（4）重復(fù)上述（2）、（3）步，直到所選的節(jié)點(diǎn)數(shù)超過(guò)上限。

圖1Fornbergamp;Flyer法節(jié)點(diǎn)生成示意圖

節(jié)點(diǎn)處殘差的絕對(duì)值越大，表示節(jié)點(diǎn)附近的區(qū)域越需要細(xì)化。由于殘差的范圍可能變化很大，首先將殘差值標(biāo)準(zhǔn)化為 [0，1） 之間的數(shù)ε。采樣半徑函數(shù)直接控制Fornbergamp;Flyer算法在每個(gè)位置的生成半徑 h

采樣半徑 h 能夠控制局部節(jié)點(diǎn)密度，而采樣參數(shù)r 和 s 決定了最大和最小生成半徑。比率參數(shù) r 控制最大與最小半徑的比率，從而確定最大與最小密度的比率，且 。縮放參數(shù) s 控制半徑的總體縮放。標(biāo)準(zhǔn)化殘差 控制采樣半徑函數(shù) h ，殘差大的位置的采樣半徑較小，節(jié)點(diǎn)密度更高；殘差小的位置處的采樣半徑較大，節(jié)點(diǎn)密度較低。

根據(jù)網(wǎng)絡(luò)輸出的殘差調(diào)整節(jié)點(diǎn)生成半徑，從而在殘差值較高的區(qū)域生成更多節(jié)點(diǎn)，在殘差值較小的區(qū)域生成較少節(jié)點(diǎn)。圖2為三種不同的采樣方式生成1000個(gè)節(jié)點(diǎn)獲得的節(jié)點(diǎn)分布，節(jié)點(diǎn)分布的密度大小由核密度估計(jì)值體現(xiàn)。核密度估計(jì)值通過(guò)高斯核平滑技術(shù)，估計(jì)每個(gè)點(diǎn)周圍的概率密度，其值越大，表示該點(diǎn)附近的樣本點(diǎn)越密集。當(dāng)然，從樣本點(diǎn)的分布也可以看出不同位置樣本點(diǎn)的疏密程度。可見(jiàn)，除了隨機(jī)節(jié)點(diǎn)分布以外，Hammersely方法生成的節(jié)點(diǎn)分布比較均勻，但仍存在局部聚集現(xiàn)象，而Fornbergamp;Flyer方法生成的節(jié)點(diǎn)不僅更加均勻，而且節(jié)點(diǎn)的生成半徑可以根據(jù)具體需求調(diào)整，是最理想的方法之一。

圖2三種方法產(chǎn)生的節(jié)點(diǎn)分布對(duì)比

1.3 自適應(yīng)節(jié)點(diǎn)調(diào)整的物理信息網(wǎng)絡(luò)

本文將自適應(yīng)節(jié)點(diǎn)調(diào)整策略應(yīng)用于求解程函方程的物理信息網(wǎng)絡(luò)，提出了一種基于殘差分布的自適應(yīng)節(jié)點(diǎn)生成物理信息網(wǎng)絡(luò)（Residual-basedAdaptive-nodePhysics-informedNeural Network，RAN-PINN）。物理信息網(wǎng)絡(luò)計(jì)算程函方程的優(yōu)勢(shì)在于訓(xùn)練完成之后可以反復(fù)利用網(wǎng)絡(luò)獲得旅行時(shí)信息。RAN-PINN在物理約束條件下訓(xùn)練網(wǎng)絡(luò)，調(diào)整網(wǎng)絡(luò)的權(quán)重和偏差，使解在給定訓(xùn)練點(diǎn)集上的損失函數(shù)最小化，能夠減小較大的殘差，改善物理信息網(wǎng)絡(luò)的旅行時(shí)計(jì)算效果。圖3為針對(duì)因式分解程函方程設(shè)計(jì)的物理信息網(wǎng)絡(luò)旅行時(shí)計(jì)算的原理示意圖。

圖3RAN-PINN旅行時(shí)計(jì)算原理示意圖

RAN-PINN首先通過(guò)自適應(yīng)節(jié)點(diǎn)方法初始化訓(xùn)練節(jié)點(diǎn)，然后將訓(xùn)練節(jié)點(diǎn)輸入全連接網(wǎng)絡(luò)，利用神經(jīng)網(wǎng)絡(luò)作為函數(shù)逼近器的能力求解程函方程，根據(jù)PDE殘差分布自適應(yīng)調(diào)整訓(xùn)練節(jié)點(diǎn)，對(duì)殘差較大的區(qū)域進(jìn)行重點(diǎn)訓(xùn)練。定義的損失函數(shù)是在選定的一組訓(xùn)練節(jié)點(diǎn)上最小化PDE殘差。式（3）的約束使得旅行時(shí)滿足邊界條件，損失函數(shù)只需要考慮程函方程。本文采用的損失函數(shù)為

式中： x_rc 是選定的配置節(jié)點(diǎn)； N_rc 是選定的配置節(jié)點(diǎn)數(shù)； q=2 時(shí)采用 L₂ 范數(shù)， q=1 時(shí)采用 L₁ 范數(shù)。節(jié)點(diǎn)生成模塊根據(jù)損失函數(shù)在每一輪訓(xùn)練開(kāi)始時(shí)生成新的計(jì)算節(jié)點(diǎn)，直到平均殘差小于給定的閾值 ε ，完成訓(xùn)練并輸出旅行時(shí) τ（x_rc） 。

2 模型測(cè)試

本文使用兩個(gè)模型驗(yàn)證程函方程依托物理信息網(wǎng)絡(luò)計(jì)算旅行時(shí)的可行性及精度。

2.1垂向速度梯度模型

首先驗(yàn)證本文方法的有效性。圖4a為一個(gè)具有恒定垂向速度梯度的模型，尺寸為 2km×2km ，垂向速度梯度為 0.5s^-1 ，橫向速度變化為0，震源位于中心處（圖4a黑星所示）。該模型的旅行時(shí)場(chǎng)解析解為

式中 g_v 和 g_h 分別表示垂直和水平方向的速度梯度。圖4b為旅行時(shí)解析解，圖4c為傳統(tǒng)FMM方法計(jì)算的旅行時(shí)，其中網(wǎng)格間距為 5m 。由圖4可以看出，對(duì)于恒定速度梯度模型，傳統(tǒng)的FMM方法能得出較高精度的旅行時(shí)解，因此本文將選用FMM方法求得的旅行時(shí)解與PINN解對(duì)比。圖5a和圖5c分別展示了隨機(jī)采樣、Hammersley采樣生成的節(jié)點(diǎn)分布，圖5b和圖5d分別為基于兩種固定節(jié)點(diǎn)分布的物理信息網(wǎng)絡(luò)解與FMM解的旅行時(shí)對(duì)比，兩種方法的結(jié)果均能夠有效地計(jì)算旅行時(shí)分布，盡管PINN的訓(xùn)練過(guò)程中只使用了2500個(gè)訓(xùn)練點(diǎn)，但能夠有效地獲取旅行時(shí)解。

圖4垂向梯度速度模型（a）的解析解（b）與傳統(tǒng)FMM方法計(jì)算結(jié)果（c）的對(duì)比

本文模型測(cè)試選取的全連接網(wǎng)絡(luò)包含4個(gè)隱藏層，每層70個(gè)神經(jīng)元，針對(duì)固定節(jié)點(diǎn)分布的情況，一旦對(duì)網(wǎng)絡(luò)進(jìn)行了訓(xùn)練，就可以在一個(gè)規(guī)則的網(wǎng)格上評(píng)估網(wǎng)絡(luò)。圖6a為RAN-PINN第1次訓(xùn)練后的PDE殘差分布、節(jié)點(diǎn)分布以及計(jì)算的旅行時(shí)，圖6b為第100O次訓(xùn)練后的PDE殘差分布、節(jié)點(diǎn)分布以及計(jì)算的旅行時(shí)。隨著RAN-PINN的訓(xùn)練，神經(jīng)網(wǎng)絡(luò)會(huì)根據(jù)式（7）計(jì)算的PDE殘差來(lái)調(diào)整節(jié)點(diǎn)分布，新生成節(jié)點(diǎn)逐漸聚集到誤差較大的區(qū)域；并且隨著訓(xùn)練次數(shù)的增加RAN-PINN的解會(huì)逐漸擬合解析解，同時(shí)二者之間的誤差也會(huì)減小。新生成的節(jié)點(diǎn)局部密度能夠根據(jù)當(dāng)前PDE殘差自適應(yīng)地變化，自適應(yīng)采樣方法收斂速度快，細(xì)化局部區(qū)域?qū)p少局部區(qū)域PDE殘差，一旦局部殘差減小則局部區(qū)域的節(jié)點(diǎn)密度也會(huì)自適應(yīng)地減小；如果局部PDE殘差在幾次迭代中再次增加，將在該區(qū)域中再次采樣更多節(jié)點(diǎn)。

本文用相對(duì)平均絕對(duì)誤差（RMAE）評(píng)價(jià)旅行時(shí)計(jì)算精度

式中： N_rc 是選定的配置節(jié)點(diǎn)總數(shù)； τ_ref（x_rc，i） 為第 i 個(gè)節(jié)點(diǎn)旅行時(shí)解析解或者采用FMM計(jì)算的數(shù)值解；τ（x_rc，i） 為第 i 個(gè)節(jié)點(diǎn)采用RAN-PINN計(jì)算的旅行時(shí)。

PINN采用隨機(jī)采樣、Hammersley采樣和本文的自適應(yīng)生成節(jié)點(diǎn)在網(wǎng)絡(luò)訓(xùn)練5000次后的輸出旅行時(shí)的RMAE分別為 0.2506%0.2505% 和0.2494% ，可見(jiàn)，采用基于殘差分布的自適應(yīng)節(jié)點(diǎn)生成物理信息網(wǎng)絡(luò)誤差更小。

圖7a給出了不同節(jié)點(diǎn)分布情況下的損失函數(shù)隨迭代次數(shù)的變化曲線。可見(jiàn)對(duì)于簡(jiǎn)單模型，三種節(jié)點(diǎn)分布都能獲得較好的訓(xùn)練結(jié)果。圖7b給出了自適應(yīng)節(jié)點(diǎn)方法在迭代過(guò)程中使用的節(jié)點(diǎn)數(shù)變化，在PINN訓(xùn)練過(guò)程中，可以于任意階段設(shè)置一組新的節(jié)點(diǎn)，而不是在訓(xùn)練期間使用相同的節(jié)點(diǎn)。本文每500次迭代根據(jù)殘差分布調(diào)整一次節(jié)點(diǎn)分布（如圖7b所示），模型測(cè)試時(shí)使用的GPU型號(hào)為3060ti，迭代5000次所耗時(shí)間為 2min16.8s ，訓(xùn)練所需時(shí)間與FMM、FSM求解程函方程所需時(shí)間相比還是可以接受的，并且隨著GPU算力的提高訓(xùn)練網(wǎng)絡(luò)的時(shí)間還會(huì)有進(jìn)一步的縮短。

圖5垂向梯度速度模型兩種固定節(jié)點(diǎn)分布物理信息網(wǎng)絡(luò)解的對(duì)比

圖6垂向梯度速度模型RAN-PINN第1（a）和第1000（b）次迭代后的PDE殘差分布（左）、節(jié)點(diǎn)分布（中）以及旅行時(shí)分布（右）

7垂向梯度速度模型三種節(jié)點(diǎn)生成方法的網(wǎng)絡(luò)損失函數(shù)變化曲線（a）及RAN-PINN的節(jié)點(diǎn)數(shù)隨迭代過(guò)程變化曲

2.2Marmousi模型

當(dāng)模型結(jié)構(gòu)復(fù)雜時(shí)，PINN的計(jì)算精度會(huì)有不同程度的下降。使用Marmousi模型（圖8a）測(cè)試RAN-PINN的計(jì)算精度，并與傳統(tǒng)FMM的結(jié)果進(jìn)行對(duì)比。模型尺寸為 3.5km×3.5km ，震源位于模型中心處。圖8b為采用FMM計(jì)算得到的旅行時(shí)。在PINN求解不同PDE的應(yīng)用中，損失函數(shù)通常采用 L₂ 范數(shù)（即 q=2， ，但 L₂ 范數(shù)適用于正態(tài)高斯分布的誤差。對(duì)于任意非均勻速度模型，誤差分布可能不是高斯的。在回歸問(wèn)題中，假設(shè)誤差具有隨機(jī)拉普拉斯分布， L₁ 范數(shù)（平均絕對(duì)誤差，即 q=1， 0可以提供對(duì)異常值魯棒的解。針對(duì)非均勻介質(zhì)本文采用 L₁ 范數(shù)自標(biāo)函數(shù)作為損失函數(shù)。

圖9給出了采用兩種固定節(jié)點(diǎn)方法生成的節(jié)點(diǎn)分布、PINN計(jì)算的旅行時(shí)以及PDE殘差分布。由圖可知，在模型復(fù)雜時(shí)，由于PINN存在譜偏差，會(huì)導(dǎo)致得到的PDE殘差相對(duì)于圖6b中平滑模型的PDE殘差大。比較FMM和PINN計(jì)算的旅行時(shí)等值線分布可知PINN計(jì)算的旅行時(shí)精度較高。

為了分析訓(xùn)練節(jié)點(diǎn)數(shù)對(duì)網(wǎng)絡(luò)性能的影響，使用不同節(jié)點(diǎn)數(shù)對(duì)PINN進(jìn)行訓(xùn)練，計(jì)算旅行時(shí)的RMAE如表1所示。由表1可知適當(dāng)?shù)脑黾庸?jié)點(diǎn)數(shù)量進(jìn)行PINN訓(xùn)練能夠減小不同節(jié)點(diǎn)分布方法的RMAE。針對(duì)復(fù)雜模型的物理信息網(wǎng)絡(luò)旅行時(shí)計(jì)算在構(gòu)造數(shù)據(jù)集時(shí)，需要使用足夠數(shù)量的訓(xùn)練樣本。對(duì)于平滑變速模型的測(cè)試，選擇總網(wǎng)格點(diǎn)可相對(duì)較少；而對(duì)于復(fù)雜模型，則需要選擇相對(duì)較多的網(wǎng)格點(diǎn)。

圖8Marmousi速度模型（a）及傳統(tǒng)FMM計(jì)算的旅行時(shí)（b）

圖10為Marmousi模型的自適應(yīng)節(jié)點(diǎn)迭代調(diào)整過(guò)程，隨著網(wǎng)絡(luò)的不斷訓(xùn)練，模型中各點(diǎn)的PDE殘差值在不斷減小，模型迭代后生成的節(jié)點(diǎn)更集中地分布在模型中速度分界面處。RAN-PINN能夠?qū)DE殘差較大的區(qū)域進(jìn)行更密集的采樣。本文模型測(cè)試采用自適應(yīng)節(jié)點(diǎn)生成方法的RMAE為0.7102% ，相較于兩種固定節(jié)點(diǎn)生成方法（見(jiàn)表1），其計(jì)算精度有所提升。當(dāng)然，與傳統(tǒng)的數(shù)值求解方法相比，該方法的實(shí)際計(jì)算優(yōu)勢(shì)取決于許多因素，包括網(wǎng)絡(luò)結(jié)構(gòu)、優(yōu)化超參數(shù)和采樣技術(shù)等。

本文采用的RAN-PINN在訓(xùn)練過(guò)程中的節(jié)點(diǎn)數(shù)目可以根據(jù)殘差的分布自適應(yīng)調(diào)整，無(wú)需人工參與，且實(shí)際采用的節(jié)點(diǎn)數(shù)目相對(duì)較少，能夠加快訓(xùn)練速度。在Marmousi模型測(cè)試中，采用隨機(jī)和Hammersley這兩種固定節(jié)點(diǎn)生成方法的PINN訓(xùn)練時(shí)間分別是190s和186s，而采用自適應(yīng)節(jié)點(diǎn)生成方法的RAN-PINN的網(wǎng)絡(luò)訓(xùn)練時(shí)間是 136s 。圖1la為不同節(jié)點(diǎn)分布情況下的損失函數(shù)曲線對(duì)比，可見(jiàn)，基于殘差分布的自適應(yīng)節(jié)點(diǎn)生成物理信息網(wǎng)絡(luò)收斂速度相對(duì)較快。圖11b為自適應(yīng)節(jié)點(diǎn)方法在迭代過(guò)程中使用的節(jié)點(diǎn)數(shù)變化曲線，開(kāi)始使用了較多的節(jié)點(diǎn)，后續(xù)節(jié)點(diǎn)數(shù)目可以適當(dāng)減少，以便節(jié)省訓(xùn)練時(shí)間。

圖9Marmousi模型隨機(jī)（a）和Hammersley（b）方法的節(jié)點(diǎn)分布（左）、PINN計(jì)算旅行時(shí)（中）及PDE殘差分布（右）

圖10Marmousi模型本文自適應(yīng)方法在第1（a）、第500（b）、第5000（c）次迭代后的PDE殘差分布（左）生成的節(jié)點(diǎn)分布（中）以及RAN-PINN計(jì)算的旅行時(shí)（右）的對(duì)比

表1兩種固定節(jié)點(diǎn)生成法節(jié)點(diǎn)數(shù)對(duì)PINN計(jì)算旅行時(shí)的影響分析

不同網(wǎng)絡(luò)架構(gòu)的旅行時(shí)計(jì)算精度不同，通過(guò)改變隱藏層的數(shù)量和每個(gè)隱藏層中神經(jīng)元的數(shù)量可以實(shí)現(xiàn)不同精度的旅行時(shí)計(jì)算。圖12為不同網(wǎng)絡(luò)架構(gòu)的旅行時(shí)計(jì)算結(jié)果對(duì)比，表2給出了物理信息網(wǎng)絡(luò)預(yù)測(cè)解和FMM解之間的RMAE。可見(jiàn)，隨著層數(shù)和神經(jīng)元數(shù)量的增加，網(wǎng)絡(luò)能夠近似更復(fù)雜的函數(shù)，計(jì)算精度也隨之提高。

圖11Marmousi模型三種節(jié)點(diǎn)生成方法網(wǎng)絡(luò)損失函數(shù)變化曲線（a）及RAN-PINN的節(jié)點(diǎn)數(shù)隨迭代過(guò)程變化曲線（b）

圖12Marmousi模型不同網(wǎng)絡(luò)架構(gòu)的計(jì)算旅行時(shí)對(duì)比

表2不同網(wǎng)絡(luò)參數(shù)的RMAE對(duì)比

3 PINN全局特征分析

大多數(shù)基于有限差分的程函方程計(jì)算方法都受到所謂的一階收斂精度的影響，需要較小的網(wǎng)格間距以提高精度和穩(wěn)定性。此外，通常需要從旅行時(shí)值及速度信息已知的點(diǎn)開(kāi)始計(jì)算，并逐漸擴(kuò)展波前以覆蓋整個(gè)區(qū)域，這意味著需要整個(gè)域的所有速度信息。而PINN具有全局特征解，只需要對(duì)計(jì)算域中總網(wǎng)格點(diǎn)的一部分進(jìn)行訓(xùn)練就可以可靠地預(yù)測(cè)旅行時(shí)。

3.1有限區(qū)域點(diǎn)模型

PINN具有外推插值能力，能夠在非間隙區(qū)域內(nèi)精確地求解旅行時(shí)。使用斜向梯度速度模型（圖13a）分析神經(jīng)網(wǎng)絡(luò)能否準(zhǔn)確計(jì)算間隙中的旅行時(shí)。模型尺寸為 2000m×2000m ，垂向梯度為 0.5s^-1 橫向梯度為 0.4s^-1 ，震源位于模型中心。設(shè)置兩種含有速度間隙的模型（該區(qū)域不選擇訓(xùn)練節(jié)點(diǎn)）如圖13b和圖13c所示。圖14給出了不同情況下的節(jié)點(diǎn)分布以及旅行時(shí)計(jì)算結(jié)果。由圖可知，PINN可以使用具有間隙的速度模型計(jì)算旅行時(shí)場(chǎng)，具有克服速度信息缺失的能力。當(dāng)然，旅行時(shí)計(jì)算的不確定性會(huì)隨著間隙的增大而增加。間隙模型1的旅行時(shí)RMAE為 0.2460% ，間隙模型2的RMAE為0.2477% 。

圖13斜向梯度速度模型（a）及其含有間隙模型1（b）和模型2（c）

圖14斜向梯度速度模型RAN-PINN的節(jié)點(diǎn)分布（a）及其計(jì)算旅行時(shí)與解析解和FMM旅行時(shí)的對(duì)比（b）左：完整模型；中：含間隙的模型1；右：含間隙的模型2

3.2起伏地表模型

在陸上地震勘測(cè)中地表通常不平坦。圖15為一個(gè)起伏地表速度模型，震源位于黑星處。傳統(tǒng)求解方法是將笛卡爾坐標(biāo)系轉(zhuǎn)換為曲線坐標(biāo)系，以數(shù)學(xué)方式展平自由表面。但是這種方法增加了額外的計(jì)算成本，而且當(dāng)?shù)匦渭眲∽兓瘯r(shí)可能導(dǎo)致計(jì)算不穩(wěn)定。圖16為本文方法第1次和第500次迭代的節(jié)點(diǎn)分布和計(jì)算旅行時(shí)分布，通過(guò)逐步迭代，起伏地表附近的節(jié)點(diǎn)分布逐漸變密，能夠更好地適應(yīng)起伏地表情況，與第1次迭代相比，第500次迭代的旅行時(shí)分布變得更加準(zhǔn)確。與傳統(tǒng)方法相比，基于PINN的方法作為一種無(wú)網(wǎng)格方法，在處理起伏地表模型旅行時(shí)計(jì)算問(wèn)題時(shí)，PINN不需要任何特殊處理，僅使用自由表面下方的節(jié)點(diǎn)訓(xùn)練和評(píng)估網(wǎng)絡(luò)即可。

圖15起伏地表模型

圖16起伏地表模型本文方法第1（a）、第500（b）次迭代時(shí)生成的節(jié)點(diǎn)分布（左）、計(jì)算的旅行時(shí)分布（右）

4結(jié)論

不依賴地形，在起伏地表模型中不需要調(diào)整坐標(biāo)系即可獲得穩(wěn)定、準(zhǔn)確的旅行時(shí)場(chǎng)。

本文提出了一種基于殘差分布的自適應(yīng)節(jié)點(diǎn)調(diào)整、使用物理信息網(wǎng)絡(luò)（RAN-PINN）求解各向同性介質(zhì)程函方程的旅行時(shí)計(jì)算方法，對(duì)于復(fù)雜速度模型能獲得穩(wěn)定和準(zhǔn)確的結(jié)果，有效減少誤差，得出的主要結(jié)論如下。

（1）基于神經(jīng)網(wǎng)絡(luò)求解程函方程能夠通過(guò)將PDE嵌入損失函數(shù)來(lái)預(yù)測(cè)旅行時(shí)，利用神經(jīng)網(wǎng)絡(luò)作為通用函數(shù)逼近器的能力，最小化在選定節(jié)點(diǎn)上的殘差分布。與純數(shù)據(jù)驅(qū)動(dòng)的學(xué)習(xí)神經(jīng)網(wǎng)絡(luò)相比，PINN在訓(xùn)練過(guò)程中施加了物理信息約束，得到的結(jié)果能夠滿足偏微分方程中包含的物理約束。

（2）基于殘差的自適應(yīng)節(jié)點(diǎn)生成方法通過(guò)引入變密度節(jié)點(diǎn)生成方法，能夠改進(jìn)訓(xùn)練點(diǎn)采樣和分布并提高旅行時(shí)的計(jì)算精度。

（3）由于PINN具有全局特征解，因此本文方法

參考文獻(xiàn)

[1]LAWTOND C.Computation of refraction static cor rections using first-break traveltime differences[J]. Geophysics，1989，54（10）：1289-1296.

[2] GUOR，LIMK，YANGF，etal.Firstarrival traveltime tomography using supervised descent learning technique[J].Inverse Problems. 2019，35（10）： 105008.

[3] 郭振波，孫鵬遠(yuǎn)，李培明，等．射線追蹤方程與程函方 程的初至旅行時(shí)層析對(duì)比[J].石油地球物理勘探， 2019，54（3）：558-564，576. GUO Zhenbo，SUN Pengyuan，LI Peiming，etal. Raytracing equation-based first-arrival traveltime tomographyand eikonal equation-based first-arrival traveltime tomography[J]. Oil Geophysical Prospecting，2019， 54（3）：558-564，576.

[4] SETHIAN JA. Fast marching methods[J]. SIAM Review，1999，41（2）：199-235

[5] ZHAO HK.A fast sweeping method for eikonal equations[J]. Mathematics of Computation，20o5，74（250）： 603-627.

[6] 謝春，劉玉柱，董良國(guó)，等．伴隨狀態(tài)法初至波走時(shí)層 析[J].石油地球物理勘探，2014，49（5）：877-883. XIEChun，LIU Yuzhu，DONG Liangguo，et al. First arrival traveltime tomography based on the adjoint state method[J]. Oil Geophysical Prospecting，2014， 49 （5） ：877-883.

[7］NAJAFABADI MM， VILLANUSTREF， KHOSHGOFTAAR T M，et al. Deep learning applications and challenges in big data analytics[J]. Journal of Big Data，2015，2（1） ：1-21.

[8] ANITESCUC，ATROSHCHENKOE，ALAJLAN N，et al. Artificial neural network methods for the solution of second order boundary value problems[J]. Computers，Materials and Continua，2019，59（1） ：345-359.

[9] TAUFIKM H，WAHEEDUB，ALKHALIFAHT A.Upwind，n more：Flexible traveltime solutions usingphysics-informed neural networks[J]. IEEE Transactions on Geoscience and Remote Sensing， 2022，60（1） ：1-12.

[10] RAISSI M，PERDIKARISP. Physics-informed neuralnetworks：adeeplearningframeworkforsolvingfor wardandinverse problems involvingnonlinearpartial differential equations[J]. Journal of Computational Physics，2019，378：686-707.

[11] SMITH JD，AZZIZZADNESHLI K，ROSS Z. Eikonet：Solving the eikonal equation with deep neural networks[J]. IEEE Transactions on Geoscience and Remote Sensing，2020，99（12） ：1-12.

[12] WAHEEDUB，HAGHIGHATE，ALKHALIFAH T，et al． PINNeik：Eikonal solution using physicsinformed neural networks[J]. Computers amp; Geosciences，2021，155：104833.

[13] GRUBAS S，ANTON D，Loginov G.Neural eikonal solver：improving accuracy of physics-informed neural networks for solving eikonal equation in case of caustics[J].Journal of Computational Physics，2O23，474： 111789.

[14] LEEH，KANG IS. Neural algorithm for solving differential equations[J].JournalofComputational Physics，1990，91（1）：110-131.

[15] 都國(guó)寧，譚軍，宋鵬，等．基于物理信息驅(qū)動(dòng)神經(jīng)網(wǎng)絡(luò) 的三維初至波旅行時(shí)計(jì)算方法[J].石油地球物理勘 探，2023，58（1）：9-20. DU Guoning，TAN Jun，SONG Peng，et al. 3D traveltime calculation of first arival wave using physics-informed neural network[J]. Oil Geophysical Prospecting，2023，58（1） ：9-20.

[16] YU J，LUL，MENGXH，et al.Gradient-enhanced physics informed neural networks for forward and inverse PDE problems[J]. Computer Methods in Applied Mechanicsand Engineering，2022，393：114823.

[17] NABIAN MA，GLADSTONERJ，MEIDANIH. Efficient training of physics-informed neural networks viaimportance sampling[J]. Computer-Aided Civil and Infrastructure Engineering，2021，36（8） ：962-977.

[18]DAS S，TESFAMARIAM S. State-of-the-art review of design of experiments for physics -informed deep learning [EB/OL]. （2022-02-16）[2024-09-17]，https： //arxiv.org/abs/2202.06416.

[19] PENGW，ZHOU W，ZHANGX，et al. RANG：A residual-based adaptive node generation method for physics-informed neural networks[EB/OL]. （2022-05- 10）[2024-09-17]，https：//arxiv.org/abs/2205.01051.

[20] QIAN J，SYMES W. An adaptive finite-difference method for traveltimes and amplitudes[J].Geophysics， 2002，67（1）：167-176.

[21] FOMEL S，LUO S T，Zhao H K. Fast sweeping method for the factored eikonal equation[J]. Journal of Computational Physics，2009，228（17） ：6440-6455.

[22] 張?jiān)疲钕＃壮ⅲ龋?D/3D起伏地表多震相地 震波走時(shí)的因式分解程函方程算法[J].石油地球物 理勘探，2023，58（4）：857-871. ZHANGYun，LIXihai，BAIChaoying，etal.Multiphase seismic traveltime computation in 2D/3D undulatedsurface model using factored eikonal equation[J]. Oil GeophysicalProspecting，2023，58（4） ：857-871.

[23] 崔寧城，黃光南，李紅星，等．求解橢圓各向異性介質(zhì) 程函方程的源點(diǎn)快速掃描算法[J].石油地球物理勘 探，2019，54（4）：796-804. CUI Ningcheng，HUANG Guangnan，LI Hongxing，et al.A fast scanning algorithm of source points for solvingelliptic anisotropic dielectric equation[J]. Oil Geophysical Prospecting，2019，54（4） ：796-804

[24] 潘奕銘，桑運(yùn)云，姚雪峰，等．一種多信息約束的初至 波走時(shí)層析反演優(yōu)化方法[J].石油物探，2023，62（2）： 271-280. PANYiming，SANG Yunyun，YAO Xuefeng，et al. Optimal first-arrival travel-time tomographywith multiinformation constraintsforcomplex near-surface scheme[J]. Geophysical Prospecting for Petroleum， 2023，62（2）：271-280.

[25]TAUFIKM H，WAHEEDUB，ALKHALIFAHT A．Aneural network based global traveltime function （GlobeNN）[J]. Scientific Reports，2023，13（1）：7179.

[26]FORNBERG B，F(xiàn)LYER N. Fast generation of 2-D node distributions for mesh-free PDE discretizations [J].Computers amp;. Mathematics with Applications， 2015，69（7） ：531-544.

（本文編輯：宜明理）

作者簡(jiǎn)介

唐杰副教授，1980年生；2003年獲中國(guó)科技大學(xué)地球物理專業(yè)學(xué)士學(xué)位；2008年獲中國(guó)科技大學(xué)地球物理專業(yè)博士學(xué)位；現(xiàn)在中國(guó)石油大學(xué)（華東）地球科學(xué)與技術(shù)學(xué)院主要從事深度學(xué)習(xí)與地球物理反演等方面的研究。