基于波動方程的地震波數值模擬研究綜述

摘要：

地震波場數值模擬在地震勘探、地震資料處理和地球構造研究等方面發揮著重要的作用。波動方程數值模擬方法充分考慮了地震波傳播的動力學特征和幾何學特征，可以為地震波傳播機理的研究和復雜地層的解釋提供強有力的理論支持，是目前應用較為廣泛的地震波場數值模擬方法之一。本文調研了五種基于波動方程的數值模擬方法：有限差分法易于理解，但數值頻散問題明顯；偽譜法精度高，但計算效率低；有限元法適用于復雜模型，但計算資源消耗大；譜元法適合高精度問題，但對計算內存需求較高；基于物理信息神經網絡的深度學習法具有較強的適應性，但訓練成本較高。并分別敘述了這五種數值模擬方法的理論基礎、適用條件和最新進展。未來，地震波場數值模擬方法應結合深度學習等最新技術，優化邊界條件模擬真實的邊界反射情況，提高模擬的精度和效率。

關鍵詞：

波場模擬；有限差分法；偽譜法；有限元法；譜元法；物理信息神經網絡

doi：10.13278/j.cnki.jjuese.20230308

中圖分類號：P631.4

文獻標志碼：A

A Comprehensive Review of Numerical Simulation of Seismic Waves Based on "Wave Equation

Li Hang¹， Sun Yuhang^{1， 2}， Li Jiahui^{1， 2}， Li Xuegui^{1， 2}， Dong Hongli^{1， 2}

1. Key Laboratory of Networked Control and Intelligent Systems， Northeast Petroleum University， Daqing 072751， "Heilongjiang， China

2. "NEPU Sanya Offshore Oil and Gas Research Institute， "Sanya 572024， Hainan， China

Abstract：

Numerical simulation of seismic wave fields is crucial for seismic exploration， seismic data processing， and the study of Earth’s structures. The wave equation numerical simulation method takes into account the dynamic and geometric characteristics of seismic wave propagation， providing a solid theoretical basis for studying the mechanism of seismic wave propagation and interpreting complex geological structures. It is currently one of the most widely used methods for simulating seismic wave fields. This article surveys five wave equation-based numerical simulation methods： The finite difference method is easy to understand， but suffers from numerical dispersion issues; The pseudo spectral method offers high accuracy but low efficiency; The finite element method is suitable for complex models but requires substantial computational resources; The spectral element method is appropriate for high-precision problems but demands significant computational memory; And the deep learning method based on physics-informed neural networks demonstrates strong adaptability， though it comes with high training costs. The theoretical foundations， applicable conditions， and latest advancements of these five numerical simulation methods are described respectively. In the future， seismic wave field numerical simulation methods should integrate cutting-edge technologies such as deep learning， optimize boundary conditions to simulate actual boundary reflections， and enhance the precision and efficiency of simulations.

Key words：

wavefield simulation; finite difference method; pseudo spectral method; finite element method; spectral element method; physics-informed neural network

0"引言

隨著地震勘探的不斷深入，目標儲藏的情況越來越復雜，勘探難度和勘探需求的精度越來越高。地震波包含了地下地質結構特征的大部分信息，地震波場正演模擬能夠驗證地震資料解釋結果的合理性，在地震勘探中的地震反演等環節發揮著重要的作用^[1]。地震波場正演模擬可以分為物理模擬和數值模擬兩種方法。相對而言，后者的計算成本更低、計算效率更高、可重復性更強，受到了學者們的廣泛關注。地震波場數值模擬方法主要包括三類：幾何射線法、積分方程法和波動方程法^[2]。

幾何射線法主要利用地震波運動學信息，以射線理論為核心，以計算機為載體實現數值模擬，得到介質中描述地震波運動特征的關鍵物理參數^[3]。幾何射線法適用范圍較廣，但由于很少利用必要的動力學信息，這類方法在處理一些復雜場景（如非均勻介質）時精度不高。積分方程法是以惠更斯原理為理論基礎的一種求解非線性問題的方法^[4]，通常用來模擬具有特定邊界的地質體產生的地震波。利用這類方法進行模擬時常常需要求解大型矩陣，因此對計算時間和計算內存的需求較高。波動方程法的本質是通過求解地震波動方程來進行數值模擬，該方法可以提供豐富的波場動力學和運動學信息。在同等條件下，波動方程法數值模擬結果的精度要高于上述其他兩種方法^[5]。隨著計算機算力的提高和算法的優化，這類方法的效率有所提高，被廣泛使用。

波動方程的數值模擬涉及多個關鍵步驟：1）構建適用的波動方程，需要針對不同的介質模型進行調整，確保方程能夠準確地描述介質中波的傳播特性；2）將求解域離散為規則或非規則的網格系統，為后續的數值計算打下基礎；3）設定恰當的邊界條件，以模擬實際情況并防止邊界反射引起誤差，常用的方法為添加人工吸收邊界條件或完全匹配層（perfectly matched layer， PML）邊界；4）確定準確的初始條件，以確保數值模型能夠反映真實的物理現象；5）采用合適的數值方法求解離散化的波動方程。

目前，基于波動方程的傳統地震波正演數值模擬方法主要包括有限差分法（finite difference method， FDM）、偽譜法（pseudo spetral method， PSM）、有限元法（finite element method， FEM）和譜元法（spectral element method， SEM）。有限差分法是地震波場數值模擬的主要方法之一，普遍存在頻散問題且需要對邊界條件進行合適的處理^[6]。Cheng等^[7]基于旋轉交錯網格，引入合適的高精度單邊仿真型有限差分算子，解決起伏地表的自由邊界條件問題，顯著提高了有限差分法的計算效率。偽譜法利用傅里葉級數展開式表征波場函數，正演結果在頻率域具有較高的精度和分辨率。劉炯等^[8]在復雜起伏地表下地震波傳播過程的研究中，引入貼體網格技術，首先對地表不規則的模型進行網格劃分，然后使用偽譜法進行數值模擬，模擬結果具有較高的精度。有限元法的網格劃分形式相對靈活，能夠適用于異常地質體界面和不規則的起伏地表。黃旭東等^[9]使用擴展有限元法模擬了地下巖體的水力壓裂過程，獲得了良好的數值模擬結果。譜元法將有限元法與偽譜法相結合，同時具有兩者的優點，速度快、精度高、適應性強，被廣泛應用在地震波數值模擬中^[10]。

隨著研究的深入，學者們嘗試結合物理信息神經網絡（physics-informed neural network， PINN）解決地震波場波動方程數值模擬問題并取得了一些有益的成果。波動方程的求解可以看作為偏微分方程的求解，神經網絡因其強大的特征提取能力和逼近能力在解決這類問題時表現優異。以PINN為主的深度學習算法及其衍生算法在不同類型的偏微分方程求解及地震波數值模擬中均發揮了較好的作用，并取得了較好的結果^[11]。

總體來看，國內外對地震波數值模擬的研究成果越來越多，也逐漸出現一些綜述性文章對該類方法進行總結^{[1， 12-13]}。隨著研究的不斷深入，如何提高地震波數值模擬的效率和模擬數值結果的精度，減少數值模擬的計算耗時等已逐步成為研究熱點。本文對近年來該領域的相關研究內容做了系統性的總結和歸納，方便讀者對該領域進行較為深入的了解。

1"有限差分法

有限差分法表達形式簡單，概念易于理解，通常可用于解決復雜問題，是目前應用范圍較為廣泛、研究較為深入的地震波數值模擬方法之一^[14]。

二維常規波動方程表達式為

²ut²=c²（²ux²+²uy²）。（1）

式中：u（x，y，t）為在二維空間中傳播的波動函數；x和y為空間坐標；c為波速；t為時間。有限差分法主要利用運用泰勒公式的展開形式進行運算，表達式為

u（x+h）=u（x）+hux+h²2！²ux²+h³3！³ux³+…（2）

式中，h為步長。以中心差分為例，將時間和空間的展開式代入其導數的差分近似中，就可以得到一個更為精確的差分表達式來近似表達波動方程中的時間和空間偏微分算子：

²ut²=

u（x，y，t+Δt）-2u（x，y，t）+u（x，y，t-Δt）Δt²;

²ux²=

u（x+Δx，y，t）-2u（x，y，t）+u（x-Δt，y，t）Δx²;

²uy²=

u（x，y+Δy，t）-2u（x，y，t）+u（x，y-Δy，t）Δy²。

（3）

將式（3）代入式（1）中轉換為差分離散波動方程并進行求解，即可實現波動方程數值模擬。

Alterman等^[15]首次將有限差分法應用于求解均勻介質的二階彈性波位移方程。在此基礎上，Boore^[16]將有限差分法應用于非均勻介質中的地震波波場模擬，提出了非均勻介質二階彈性波有限差分法。Madariaga^[17]采用交錯網格技術求解一階速度應力方程組。Dablain^[18]將高階有限差分算子應用于標量波場的數值模擬。這些研究進一步推動了有限差分法在地震波數值模擬領域的發展。

有限差分法雖然計算效率較高，但因利用離散形式代替連續形式進行運算，通常會產生數值頻散問題。低階有限差分在數值模擬中普遍存在較為嚴重的數值頻散情況，影響了數值模擬結果的精度。此外，有限差分法的網格剖分情況也會直接影響數值模擬結果的精度。常見的網格類型主要包括三種：規則網格、交錯網格和旋轉交錯網格。它們的不同之處在于網格排列方式和內部結構不同。規則網格主要用于處理常規介質；交錯網格主要用于高階差分地震波場的數值模擬中，具有更高的精度^[19]。為了有效控制頻散誤差，學者們研究了高階交錯網格有限差分的優化技術。唐超等^[20]提出利用交替方向乘子法進行一階應力-速度聲波方程交錯網格有限差分數值模擬，對交錯網格有限差分系數進行了較為精確的計算，并通過長時程數值實驗，驗證了基于L₁范數的優化方法在減小誤差積累方面具有一定的優越性。Zhou等^[21]提出了一種低成本的混合交錯網格有限差分形式，相比于傳統的交錯網格有限差分，利用該差分形式進行數值模擬時具有更好的精度和穩定性。相比于交錯網格，更為復雜的旋轉交錯網格能夠進一步提高數值模擬結果的精度^[22]。

除了常見的幾種網格類型外，混合網格、變網格及自適應網格等也是研究的重點。其中，自適應網格憑借其靈活的結構可以在低速區的不規則區間使用密網格，在其他區間使用疏網格，具有更高的計算效率。張志佳等^[23]提出一種起伏多重變加密網格有限差分波動方程數值模擬方法，能兼顧計算精度和效率，且能很好地適應實際復雜近地表介質。Zhang等^[24]為了高效模擬強異質介質中的聲波傳播，采用了自適應網格細化技術，通過模擬區域中的介質變化和波場傳播情況動態調整網格分辨率，在關鍵區域增加網格密度，更加準確地描述和模擬聲波在復雜介質中的傳播行為，顯著提高了計算速度。于明浩等^[25]采用不等距有限差分法網格技術，對起伏地表條件下的地震波場進行了數值模擬研究。

根據實現域的不同，有限差分法可以分為時間域有限差分^[26]和頻率域有限差分^[27]兩種。其中，頻率域有限差分法常應用于與頻率相關的介質模型中。相對而言，時間域有限差分法更為直觀、效率更高且應用更為廣泛。此外，根據求取差分系數所采用的頻散關系類型，有限差分法可被分為空間域有限差分和時空域有限差分兩種^[28-29]。在地震波動方程數值模擬中，時空域有限差分法因同時考慮時間和空間精度方面的優勢而備受青睞。為提高地震波動方程數值模擬精度，彭更新等^[30]采用了時間二階、空間2M階的交錯網格有限差分法，通過引入時空域的頻散關系和泰勒級數展開方法，提出了一種改進的時空域交錯網格有限差分法，在保持差分格式不變性的同時，提高了數值模擬精度。此外，他們對空間域與時空域交錯網格有限差分法數值模擬結果進行了比較，在3.0 s時壓力場波場快照與炮集記錄的模擬結果中：空間域交錯網格有限差分法數值模擬均呈現出明顯的數值頻散現象，表現出較低的模擬精度；相比之下，時空域交錯網格有限差分法的數值模擬結果在相同時間點沒有顯示明顯的數值頻散，精度較高（圖1）。根據時間和空間離散化的方式，有限差分法可以分為顯式和隱式兩種格式。顯式有限差分法使用當前時間步的信息來更新下一個時間步，計算速度較快，但可能會存在數值穩定性較差的問題；而隱式有限差分法利用當前和未來時間步的信息進行迭代，計算速度較慢但更加穩定。在實際應用中，通常會根據問題的特點和數值穩定性的要求選擇合適的方法。徐世剛等^[31]提出一種改進的時間高精度-空間隱式交錯網格差分法，并將其與幾種現有的差分方法進行了對比，彈性波數值模擬結果表明該方法能有效壓制數值頻散干擾，產生高精度的模擬結果。

廣義有限差分法是在傳統有限差分法基礎上的擴展和改進。傳統有限差分法將計算區域離散為規則的網格，并利用差分格式對空間和時間進行離散化，從而將偏微分方程轉化為代數方程組。然而，在復雜介質或非均勻網格的情況下，傳統有限差分法的精度和穩定性可能受到限制。而廣義有限差分法通過引入非均勻網格、高階差分格式和高維插值等技術，克服了傳統有限差分法的局限性。在全波形反演中，Li等^[32]采用了不規則節點分布的無網格廣義有限差分法建立聲彈耦合介質中的波動方程，在空間離散化過程中引入了更高階的差分算子，提高了數值模擬結果的精度和穩定性，更精確地模擬了波場的傳播過程。

緊湊型有限差分法相比于傳統有限差分法可在較少的離散點下實現相同或更高的精度，且可以減少數值耗散或數值擴散現象。Kosloff等^[33]提出了一種通用的緊湊型有限差分方案，該方案可以使用任意個點的函數值計算當前點的導數值，并用于二維聲波方程的數值解以及各向同性介質中二維動態彈性方程的求解。Li等^[34]提出了一種新的緊湊型有限差分方案求解三維聲波方程。Chen等^[35]提出優化的緊湊型有限差分方案和二階中心有限差分方案，分別對二維聲波方程的空間和時間導數進行近似，數值模擬結果表明，優化的緊湊型有限差分方案優于同階的中心有限差分方案。整體而言，緊湊型有限差分方案因其較高的空間效率、較好的邊界條件處理能力以及較高的精度，在求解勘探地震學中的

地震波動方程方面顯示出巨大的潛力。

目前，關于有限差分法在地震波場數值模擬中的研究主要集中于網格類型的選擇和差分系數的求取等方面。隨著計算機技術的進步，為增強數值模擬的穩定性并提高數值模擬結果的精度，學者們在先前的研究基礎上提出一系列優化策略。Guo等^[36]為了有效地獲得最優交錯網格差分階數和差分系數，結合人工神經網絡進行參數估計。韓復興等^[37]將二階CE（Clayton-Enquist）吸收邊界條件和PML衰減邊界條件進行組合，有效提高了數值模擬效率。趙明哲等^[38]利用GPU（graphics processing unit）并行算法提高了高階旋轉交錯網格有限差分法的計算效率，模型數據實驗結果證明了其有效性。

綜上所述，有限差分法具有表達形式簡單、概念易于理解的優點，在地震波動方程數值模擬中得到了廣泛應用。未來，研究者們可開發更多高精度的

有限差分算法，以減少數值頻散、數值耗散等常見問題。同時，自適應網格技術和多尺度模擬方法有望進一步發展，從而更準確地捕捉復雜的地質結構，模擬地震波傳播過程。此外，計算機技術的持續進步，尤其是在并行處理和算力等方面技術的突破，利用GPU加速等資源，將大幅提升計算效率，增強處理復雜科學問題的能力。

2"偽譜法

偽譜法首次提出于上世紀70年代^[39]，最初應用于大氣預測以及非線性波動模擬等領域，20世紀80年代被引入到地震波數值模擬的相關研究中。Kosloff等^[40]利用傅里葉偽譜法求解二維聲波方程，并證明該方法解決此問題時比有限差分法更有效。

對于二維波動方程來說，偽譜法首先將二維空間離散化為n_x×n_y個網格點，然后對每個網格點（x_i，y_i）使用插值多項式來逼近待求解函數u（x，y，t），最后利用數學領域中的傅里葉變換將導數求解過程從空間域轉換到波數域，表達式為

F（k_x，k_y，t）=1L_xL_y∫－∫－u（x，y，t）e^－ik_x^xe^－ik_y^ydxdy。（4）

式中：F（k_x，k_y，t）為u（x，y，t）的傅里葉變換；L_x和L_y分別為x方向和y方向的空間域尺寸；k_x和k_y分別為x方向和y方向的波數。傅里葉變換將導數的求解過程從空間域轉換到波數域，波數域中的導數運算可以轉換為相乘運算，可以相對簡單地求解方程。之后，將其結果通過傅里葉逆變換轉換到時間域即可完成二維波動方程的求解。這個傅里葉正變換再逆變換的過程能夠消除直接求導產生的頻散現象，提高數值模擬的效率和精度。

偽譜法可以看作為一種無窮階極限情況下的空間域有限差分法。在同等模擬精度下，偽譜法比有限差分法采用更大的網格間距和更少的采樣點數，需要更小的存儲空間^[41]。在網格間距和時間步長相同的情況下，孫獻果等^[42]對比了有限差分法和偽譜法不同時刻的波場快照，發現有限差分法會出現頻散現象，主要表現為主圓環內部存在較多不明顯的小圓環（圖2）。表1展示了有限差分法和偽譜法的綜合性能對比情況，當兩者正演結果的精度相似時，有限差分法雖然計算耗時小，但所需要的計算內存明顯高于偽譜法。

近年來，利用偽譜法進行地震波動方程數值模擬的研究日漸增多^[43-44]。通過網格劃分的方式可以

提高偽譜法數值模擬結果的精度。相比于常規網格偽譜法和高階有限差分法，譚文卓等^[45]采用梯形網格偽譜法進行地震波動方程數值模擬，有效地減少了計算時間，降低了數值頻散。李元燮等^[46]利用時間分裂法和旋轉交錯網格偽譜法相結合的數值模擬方法實現了波場的數值模擬，能夠得到穩定的、高精度的模擬結果。

在地震波場數值模擬中，由于計算資源的限制，無法模擬無限大的區域，因此必須在某個有限的計算區域上設置邊界。這個邊界是人為引入的，稱為人工截斷邊界。但人工截斷邊界會在計算區域邊界產生強邊界反射，很大程度上降低了地震波場模擬的精度。李青陽等^[47]基于PML邊界，采用高階偽譜法進行數值模擬，有效壓制了時間頻散，提高了計算效率。為了降低數值噪聲，提高模擬結果的精度，馬銳等^[48]首先利用分裂格式構建PML吸收邊界和海綿邊界，然后組合形成復合吸收邊界，最后基于此使用偽譜法進行數值模擬，獲得了較好的數值模擬結果。陳漢明等^[49]在對分數階黏滯聲波方程進行數值模擬時，采用近似卷積型PML作為吸收邊界條件的交錯網格偽譜法，得到了具有較高精度的模擬結果。隨著計算機算力的提高，基于GPU加速的偽譜法的計算耗時得到了大幅度減少^[50]。

綜上所述，在地震波動方程數值模擬中，偽譜法具有高精度和低內存需求的優勢。然而，該方法數值模擬

計算效率較低，且在空間分辨率和

邊界條件的處理上存在一定的限制。因此，在應用偽譜法進行數值模擬時，需要綜合考慮其優缺點，以獲得準確可靠的模擬結果。未來，偽譜法可以引入高階離散化方法進一步提高精度，并借鑒動態網格技術，智能地調整網格分辨率，通過改進邊界條件處理方法，增強其適用性；也可以開發對應的優化算法，提升數值模擬效率，并與深度學習相結合，優化模擬參數、提高地震波數值模擬結果的準確性和可靠性。

3"有限元法

有限元法在求解地震波動方程時，首先根據數值模擬要求選擇單元網格類型，將連續的求解域Ω離散化成有限數量的網格小單元，在每個網格小單元Ω_i內部定義節點，這些節點通常位于單元的頂點或適當的位置，以滿足數值精度和計算效率的要求。然后根據求解精度要求和單元類型選擇單元內的插值函數，如線性插值、二次插值等，用于在單元的節點之間插值，以得到單元內任意點上的場變量值。以二維波動方程為例，位移近似解u_h（x，y，t）可以表示為

u（x，y，t）=∑ni=1φ_i（x，y）ψ_i（t）。（5）

式中：φ_i（x，y）為插值函數；ψ_i（t）為在每個時間步長上位移的近似值；n的具體值取決于所選擇的插值函數和單元類型。接著將每個單元的局部位移解u_h組裝成全局位移解U，再將每個小單元上計算的剛度矩陣和質量矩陣組裝為全局剛度矩陣K和全局質量矩陣M。加載邊界條件向量F_b和震源項F_s（t），將整個系統的方程組裝成全局的代數方程組：

M+KU=Fb+Fs（t）。（6）

式中，為空間變量上位移的二階時間導數向量。最后，將這個方程組（式（6））在每個時間步長上求解，以得到求解域內各個時刻的位移場。

20世紀70年代，Lysmer等^[51]和Drake^[52]將有限元法應用于地震波場模擬的研究。在二維波動方程數值模擬中，Mullen等^[53]對采用不同時間離散方法的低階有限元方法進行了評估，研究其精度、穩定性和數值頻散性能。Chen^[54]則在各向異性介質下應用有限元法求解波動方程。此外，Serón等^[55]分析了有限元法在計算彈性波動方程方面的性能。總體而言，相較于有限差分法，在邊界處和奇點附近，有限元法能夠更好地處理精度損失和振蕩現象，這得益于有限元法采用逐單元的離散化方法^[56]。

有限元法在剖分非結構化網格、處理自然邊界條件方面表現出色，更適合模擬復雜地下介質中波的傳播；然而該方法計算量較大，且低階有限元法地震波動方程數值模擬結果精度不高。傳統有限元法在地震波場數值模擬中通常需要大量的計算內存，且效率不高^[57]。為提升有限元法的計算效率，蘇波等^[58]采用一種改進的核矩陣存儲策略，在網格剖分時采取了任意非規則單元的形式，能夠在保證數值模擬結果精度的同時，降低計算所需的內存，提高計算效率。他們在此基礎上發展了修正有限元法，該新方法在復雜模型的高精度波場數值模擬中表現出色，且在大尺度及復雜模型的地震波傳播數值模擬中具有一定優勢。在使用有限元法進行數值模擬時，選擇合適的邊界條件對于確保模擬結果的準確性和可靠性至關重要。Zhao等^[59]使用有限元方法進行數值模擬時，提出了一種非分裂式的完全匹配層（non-split perfectly matched layer， NPML）邊界條件，以吸收由人工邊界產生的反射，該研究主要針對擴散黏性波進行數值模擬，結果表明NPML邊界條件可以有效吸收人工邊界反射（圖3）。

間斷有限元法（discontinuous galerkin method， DGM）是以有限元法為基礎的一種數值模擬方法，該方法由Reed等^[60]在求解中子輸運方程時提出的。間斷有限元法在進行數值模擬時能夠有效壓制數值頻散，且網格剖分形式更加靈活，能夠減少計算量并提高計算精度。賀茜君等^[61]在間斷有限元法的基礎上提出一種新的加權Runge-Kutta間斷有限元法（Runge-Kutta DGM， RKDG），數值模擬結果顯示，該方法可有效模擬包含球體在內的非規則模型以及非均勻Marmousi模型中衰減的聲波波場。He等^[62]將加權的RKDG擴展到求解二維橫向各向同性介質中的彈性波動方程，并使用非結構化三角網格行模擬，證明該方法具有較好的靈活性。

現有的間斷有限元法通常對二維彈性波進行數值模擬，Hong等^[63]提出了一種用于三維彈性波模擬的間斷有限元法，新方法中采用的自適應四面體分裂法可以有效地模擬任意復雜起伏地表中彈性波的傳播過程，具有較高的模擬精度和計算效率。韓德超等^[64]通過構建不同的三角形單元來組成周期性網格進行數值模擬，并提供了基于三角形網格的間斷有限元法彈性波模擬的參數選擇指南，研究表明數值各向異性最小的是等邊三角形網格，而直角三角形和鈍角三角形網格則表現出較強的數值各向異性。間斷有限元法相比傳統的有限元法和有限差分法具有許多優點：RKDG能夠更準確地捕捉解的局部特征和不連續性，適用于處理具有強烈變化或不光滑解的問題；其次，允許使用更高階的多項式來逼近解，從而獲得更高的精度。學者們一直關注間斷有限元的優化問題，以期在求解波動方程時達到較小的誤差和較好的收斂性^[65]。除此之外，Yang等^[66]提出一種有限元-有限差分法，能夠較好地處理復雜地下結構的數值模擬問題，其中：有限元法可以更好地表示復雜地下結構的幾何特征，適應不規則邊界；有限差分法保持了計算的效率。

綜上所述，有限元法是一種常用的地震波動方程數值模擬方法，在理論和實踐上都得到了驗證，被廣泛接受并應用。首先，有限元法可以對復雜的地震波動問題進行較為精確的模擬。其次，有限元法靈活性強，可以適應各種實際情況，包括不規則地質結構。再次，間斷有限元法通過引入間斷近似空間和數值通量來處理解的間斷問題，能夠處理復雜地質結構、適應解的間斷特性、具有高精度和數值穩定性等優勢。但這類方法的計算成本通常會高于傳統的連續有限元法，特別是在高階逼近和復雜幾何形狀的情況。

然而，由于要將區域和結構劃分成許多小單元進行計算，有限元法計算量較大，尤其是在模擬大型結構或復雜地震場地時。此外，有限元法建立在一定的假設和簡化條件下，其模擬結果受建模誤差的影響，可能與實際情況存在一定差異。有限元法在地震波數值模擬方面可以朝著更加精細化的模型發展，以在考慮地下結構的復雜性時，包括地下的不規則邊界、地層的變化等，更好地模擬實際情況。

4"譜元法

譜元法是一種融合有限元法和譜方法的數值模擬技術。該方法兼具了有限元法處理復雜幾何邊界問題上的靈活性以及譜方法在高精度解決問題時的優點^[67]。利用譜元法進行求解時首先需要選擇合適的基函數（例如Chebyshev多項式、Legendre多項式等）以構建高階插值多項式，這些基函數在整個求解域內滿足正交性和完備性；然后在每個小單元內逼近待求解的函數。這與偽譜法中利用插值方法在離散點處逼近待求解函數的步驟類似。

與偽譜法不同的是，譜元法得到逼近解之后，不需要進行傅里葉變換，相反，它更類似于有限元法，將插值得到的近似解代入偏微分方程中，并將局部剛度矩陣疊加形成全局剛度矩陣。需要注意的是，譜元法與普通有限元法在處理全局質量矩陣時方法不同，譜元法首先通過數學處理，使全局質量矩陣成為對角矩陣，極大地提高了計算效率。然后按照需求添加震源項。與有限元法處理方式類似，整理得到

M+KU=Fs。（7）

最后選擇合適的數值求解方法求解全局方程組，得到各個時間步長上的位移場。由于結合了有限元法處理復雜結構的幾何靈活性和偽譜法高精度、快速收斂的特性，目前譜元法已成為地震波動方程數值模擬的重要工具。

Patera^[68]最初提出了譜元法，并將該方法應用于解決流體力學問題。在此基礎上，Seriani等^[69]應用譜元法研究聲波在復雜地質結構中的傳播問題。相比于一些傳統方法，譜元法的數值模擬結果更準確，可以較好地描述不同介質之間的不規則界面。Komatitsch等^[70]基于二維、三維模型，考慮復雜的幾何形狀、自由表面和材料界面，運用譜元法對真實地質結構中的彈性波傳播規律進行了模擬研究。

早期的譜元法在進行地震波動方程數值模擬時往往采用Chebyshev正交多項式對波動方程中的空間變量進行離散化處理。Lagrange譜元法是一種常用的譜元法，使用高次多項式作為插值函數，能夠適應各種復雜的幾何形狀，并且可以提供較高的數值模擬精度。但常規Lagendre譜元法在求解地震波運動方程時存在精度損失的情況。針對這一問題，劉少林等^[71]提出了一種提升譜元法模擬精度的優化算法，數值模擬結果驗證了其可行性和有效性。Liu等^[72]對有限元法和譜元法在地震波場建模中的應用效果進行了對比和分析，分別列舉了二階基于四邊形單元的有限元法和譜元法，以及三階基于三角形單元的有限元法和譜元法四種方法在求解特定接收器位置的計算誤差（包括水平方向和垂直方向分量的最大范數），還比較了這四種方法在半空間模型下CPU時間消耗、總時間步數、平均每步耗時、計算機內存占用（表2）。由2表可知：基于有限元法得到的計算結果在垂直方向上的誤差多大于水平方向上的誤差，然而總體上基于有限元法得到的結果與基于譜元法相當。在使用半空間模型進行數值模擬計算時，采用三角形單元三階多項式的方法需要更多的計算時間，并且計算機內存占用量與節點數量和多項式階數成正比。

近年來，針對譜元法優化方法的研究逐漸增多。Lyu等^[73]采用更高階的譜元方法，準確地模擬了彈性波傳播的物理過程，比低階譜元法更顯著地節省計算內存。Zhang等^[74]提出一種結合譜元法與有限

元法的SEM-FEM數值模擬方法，先使用譜元法在粗網格上模擬整個區域的地震波傳播過程，在此結果上利用有限元法在細網格上對局部地質和地形條件下的地震波場進行分析。劉立民等^[75]提出了適用于起伏地表復雜構造成像的Chebyshev譜元逆時偏移成像算法，取得了較好的應用效果。邢浩潔等^[76]將一種多人工波速優化透射邊界應用于高精度譜元法的地震波動模擬中，將該邊界條件與其他常見的邊界條件進行了比較，結果表明這是一種便捷、高效的人工邊界條件。孟雪莉等^[77]提出了一種基于優化數值積分的譜元法（optimal SEM， OSEM），該方法不僅能夠適應各類邊界條件，而且在壓制數值頻散、提高計算精度等方面具備一定的優勢。隨著研究不斷深入，譜元法在地震波動方程數值模擬應用中表現越來越出色，其模擬精度更高，計算速度更快。

綜上所述，譜元法在地震波動方程數值模擬中具有顯著的優勢。通過使用高次多項式作為插值函數，譜元法不僅能夠提供高精度的地震波傳播模擬結果，而且能夠準確地捕捉地震波傳播過程。譜元法能夠適應各種復雜的地形和非均勻性介質，在地震波模擬中具有廣泛的適用性。此外，譜元法基于頻域分析，能夠直接獲得地震波的頻譜信息，有助于進一步研究地震波的頻率特性。然而，譜元法也面臨一些挑戰。首先，由于需要進行頻域計算和求解大規模矩陣特征值，譜元法的計算復雜度較高，可能需要更多的計算資源和時間。其次，譜元法在數值模擬中高頻區域可能會存在的耗散誤差，會影響地震波傳播的高頻部分。最后，對于邊界條件的處理更為復雜，特別是對于非完全吸收邊界條件，可能需要采用額外的技術來減小邊界反射和影響。未來，譜元法將繼續發展并在地震波數值模擬中發揮重要作用。通過

開發更高階的插值函數減少數值耗散和頻散，開發更高效的算法減少計算復雜度，提高計算效率。

5"基于物理信息神經網絡的深度學習法

PINN是一種基于神經網絡，并引入物理信息約束以求取偏微分方程近似解的方法。Lagaris等^[78]在1998年使用人工智能求解常微分方程和偏微分方程，但受限于當時的計算水平，并未發展起來。直到Raissi等^[79]使用深度神經網絡來逼近函數，并根據物理方程和定解條件對網絡的預測結果進行約束，該方法才得到了顯著的進步和應用。傳統的數值求解方法只能在給定初始條件和邊界條件時才能進行求解，如果初始條件和邊界條件發生變化，就要重新進行求解^[80]。PINN最基礎的設計思想是將物理方程中的微分形式約束條件融入神經網絡的損失函數設計中。以二維波動方程為例，將空間坐標（x，y）和時間t作為輸入，將波函數u（x，y，t）的近似值作為輸出，損失函數的總體L=L_data+L_physics。為確保神經網絡模型能夠逼近已知的數據點，包括初始條件和邊界條件，損失函數通常采用歐氏距離的平方來衡量預測值與實際值之間的誤差。數據驅動部分的損失函數可以表示為

L_data=∑[u_predicted（x，y，t）－u_actual（x，y，t）]。（8）

式中：u_precticed為神經網絡預測的波函數；u_actual為真實的波函數。物理模型驅動部分的損失函數為

L_physics=∑2u_predictedt²－c²²u_predicted²。 （9）

式中，²表示拉普拉斯算子，目的是使預測的波函數二階時間導數和空間拉普拉斯算子與物理方程對應部分之間的差異最小化，確保神經網絡在訓練過程中逐漸調整自身參數，預測出符合物理方程的解，從而得到符合物理規律的結果。

PINN遵守訓練數據的分布規律以及偏微分方程所描述的物理定律。結合深度學習和物理模型的優勢，只需少量的訓練樣本就能訓練出滿足物理約束條件的神經網絡模型，模擬結果更逼近真實的地震波場^[81]。PINN作為一種求解波動方程的無網格方法，具有較高的精度，但計算效率較低、穩定性不高、適用性較低。隨著研究的不斷深入，越來越多學者使用PINN進行地震波場數值模擬，并針對存在的不同問題進行優化^[82]。

Guo等^[83]改進傳統的PINN算法，使神經網絡模型的性能進行進一步提升，成功提升了一維聲波方程的求解結果精度。Ovadia等^[84]考慮了聲波問題不滿足數值穩定性條件的情況，提出了一種利用PINN的方法以提高聲波方程求解的穩定性，并通過實驗驗證了這一方法的有效性。Song等^[85]提出了一種利用PINN求解垂直各向異性介質聲波方程散射波場的方法，采用全連接深度神經網絡來構建模型，并通過訓練來優化網絡參數，對空間中每個位置散射波場的實部和虛部進行預測，結果表明，在大部分區域，數值解和PINN預測解的散射波場的總體形狀非常接近（圖4）。其中，較小的散射波場差異主要來自缺失的散射分類。數值模擬結果表明，在節省計算資源的情況下，該方法能夠得到沒有數值頻散偽影的聲波波場。Rasht-Behesht等^[86]提出了一種利用PINN求解二維聲波方程的算法，

并使用模擬數據進行測試;結果表明，PINN能夠處理不同的復雜構造。針對PINN在近似高頻分量時存在頻譜偏置的問題，Song等^[87]采用具有嵌入式傅里葉特征的PINN模擬地震多頻波場，數值模擬結果表明該方法與簡單PINN相比，在生成多頻波場的準確性方面有顯著改進;此外，該方法能有效模擬具有不同頻率和源位置的波場。Zhang等^[88]運用PINN進行地震反演，將一階聲波方程嵌入到損失函數中，將此作為訓練神經網絡的正則化項，以估計地下速度和密度場，結果表明，該方法能有效預測地震波在地下傳播的過程及對應的速度和密度，并在未添加吸收邊界條件下成功避免了地震波在遇到邊界時產生的反射現象，從而提高了模擬的真實性和可靠性。

在地震波動方程數值模擬中，模擬地震波的彈性行為非常重要。與求解聲波方程相比，求解彈性波動方程往往需要更大的計算成本。傳統數值模擬方法需要在有限點的網格上進行建模，這些方法的適用性隨著維度數量的增加而降低。相對而言，采用PINN進行數值模擬能有效解決這一問題。Rao等^[89]提出了一種具有混合變量輸出且不依賴于訓練樣本的PINN方法求解彈性動力學問題，該方法將位移和應力分量同時作為網絡輸出，極大地提高了網絡模型的可訓練性和預測結果的準確性。

PINN能以較高的效率提供微分方程的無網格解，但難以有效精準地解決大規模、多尺度的問題。對此，Moseley等^[90]以經典有限元法為基礎，提出了一種基于有限元和物理神經網絡的深度學習方法（finite basis PINN， FBPINN），在每個子域上將輸入進行單獨的歸一化來解決神經網絡的譜偏差問題，通過并行分治的方法來降低基礎優化問題的復雜性；數值模擬結果表明，該方法的目標函數在訓練過程中能夠較快收斂，且其預測結果的精度較高。陳蘇等^[91]使用物理驅動的深度學習方法結合譜元法形成的稀疏初始波場數據對二維波動問題進行數值模擬，成功處理了具有自由邊界和復雜地形的波場計算。他們改變不同初始條件進行數值模擬檢驗神經網絡的泛化精度，并提出了一種遷移學習策略，該策略可以高精度地預測在無限介質中不同源位置的波場。與譜元法數值模擬結果對比，驗證了該方法在模擬均質介質、空間非均勻性以及復雜地形下的波傳播過程中具有可靠性。

總體來說，基于PINN的深度學習法為地震波場數值模擬提供了一種新的思路和技術手段。神經網絡具備強大的學習能力，可以通過訓練來自動調整模型參數，適應不同的地震事件和地質條件。但該方法也存在一些限制。首先，獲取大量的地震數據用于訓練比較困難，這包括收集地震波場觀測數據和相關的地質參數。此外，深度學習方法對計算資源的需求很龐大，需要充足的計算設備和存儲空間。還需要對其參數進行細致的調整和優化，這通常需要專家具有豐富的專業知識和經驗。

6"討論

上述五種方法的適用范圍、計算精度、計算效率、計算內存、數值穩定性和邊界條件處理能力如表3所示。有限差分法適用于均勻介質和常見邊界條件，在求解波動方程時具有較高的穩定性，但其計算效率和精度受網格分辨率和時間步長的影響；偽譜法

具有高精度和低內存需求的特點

，但計算效率較低，數值模擬結果受空間分辨率影響，在高頻區域可能不穩定

且依賴于合適的邊界條件；有限元法適合處理復雜地質結構，計算精度較高，但其計算效率較低和且內存需求較大；譜元法適用于解決高精度地震波傳播模擬問題，靈活性較好，但其數值模擬過程較為復雜；基于PINN的深度學習法擅長處理復雜非線性問題，具有自適應性，但其訓練過程耗時高且依賴網絡結構與參數調優。

五種方法的優點、挑戰和未來的發展方向如表4所示。未來：有限差分法需要發展新的高精度算法以解決數值頻散和數值耗散問題，提高模擬精度；偽譜法可以改進邊界條件處理方法以增強適用性，并智能調整網格分辨率來提高模擬結果的精度；有限元法著重提高計算效率，并優化計算資源利用，以適應更復雜的地震波動問題；譜元法在保證精度的前提下，開發更高效的算法以減少計算復雜度，提高計算效率；基于PINN的深度學習法需獲取大量地震數據進行訓練、提升計算資源，進行參數調整和優化，并結合最新網絡結構及訓練算法，提高計算效率以及模擬結果的準確性和可靠性。

通過對Web of Science核心期刊中地球物理學領域五種數值求解方法的年發表量進行深入的統計分析（圖5），可以得出一些有益結論：1）有限差分法是五種方法中研究頻率最高的數值模擬方法，在學術界和工業界得到了廣泛的認可和應用。這可能是由于有限差分法易于理解且計算效率較高，在處理復雜地質問題時靈活性和有效性較高。2）有限元法與譜元法的研究頻率緊隨其后。這可能是由于在實際生產中，數值求解方法面臨復雜性和計算成本的挑戰。有限元法能夠靈活地適應復雜的地形和地質邊界條件，譜元法可以實現對波動方程的高效求解。這些優勢導致它們受到學者和專家們的持續關注，其應用和研究頻率逐年增加。3）偽譜法的研究頻率較低，這可能是由于其在解決某些問題時的收斂性和精度存在局限性。盡管如此，偽譜法理論上仍具有較高的模擬精度，其發展潛力和應用價值仍值得期待。4）近年來，隨著計算機技術的飛速發展，基于PINN的深度學習法開始嶄露頭角。盡管目前相關文章的發行量較少，但其增長趨勢表明了該方法在地球物理學領域的潛力。基于PINN的深度學習法有望為地震波動方程數值模擬提供新的視角和解決方法。總體來說，這些數值求解方法在地球物理學領域都有獨特的地位和價值。隨著科技的進步和研究需求的不斷變化，它們將會在這個領域發揮越來越重要的作用。

7"結論和展望

波動方程包含豐富的動力學和運動學信息，基于波動方程的地震波數值模擬對地震資料處理、解釋和地震勘探至關重要。隨著研究的深入和計算機算力的提高，地震波數值模擬方法在理論上不斷創新，在應用上也取得了較多成果，但這些方法也存在著不足之處。如有限差分法雖然易于理解且計算效率較高，但通常會產生數值頻散問題，影響模擬精度；偽譜法雖然計算精度高且內存需求低，但在空間分辨率和邊界條件處理上存在一定的限制，且在高頻情況下可能出現振蕩現象，影響結果的可靠性；有限元法在復雜地下介質中能獲得較高精度的波場模擬結果，但其計算量較大，且低階數值模擬精度較低；譜元法雖然結合了有限元法和偽譜法的優點，能夠較好地處理復雜介質數值模擬問題，但在求解大規模矩陣特征值時計算復雜度較高，計算效率較低；基于物理信息神經網絡的深度學習法在數值模擬中具有高效、靈活等優點，能夠處理復雜介質的數值模擬問題，但需獲取大量地震數據進行訓練，并需要較高的計算資源，此外，該方法難以有效地處理振幅較小的數據，且其損失函數的適用性有限，限制了其在實際問題中的應用。

基于上述分析，本文對地震波動方程數值模擬方法的發展趨勢作如下展望：

1）在數值模擬中使用更為真實的地層模型和介質參數，以提高模擬精度；2）優化邊界條件，更為精確地模擬邊界反射情況；3）結合多物理場耦合理論，更為真實地模擬復雜介質中地震波的傳播情況；4）利用并行算法等高效計算機技術，提高數值模擬的效率；5）更為深入地結合PINN等神經網絡模型，提升數值模擬的性能。

參考文獻（References）：

[1] 李信富， 李小凡， 張美根. 地震波數值模擬方法研究綜述[J].防災減災工程學報， 2007， 27（2）： 241-248.

Li Xinfu， Li Xiaofan， Zhang Meigen.Review of Seismic Wave Numerical Modeling Methods[J]. Journal of Disaster Prevention and Mitigation Engineering， 2007， 27（2）： 241-248.

[2] 皮紅梅， 蔣先藝， 劉財， 等. 波動方程數值模擬的三種方法及對比[J].地球物理學進展， 2009， 24（2）： 391-397.

Pi Hongmei， Jiang Xianyi， Liu Cai， et al.Three Methods of Wave Equation Forward Modeling and Comparison[J]. Progress in Geophysics， 2009， 24（2）： 391-397.

[3] 邴琦， 孫章慶， 韓復興， 等. 地震波射線追蹤方法綜述：方法、分類、發展現狀與趨勢[J].地球物理學進展， 2020， 35（2）： 536-547.

Bing Qi， Sun Zhangqing， Han Fuxing， et al. Review on the Seismic Wave Ray Tracing： Methods， Classification， Development Present and Trends[J]. Progress in Geophysics， 2020， 35（2）： 536-547.

[4] Ouyang F， Chen L. Iterative Magnetic Forward Modeling for High Susceptibility Based on Integral Equation and Gauss-Fast Fourier Transform[J]. Geophysics， 2020， 85（1）： J1-J13.

[5] Li C， Liu J， Chen B， et al. Optimal Third‐Order Symplectic Integration Modeling of Seismic Acoustic Wave Propagation[J]. Bulletin of the Seismological Society of America， 2020， 110（2）： 754-762.

[6] Chen L， Huang J， Fu L Y， et al. A Compact High-Order Finite-Difference Method with Optimized Coefficients for 2D Acoustic Wave Equation[J]. Remote Sensing， 2023， 15（3）： 604.

[7] Cheng J W， Fan N， Zhang Y Y， et al. Irregular Surface Seismic Forward Modeling by a Body-Fitted Rotated-Staggered-Grid Finite-Difference Method[J]. Applied Geophysics， 2018， 15（3/4）： 420-431.

[8] 劉炯， 趙殿棟， 劉喜武. 起伏地表下的地震波偽譜法模擬 [C]//中國石油學會石油物探專業委員會， 中國地球物理學會勘探地球物理委員會. 2022年中國石油物探學術年會論文集（下）. 海口： 中國石油學會石油物探專業委員會， 中國地球物理學會勘探地球物理委員會， 2022： 377-380.

Liu Jiong， Zhao Diandong， Liu Xiwu. Pseudo Spectral Simulation of Seismic Waves Under Undulating Surface[C]//Petroleum Geophysical Exploration Committee of "Chinese Petroleum Society， Exploration Geophysical Committee of "Chinese Geophysical Society. Proceedings of the 2022 Annual Academic Conference on Petroleum Exploration and Production in China （Part II）. Haikou：

Petroleum Geophysical Exploration Committee of "Chinese Petroleum Society， Exploration Geophysical Committee of "Chinese Geophysical Society， 2022： 377-380.

[9] 黃旭東， 田源. 基于擴展有限元法的水力壓裂數值模擬研究[J].水利技術監督， 2021（6）： 149-155.

Huang Xudong， Tian Yuan. Numerical Simulation of Hydraulic Fracturing Based on Extended Finite Element Method[J].Water Resources Technical Supervision， 2021， 2021（6）： 149-155.

[10] 董興朋， 楊頂輝. 球坐標系下譜元法三維地震波場模擬[J].地球物理學報， 2017， 60（12）： 4671-4680.

Dong Xingpeng， Yang Dinghui.Numerical Modeling of the 3-D Seismic Wavefield with the Spectral Element Method in Spherical Coordinates[J]. Chinese Journal of Geophysics， 2017， 60（12）： 4671-4680.

[11]

Cai S， Mao Z， Wang Z， et al. Physics-Informed Neural Networks （PINNs） "for Fluid Mechanics： A Review[J]. Acta Mechanica Sinica， 2021，37（12）： 1727-1738.

[12] 朱多林， 白超英. 基于波動方程理論的地震波場數值模擬方法綜述[J].地球物理學進展， 2011， 26（5）： 1588-1599.

Zhu Duolin， Bai Chaoying. Review on the Seismic Wavefield Forward Modeling[J].Progress in Geophysics， 2011， 26（5）： 1588-1599.

[13] Virieux J， Calandra H， Plessix R E. A Review of the Spectral， Pseudo-Spectral， Finite-Difference， and Finite-Element Modelling Techniques for Geophysical Imaging[J]. Geophysical Prospecting， 2011， 59（5）： 794-813.

[14] Bouslous H， Boujaoui E H， Maniar L. Uniform Boundary Stabilization for the Finite Difference Semi-Discretization of 2-D Wave Equation[J]. Afrika Matematika， 2014， 25（3）： 623-643.

[15] Alterman Z， Karal F C. Propagation of Elastic Waves in Layered Media by Finite Difference Methods[J]. Bulletin of the Seismological Society of America， 1968， 58（1）： 367-398.

[16] Boore D M. Finite Difference Methods for Seismic Wave Propagation in Heterogeneous Materials[J]. Methods in Computational Physics： Advances in Research and Applications， 1972， 11： 1-37.

[17] Madariaga R.Dynamics of an Expanding Circular Fault[J]. Bulletin of the Seismological Society of America， 1976， 66（3）： 639-666.

[18] Dablain M A.The Application of High-Order Differencing to the Scalar Wave Equation[J]. Geophysics， 1986， 51： 54-66.

[19] Liang W， Wang Y， Cao J， et al. A Hybrid Explicit Implicit Staggered Grid Finite-Difference Scheme for the First-Order Acoustic Wave Equation Modeling[J]. Scientific Reports， 2022， 12（1）：10967.

[20] 唐超， 文曉濤， 王文化. 基于最小范數優化交錯網格有限差分系數的波動方程數值模擬[J].石油地球物理勘探， 2021， 56（5）： 1039-1047.

Tang Chao， Wen Xiaotao， Wang Wenhua. Numerical Simulation of Wave Equations Based on Minimum-Norm Optimization of Staggered-Grid Finite-Difference Coefficients[J]. Oil Geophysical Prospecting， 2021， 56（5）： 1039-1047.

[21] Zhou H， Liu Y， Wang J. Time-Space Domain Scalar Wave Modeling by a Novel Hybrid Staggered-Grid Finite-Difference Method with High Temporal and Spatial Accuracies[J]. Journal of Computational Physics， 2022， 455： 111004.

[22] 馮德山，王向宇．基于卷積完全匹配層的旋轉交錯網格高階差分法模擬彈性波傳播[J]．物探與化探，2018，42（4）： 766-776．

Feng Deshan， Wang Xiangyu. Elastic Wave Propagation Simulation in Anisotropic Media and Random Media Using High-Order Difference Method of Rotation Staggered Grids Based on Convolutional Perfectly Matched Layer[J]. Geophysical and Geochemical Exploration， 2018， 42（4）： 766-776.

[23] 張志佳， 孫章慶， 王瑞湖，等. 復雜地表起伏多重變加密網格有限差分法波動方程數值模擬[J].吉林大學學報（地球科學版）， 2023， 53（2）： 598-608.

Zhang Zhijia， Sun Zhangqing， Wang Ruihu， et al.Numerical Simulation of Wave Equation Using Finite Difference Method with Rugged Variable Refined Multigrid in Complex Surface[J]. Journal of Jilin University （Earth Science Edition）， 2023， 53（2）： 598-608.

[24] Zhang C， Zhang W.Efficient 2D Acoustic Wave Finite-Difference Numerical Simulation in Strongly Heterogeneous Media Using the Adaptive Mesh Refinement Technique[J]. Geophysics， 2022， 87（1）： T29-T42.

[25] 于明浩， 鞏向博， 萬曉杰. 復雜地表相控震源正演模擬與照明分析[J]. 吉林大學學報（地球科學版）， 2022， 52（6）： 2034-2041.

Yu Minghao， Gong Xiangbo， Wan Xiaojie. Seismic Simulation and Illumination Analysis of Phased-Array Vibrator System on Complex Surface[J]. Journal of Jilin University （Earth Science Edition）， 2022， 52 （6）： 2034-2041.

[26] 楊尚倍， 白超英， Zhou Bing. 時間域廣義2.5D一階波動方程曲網格有限差分法數值模擬[J].石油地球物理勘探， 2021， 56（6）： 1262-1278.

Yang Shangbei， Bai Chaoying， Zhou Bing.Curvillinear-Grid Finite-Difference Numerical Simulation Method for Generalized First-Order 2.5D Time Domain Wave Equation[J]. Oil Geophysical Prospecting， 2021， 56（6）： 1262-1278.

[27] 范娜， 成景旺， 秦雷， 等. 一種優化的頻率域三維聲波有限差分模擬方法[J].地球物理學報， 2018， 61（3）： 1095-1108.

Fan Na， Cheng Jingwang， Qin Lei， et al.An Optimized Method for Frequency-Domain Finite-Difference Solution of 3D Scalar Wave Equation[J]. Chinese Journal of Geophysics， 2018， 61（3）： 1095-1108.

[28] Poorya H， Mark G， Chorsi H T， et al.Using Eccentricity to Locate Ionospheric Exit Points of Magnetospheric Whistler Mode Waves[J]. IEEE Transactions on Geoscience and Remote Sensing， 2018， 56（12）： 1-13.

[29] 黃建平， 婁璐烽， 彭煒颋， 等. 一種基于拉格朗日乘子的空間域差分系數優化方法[J].中國石油大學學報（自然科學版）， 2022， 46（4）： 30-40.

Huang Jianping， Lou Lufeng， Peng Weiting， et al. An Optimization Method of Finite Difference Coefficient in Spatial Domain Based on Lagrange Multipliers[J]. Journal of China University of Petroleum （Edition of Natural Science）， 2022， 46（4）： 30-40.

[30] 彭更新， 劉威， 郭念民， 等. 基于時空域交錯網格有限差分法的應力速度聲波方程數值模擬[J].石油物探， 2022， 61（1）： 156-165.

Peng Gengxin， Liu Wei， Guo Nianmin， et al.A Time-Space Domain Dispersion-Relationship-Based Staggered-Grid Finite-Difference Scheme for Modeling the Stress-Velocity Acoustic Wave Equation[J]. Geophysical Prospecting for Petroleum， 2022， 61（1）： 156-165.

[31] 徐世剛， 包乾宗， 任志明， 等. 基于時間高精度-空間隱式矩形交錯網格有限差分的彈性波數值模擬[J].地球物理學報， 2022， 65（4）： 1389-1401.

Xu Shigang， Bao Qianzong， Ren Zhiming， et al.Simulation Elastic Wave Using Temporal High Accuracy and Implicit Spatial Rectangular Staggered-Grid Finite-Difference Approaches[J]. Chinese Journal of Geophysics， 2022， 65（4）： 1389-1401.

[32] Li Q， Wu G， Jia Z， et al. Full-Waveform Inversion in Acoustic-Elastic Coupled Media with Irregular Seafloor Based on the Generalized Finite-Difference Method[J]. Geophysics， 2023， 88（2）： T83-T100.

[33] Kosloff D， Pestana R C， Tal-Ezer H.Acoustic and Elastic Numerical Wave Simulations by Recursive Spatial Derivative Operators[J]. Geophysics， 2010， 75（6）： T167-T174.

[34] Li K， Liao W. An Efficient and High Accuracy Finite-Difference Scheme for the Acoustic Wave Equation in 3D Heterogeneous Media[J]. Journal of Computational Science， 2020， 40： 101063.

[35] Chen L， Huang J， Fu L Y， et al. A Compact High-Order Finite-Difference Method with Optimized Coefficients for 2D Acoustic Wave Equation[J]. Remote Sensing， 2023， 15（3）： 604.

[36] Guo X， Wang J， Yang S， et al. Optimal Staggered-Grid Finite-Difference Method for Wave Modeling Based on Artificial Neural Networks[J]. Computers amp; Mathematics with Applications， 2022， 108： 141-158.

[37] 韓復興， 王若雯， 孫章慶， 等. 地震聲波數值模擬中人工邊界條件的差別與組合[J]. 吉林大學學報（地球科學版）， 2022， 52（1）： 261-269.

Han Fuxing， Wang Ruowen， Sun Zhangqing， et al. Difference and Combination of Artificial Boundary Conditions in Seismic Acoustic Numerical Simulation[J]. Journal of Jilin University （Earth Science Edition）， 2022， 52 （1）： 261-269.

[38] 趙明哲， 楊軍， 張陸軍， 等. 基于CUDA的高階旋轉交錯網格有限差分法的彈性波正演模擬[J].地球物理學進展， 2022， 37（4）： 1697-1703.

Zhao Mingzhe， Yang Jun， Zhang Lujun， et al. Forward Modeling of Elastic Wave Based on CUDA Finite Difference Method with High Order Rotating Staggered Grid[J]. Progress in Geophysics， 2022， 37（4）： 1697-1703.

[39] Kreiss H， Oliger J. Comparison of Accurate Methods for the Integration of Hyperbolic Equation[J]. Tellus， 1972， 24（3）： 199-215.

[40] Kosloff D D， Baysal E.Forward Modeling by a Fourier Method[J]. Geophysics， 1982， 47（10）： 1402-1412.

[41] Reshef M， Kosloff D， Edwards M， et al.Three-Dimensional Elastic Modeling by the Fourier Method[J]. Geophysics， 1988， 53（9）： 1184-1193.

[42] 孫獻果， 張東. 地震波模擬中有限差分法與偽譜法的對比研究[J].中國科技論文， 2018， 13（17）： 2005-2008.

Sun Xianguo， Zhang Dong. Study on the Differences Between Finite Difference Method and Pseudo-Spectral Algorithm in Seismic Numerical Modeling[J]. China Science and Technology Paper， 2018， 13（17）： 2005-2008.

[43] 張慶朝， 朱國維， 周俊杰， 等. TTI介質qP波偽譜法正演模擬[J].石油地球物理勘探， 2019， 54（2）： 302-311， 236.

Zhang Qingchao， Zhu Guowei， Zhou Junjie， et al. qP-Wave Numerical Simulation in TTI Media Using Pseudo-Spectral Method[J]. Oil Geophysical Prospecting， 2019， 54（2）： 302-311， 236.

[44] Gunes A， Saydam T， Aksoy S. Stability and Accuracy of Engquist-Majda Absorbing Boundary Condition for Pseudo Spectral Time Domain Method[J]. International Journal of Numerical Modelling： Electronic Networks， Devices and Fields， 2022， 35（3）： e2971.

[45] 譚文卓， 吳幫玉， 李博， 等. 梯形網格偽譜法地震波場模擬[J].石油地球物理勘探， 2020， 55（6）： 1282-1291.

Tan Wenzhuo， Wu Bangyu， Li Bo， et al. Seismic Wave Simulation Using a Trapezoidal Grid Pseudo-Spectral Method[J]. Oil Geophysical Prospecting， 2020， 55（6）： 1282-1291.

[46] 李元燮， 劉財， 郭智奇， 等. 基于旋轉交錯網格偽譜法和L-S理論的熱彈性介質波場模擬與特征分析[J].中國科學： 地球科學， 2021， 51（8）： 1355-1374.

Li Yuanxie， Liu Cai， Guo Zhiqi， et al. Simulation and Characteristics Analysis of a Wavefield in a Thermoelastic Medium Adopting the Rotated Staggered-Grid Pseudo-Spectral Method and L-S Theory[J]. Science China Earth Sciences， 2021， 64（8）： 1390-1408.

[47] 李青陽， 吳國忱， 梁展源. 基于PML邊界的時間高階偽譜法彈性波場模擬[J].地球物理學進展， 2018， 33（1）： 228-235.

Li Qingyang， Wu Guochen， Liang Zhanyuan.Time Domain High-Order Pseudo Spectral Method Based on PML Boundary for Elastic Wavefield Simulation[J]. Progress in Geophysics， 2018， 33（1）： 228-235.

[48] 馬銳， 鄒志輝， 芮擁軍， 等. 基于SPML和海綿邊界的偽譜法彈性波模擬復合吸收邊界條件[J].石油物探， 2018， 57（1）： 94-103.

Ma Rui， Zou Zhihui， Rui Yongjun， et al.A Composite Absorbing Boundary Based on the SPML and Sponge Absorbing Boundary for Pseudo-Spectral Elastic Wave Modeling[J]. Geophysical Prospecting for Petroleum， 2018， 57（1）： 94-103.

[49] 陳漢明， 汪燚林， 周輝. 一階速度—壓力常分數階黏滯聲波方程及其數值模擬[J].石油地球物理勘探， 2020， 55（2）： 302-310.

Chen Hanming， Wang Yilin， Zhou Hui.A Novel Constant Fractional-Order Laplacians Viscoacoustic Wave Equation and Its Numerical Simulation Method[J]. Oil Geophysical Prospecting， 2020， 55（2）： 302-310.

[50] 李展輝， 曹學峰， 朱自強. 基于GPU并行加速的探地雷達偽譜法三維正演模擬[J].工程地球物理學報， 2021， 18（1）： 67-74.

Li Zhanhui， Cao Xuefeng， Zhu Ziqiang.A Three-Dimensional Forward Modeling with Pseudo-Spectral Method for Ground Penetrating Radar Accelerated by GPU[J]. Chinese Journal of Engineering Geophysics， 2021， 18（1）： 67-74.

[51] Lysmer J， Drake L A.A Finite Element Method for Seismology[J]. Methods in Computational Physics， 1972， 11（4）： 181-216.

[52] Drake L A. "Rayleigh Waves at a Continental Boundary by the Finite Element Method[J]. Bulletin of the Seismological Society of America， 1972， 62（5）： 1259-1268.

[53] Mullen R， Belytschko T. Dispersion Analysis of Finite Element Semidiscretizations of the Two-Dimensional Wave Equation[J]. International Journal for Numerical Methods in Engineering， 1982， 18（1）： 11-29.

[54] Chen K H. Propagating Numerical Model of Elastic Wave in Anisotropic in Homogeneous Media-Finite Element Method[C]// Symposium of 54th SEG， Annual Meeting Expanded Abstracts. Houston： [s. n.]， 1984： 631-632.

[55] Serón F J， Sanz F J， Kindelan M， et al. Finite-Element Method for Elastic Wave Propagation[J]. Communications in Applied Numerical Methods， 1990， 6（5）： 359-368.

[56] Jagota V， Sethi A P S， Kumar K.Finite Element Method： An Overview[J]. Walailak Journal of Science and Technology， 2013， 10（1）： 1-8.

[57] Wu Y， Meng B， Wang L， et al. Seismic Vulnerability Analysis of Buried Polyethylene Pipeline Based on Finite Element Method[J]. International Journal of Pressure Vessels and Piping， 2020， 187： 104167.

[58] 蘇波， 李懷良， 劉少林， 等. 修正辛格式有限元法的地震波場模擬[J].地球物理學報， 2019， 62（4）： 1440-1452.

Su Bo， Li Huailiang， Liu Shaolin， et al. Modified Symplectic Scheme with Finite Element Method for Seismic Wavefield Modeling[J]. Chinese Journal of Geophysics， 2019， 62（4）： 1440-1452.

[59] Zhao H， Xu W， Gao J， et al. A Finite-Element Algorithm with a Perfectly Matched Layer Boundary Condition for Seismic Modelling in a Diffusive-Viscous Medium[J]. Journal of Geophysics and Engineering， 2022， 19（1）： 51-66.

[60] Reed W H， Hill T R.Triangular Mesh Methods for the Neutron Transport Equation[R]. Los Alamos： Los Alamos Scientific Laboratory， 1973.

[61] 賀茜君， 楊頂輝， 仇楚鈞， 等. 基于非結構網格求解三維D’Alembert介質中聲波方程的并行加權Runge-Kutta間斷有限元方法[J].地球物理學報， 2021， 64（3）： 876-895.

He Qianjun， Yang Dinghui， Qiu Chujun， et al. A Parallel Weighted Runge-Kutta Discontinuous Galerkin Method for Solving Acoustic Wave Equation in 3D D’Alembert Media Based on Unstructured Meshes[J]. Chinese Journal of Geophysics， 2021， 64（3）： 876-895.

[62] He X J， Li J S， Huang X Y， et al. Solving Elastic Wave Equations in 2D Transversely Isotropic Media by a Weighted Runge-Kutta Discontinuous Galerkin Method[J]. Petroleum Science， 2023， 20（2）： 827-839.

[63] Hong C， Wang X， Zhao G， et al. Discontinuous Finite Element Method for Efficient Three-Dimensional Elastic Wave Simulation[J]. Journal of Geophysics and Engineering， 2021， 18（1）： 98-112.

[64] 韓德超， 劉衛華， 司文朋. 三角網格間斷有限元法彈性波模擬精度分析[J].石油地球物理勘探， 2021， 56（4）： 758-770.

Han Dechao， Liu Weihua， Si Wenpeng.Analysis of Elastic Wave Simulation Accuracy with Discontinuous Galerkin Finite Element Method Based on Triangular Meshes[J]. Oil Geophysical Prospecting， 2021， 56（4）： 758-770.

[65] 廉西猛， 張睿璇. 地震波動方程的局部間斷有限元方法數值模擬[J].地球物理學報， 2013， 56（10）： 3507-3513.

Lian Ximeng， Zhang Ruixuan. Numerical Simulation of Seismic Wave Equation by Local Discontinuous Galerkin Method[J].Chinese Journal of Geophysics， 2013， 56（10）： 3507-3513.

[66] Yang J， Liu T， Tang G， et al. Modeling Seismic Wave Propagation Within Complex Structures[J]. Applied Geophysics， 2009， 6（1）： 30-41.

[67] Shi L， Zhou Y， Wang J M， et al. Spectral Element Method for Elastic and Acoustic Waves in Frequency Domain[J]. Journal of Computational Physics， 2016， 327： 19-38.

[68] Patera A T. Spectral Element Method for Fluids Dynamics： Laminar Flow in a Channel Expansion[J]. Journal of Computational Physics， 1984， 54（3）： 468-488.

[69] Seriani G， Priolo E. Spectral Element Method for Acoustic Wave Simulation in Heterogeneous Media[J]. Finite Elements in Analysis and Design， 1994， 16（3/4）： 337-348.

[70] Komatitsch D， Vilotte J P.The Spectral Element Method： An Efficient Tool to Simulate the Seismic Response of 2D and 3D Geological Structures[J]. Bulletin of the Seismological Society of America， 1998， 88（2）： 368-392.

[71] 劉少林， 楊頂輝， 徐錫偉， 等. 模擬地震波傳播的三維逐元并行譜元法[J].地球物理學報， 2021， 64（3）： 993-1005.

Liu Shaolin， Yang Dinghui， Xu Xiwei， et al.Three-Dimensional Element-by-Element Parallel Spectral-Element Method for Seismic Wave Modeling[J]. Chinese Journal of Geophysics， 2021， 64（3）： 993-1005.

[72] Liu Y， Teng J， Lan H， et al. A Comparative Study of Finite Element and Spectral Element Methods in Seismic Wavefield Modeling[J]. Geophysics， 2014， 79（2）： T91-T104.

[73] Lyu C， Capdeville Y， Zhao L.Efficiency of the Spectral Element Method with Very High Polynomial Degree to Solve the Elastic Wave Equation[J]. Geophysics， 2020， 85（1）： T33-T43.

[74] Zhang L， Wang J T， Xu Y J， et al. A Procedure for 3D Seismic Simulation from Rupture to Structures by Coupling SEM and FEM[J]. Bulletin of the Seismological Society of America， 2020， 110（3）： 1134-1148.

[75] 劉立民， 劉定進， 李博. 起伏地表譜元逆時偏移方法[J].石油物探， 2021， 60（2）： 312-322.

Liu Limin， Liu Dingjin， Li Bo. Reversetime Migration for Irregular Topography Based on Spectral Element Method[J]. Geophysical Prospecting for Petroleum， 2021， 60（2）： 312-322.

[76] 邢浩潔， 劉愛文， 李小軍， 等. 多人工波速優化透射邊界在譜元法地震波動模擬中的應用[J].地震學報， 2022， 44（1）： 26-39.

Xing Haojie， Liu Aiwen， Li Xiaojun， et al. Application of an Optimized Transmitting Boundary with Multiple Artificial Wave Velocities in Spectral-Element Simulation of Seismic Wave Propagation[J]. Acta Seismologica Sinica， 2022， 44（1）： 26-39.

[77] 孟雪莉， 劉少林， 楊頂輝， 等. 基于優化數值積分的譜元法模擬地震波傳播[J].石油地球物理勘探， 2022， 57（3）： 602-612.

Meng Xueli， Liu Shaolin， Yang Dinghui， et al.Spectral-Element Method Based on Optimal Numerical Integration for Seismic Waveform Modeling[J]. Oil Geophysical Prospecting， 2022， 57（3）： 602-612.

[78] Lagaris I E， Likas A， Fotiadis D I.Artificial Neural Networks for Solving Ordinary and Partial Differential Equations[J]. IEEE Transactions on Neural Networks， 1998， 9（5）： 987-1000.

[79] Raissi M， Perdikaris P， Karniadakis G E.Physics-Informed Neural Networks： A Deep Learning Framework for Solving Forward and Inverse Problems Involving Nonlinear Partial Differential Equations[J]. Journal of Computational Physics， 2019， 378： 686-707.

[80] Krishnapriyan A， Gholami A， Zhe S， et al.Characterizing Possible Failure Modes in Physics-Informed Neural Networks[J]. Advances in Neural Information Processing Systems， 2021， 34： 26548-26560.

[81] 李野， 陳松燦. 基于物理信息的神經網絡：最新進展與展望[J].計算機科學， 2022， 49（4）： 254-262.

Li Ye， Chen Songcan.Physics-Informed Neural Networks： Recent Advances and Prospects[J]. Computer Science， 2022， 49（4）： 254-262.

[82] Rao C， Sun H， Liu Y.Physics-Informed Deep Learning for Incompressible Laminar Flows[J]. Theoretical and Applied Mechanics Letters， 2020， 10（3）： 207-212.

[83] Guo Y， Cao X， Liu B， et al. Solving Partial Differential Equations Using Deep Learning and Physical Constraints[J]. Applied Sciences， 2020， 10（17）： 5917.

[84] Ovadia O， Kahana A， Turkel E， et al.Beyond the Courant-Friedrichs-Lewy Condition： Numerical Methods for the Wave Problem Using Deep Learning[J]. Journal of Computational Physics， 2021， 442： 110493.

[85] Song C， Alkhalifah T， Waheed U B.Solving the Frequency-Domain Acoustic VTI Wave Equation Using Physics-Informed Neural Networks[J]. Geophysical Journal International， 2021， 225（2）： 846-859.

[86] Rasht-Behesht M， Huber C， Shukla K， et al. Physics-Informed Neural Networks （PINNs） for Wave Propagation and Full Waveform Inversions[J]. Journal of Geophysical Research： Solid Earth， 2022， 127（5）： e2021JB023120.

[87] Song C， Wang Y. Simulating Seismic Multifrequency Wavefields with the Fourier Feature Physics-Informed Neural Network[J]. Geophysical Journal International， 2023， 232（3）： 1503-1514.

[88] Zhang Y， Zhu X， Gao J.Seismic Inversion Based on Acoustic Wave Equations Using Physics-Informed Neural Network[J]. IEEE Transactions on Geoscience and Remote Sensing， 2023， 61： 1-11.

[89] Rao C， Sun H， Liu Y.Physics-Informed Deep Learning for Computational Elastodynamics Without Labeled Data[J]. Journal of Engineering Mechanics， 2021， 147（8）： 04021043.

[90] Moseley B， Markham A， Nissen-Meyer T.Finite Basis Physics-Informed Neural Networks （FBPINNs）： A Scalable Domain Decomposition Approach for Solving Differential Equations[J]. Advances in Computational Mathematics， 2023， 49（4）： 62.

[91] 陳蘇， 丁毅， 孫浩， 等. 物理驅動深度學習波動數值模擬方法及應用[J]. 力學學報， 2023， 55（1）： 272-282.

Chen Su， Ding Yi， Sun Hao， et al. Methods and Applications of Physical Information Deep Learning in Wave Numerical Simulation[J]. Chinese Journal of Theoretical and Applied Mechanics， 2023， 55（1）： 272-282.