交錯網格任意階導數有限差分格式及差分系數推導

楊慶節，劉 財，耿美霞，馮 晅，郭智奇，劉 洋

吉林大學地球探測科學與技術學院,長春 130026

交錯網格任意階導數有限差分格式及差分系數推導

楊慶節，劉 財，耿美霞，馮 晅，郭智奇，劉 洋

吉林大學地球探測科學與技術學院,長春 130026

交錯網格有限差分算法以其高效、精確、實用等優點在地震波數值模擬中得到廣泛應用。目前交錯網格有限差分的精度已達到時間4階、空間2N階；然而在求空間三次導數時，差分格式實際上并未達到所謂的2N階精度，而是采用了低階的差分格式及差分系數,這樣有利于提高大尺度空間正演時的計算效率；但從計算精度的角度考慮，有必要推導出準確的滿足2N階精度的交錯網格有限差分格式及差分系數，以得到更高精度的正演結果。筆者利用Taylor公式展開首次推導出了可導函數任意次導數的任意偶數階精度的差分近似式及相應的差分系數，從而完善了常規高精度交錯網格有限差分算法。采用新推導的交錯網格有限差分格式得到的正演波形與解析解進行了對比，證明了新推導的差分格式的正確性，并與常規差分格式的正演波形進行了比較，結果顯示,新推導出的交錯網格有限差分格式模擬結果穩定性好，精度更高。

交錯網格；差分格式；差分系數；高精度

0 引言

隨著我國石油天然氣勘探開發工作的不斷發展，研究人員面臨的勘探對象和開發條件越來越復雜，越來越困難。尋找復雜構造油氣藏、巖性油氣藏和裂縫油氣藏等“剩余油”的任務艱巨[1-2]。為了解決這些復雜油氣藏的勘探、開發問題，需要對地下復雜介質的地震響應進行高精度數值模擬研究。對地震彈性波方程進行數值模擬的方法主要包括有限差分法、有限元法和偽譜法[3]，而有限差分法由于計算速度快、精度高等特點使用最為廣泛[4-5]。

自從Alterman和Karal首次將有限差分法應用到各向同性介質彈性波的模擬中后，由于其自身的優點，很快被用于各種地震勘探學的數值問題上，且在應用中不斷發展[6],先后出現了變網格有限差分[7-8]、非連續網格有限差分[9]、不規則網格有限差分[10]、交錯網格有限差分[11-13]、旋轉交錯網格有限差分[14-15]、可變時間步長有限差分[16]、自適應可變空間步長網格有限差分以及隱式有限差分[17-18]。其中，交錯網格由Madariage[11]最早提出，Virieux[12]首先將其使用到一階速度-應力方程中。

為了避免對彈性常數進行空間微分，在彈性波正演模擬時，采用一階速度-應力彈性波方程。這樣可以在不同的時間層上使用不同的網格，分別進行應力和速度的計算和傳播。Virieux發展的交錯網格精度為O(Δt2+Δx2)(時間2階、空間2階)，與常規網格相比，在沒有增加計算量和存儲空間的情況下，局部精度提高了4倍，收斂速度也有所提高。Levander[19]又將交錯網格有限差分的精度提高到O(Δt2+Δx4)。隨后董良國等[20-21]發展了更高精度的交錯網格有限方法，精度達到了O(Δt4+Δx2N)。為了使用較少的時間層，不增加計算存儲空間，董良國等[20-21]將速度(應力)對時間的奇數次高階導數轉化為應力(速度)對空間的導數，把交錯網格和高階差分法成功地結合在一起。然而在高階交錯網格有限差分算法中，當變量對空間求三次導數時，采用的是低階的差分格式和差分系數[22-25]，沒有充分考慮“空間2N階精度”這一事實，這樣做雖然有助于提高有限差分的計算效率，但要想獲得更高精度的正演結果就有必要推導出準確的時間4階、空間2N階精度的交錯網格有限差分格式及相應的差分系數。

差分系數是決定交錯網格有限差分算法精度的關鍵，Taylor公式展開和最優化方法是求取差分系數的主要方式[1,17,26,28]。本文以Taylor公式展開為基礎，推導了可導函數任意次導數的任意偶數階精度的差分近似式以及差分系數，完善O(Δt4+Δx2N)一階速度-應力方程組差分格式。

1 常規交錯網格差分格式

在二維各向同性介質xoz平面內，假定體力為0，一階速度-應力彈性波方程為

(1-1)

(1-2)

(1-3)

(1-4)

(1-5)

其中：vx,vz為速度分量；τxx,τzz,τxz為應力分量；ρ為密度；c11，c13，c33，c44為介質的彈性常數。

1.1 時間2M階差分近似

(2)

式中，Δt為時間步長。令M=2，式(2)就是常規的時間4階精度差分近似。

為了減少計算內存，利用速度和應力的耦合關系，得到方程組式(1)的時間4階精度差分近似，以式(1-1)為例

(3)

1.2 空間2N階差分近似以及差分系數

在常規交錯網格算法中，速度(應力)分量的空間一次導數是由相錯半網格點的應力(速度)分量計算的。設函數f(x)具有2N-1階導數，令x=x0±[(2n-1)/2]Δx，則由f(x)在x處的2N-1階Taylor展開可以得到式(3)中空間一次導數的差分近似式

(4)

(5)

實心球表示函數值在半網格點上，空心球表示函數值在整網格點上。圖1 差分形式示意圖Fig.1 Schematic differential form

1.3 常規時間4階空間2N階差分格式

取x=iΔx，z=jΔz，t=kΔt，i、j和k分別表示空間和時間網格點，U、V分別代表速度分量vx、vz的離散值，R、T、H分別代表應力τxx、τzz、τxz的離散值。則方程(3)的精度為O(Δt4+Δx2N)的常規差分格式如下(Δx=Δz，下同)

(6)

可以明顯地看出，式(6)中P11、P12、P13、P14及P15的差分格式不是空間2N階精度的，而是用低階差分格式代替。

2 任意階導數有限差分近似式及差分系數推導

為了得到交錯網格有限差分的時間4階、空間2N階精度準確的差分格式及差分系數，首先需要推導出函數對空間二次、三次導數的差分格式及差分系數，然后將其代入到一階速度-應力波動方程組中，從而推導出該方程組精度為時間4階、空間2N階的準確差分格式。

2.1 函數任意次導數差分近似的差分系數

函數對空間任意次(奇數、偶數次)導數的差分近似式與空間一次導數的差分近似式形式相同，只是差分系數不同；所以推導函數任意次導數的準確的差分近似式的關鍵就是推導差分系數。交錯網格有限差分算法中，求離散化函數f在x=x0的導數時，假設x0為整網格點，差分形式有2種情況，一種是參與計算的函數值在半網格點上，另一種是參與計算的函數值在整網格點上，如圖1所示。

首先推導離散化函數值在半網格點上時，差分近似式的差分系數。令函數f(x)在x=x0±[(2n-1)/2]Δx兩處的Taylor展開式相減，得

(7)

寫成矩陣的形式為

(8)

其中

令

即X=B·D，所以F=D-1·B-1·S。若設F=C·S，則有

(9)

式中，C就是差分系數矩陣。由于D是初等矩陣，所以

則式(9)可寫為

(10)

式中，j=1,2,…,N。

為方便起見，只觀察式(10)中的第一個式子，兩邊取轉置有

整理有

即

(11)

通過求解式(11)可以得到變量一次導數的差分近似式的準確差分系數，其精度達到2N階。同理式(10)中其他的N-1個式子都可以得到類似式(11)的方程，即

(12)

(13)

再令函數f(x)在x=x0±[(2n-1)/2]Δx處的Taylor展開式相加，得

(14)

使用推導任意奇數次導數差分近似式和差分系數相同的方法，從式(14)可得任意偶數次導數差分近似式為

(15)

(16)

下面推導離散化函數值在整網格點上時，差分近似式的差分系數。令f(x)在x=x0±nΔx處進行Taylor展開，與函數值在半網格點上時的推導過程類似，該情況下的差分近似式及差分系數求解方程為

(17)

(18)

(19)

(20)

式(17)和式(19)分別表示離散化函數值在整網格點上時的任意奇數次和任意偶數次導數的差分近似式，其差分精度與半網格點上的一致。同時給出了差分系數的求解矩陣方程即(18)和(20)式。這樣函數f任意次導數的準確的差分近似式以及相應的差分系數全部推導完畢。

2.2 準確的時間4階、空間2N階差分格式

圖2 τxz對x一次偏導、對z二次偏導的差分示意圖Fig.2 τxz on x of a partial derivative and on z of secondary partial derivative

在確定了差分系數之后，將對應的變量對空間二次、三次導數差分近似式代入到一階速度-應力方程式(1-1)中，經過整理，得到其精度為O(Δt4+Δx2N)的準確差分格式：

(21)

方程(21)是一階速度-應力方程式(1-1)的精度為時間4階、空間2N階的通式，達到真正的空間2N階精度。其中的差分系數如表1所示。

3 穩定性以及邊界條件

通過傅里葉分析方法，得到二維各向同性介質xoz面內O(Δt4+Δx2N)精度的穩定性條件[21]：

式中，vP為縱波速度。

在模型試算時，筆者采用Collino[27]提出的完全匹配層(PML)吸收邊界條件對模型邊界進行處理，該完全匹配層(PML)能夠很好地吸收進入邊界的地震彈性波，消除了邊界反射對正演結果的影響。

表1 精度為O(Δt4+Δx2N)的差分系數

表2 常用差分精度的穩定性條件

Table 2 Stability conditions of common difference accuracy

差分精度WO(Δt4+Δx6)≤0.9068O(Δt4+Δx8)≤0.8753O(Δt4+Δx10)≤0.8553O(Δt4+Δx12)≤0.8410

4 模型試算

4.1 正確性驗證與精度對比

正演時采用各向同性均勻介質模型如圖3所示。物性參數為縱波速度vP=3 000 m/s，橫波速度vS=2 000 m/s，介質密度ρ=2 000 kg/m3。模型網格大小為200×200，震源坐標為(1 000 m，200 m)，空間步長為10 m，時間采樣間隔為0.001 s，用主頻為25 Hz的Richer子波激發。模型邊界采用PML吸收邊界。2個接收器分別位于A(800 m，600 m)和B(1 000 m，600 m)。

五角星表示震源位置，A、B點表示接收器位置。圖3 各向同性均勻介質模型Fig.3 Isotropic homogeneous model with its properties

圖4a和圖4b分別對應A點和B點的地震記錄，其中新推導出的交錯網格有限差分計算時采用時間4階、空間8階精度。從圖4中可以看出,無論在波初至還是在震相波形上都基本一致。從而證明了本文推導出的差分格式及差分系數是正確的。

a.A點；b.B點。圖4 新推導交錯網格差分格式的正演波形與解析解波形比較Fig.4 Comparisons between improved staggered-grid finite difference format solution waveform and analytic solutions waveform

為了進一步研究新推導的交錯網格有限差分的計算精度，對其與常規交錯網格有限差分正演波形分別和解析解波形的誤差進行了比較。圖5a和圖5b分別是A點和B點的地震記錄誤差。從圖5可以看出，新推導的交錯網格有限差分的正演波形誤差更小。因此，新推導的差分格式及差分系數比常規的差分格式及差分系數具有更高的精度。

a.A點；b.B點。圖5 常規交錯網格差分格式數值解誤差與新推導的差分格式數值解誤差對比Fig.5 Contrast between conventional staggered grid finite difference format solution error and improved difference format solutions error with analytic solution, respectively

4.2 水平層狀模型試算

為了驗證新推導出的差分格式對層狀介質中的正演效果，建立了圖6中的層狀地質模型，每層的物性參數如圖中所示。網格大小為250×250，空間步長為10 m，2個介質分界面分別在500 m和1 000 m處，震源坐標為(1 250 m，200 m)。在計算時采用新推導的時間4階、空間8階精度的差分格式，用主頻為25 Hz的Richer子波激發。

圖6 分層介質模型Fig.6 Homogeneous layered model with the layer properties

a．x分量；b．z分量。圖7 400 ms時模型波場快照Fig.7 Snapshots of Fig. 6 using the proposed format at time is 400 ms

圖7是精度為O(Δt4+Δx8)新推導的差分格式的400 ms波場快照，計算用時為207.9 s。從圖7中可以看出:新推導出的高精度交錯網格有限差分的差分格式能精確地模擬彈性波在地下層狀介質中的傳播特性；彈性波在地下500 m處和1 000 m處發生反射和透射，且出現轉換橫波。在其他計算條件完全一致的情況下，用常規的交錯網格有限差分公式模擬圖6中的模型計算用時也需要198.5 s，新推導的差分技術計算效率上只相差4.5%，但計算精度提高了，達到了真正的空間2N階，在本文4.1節已經做過詳細的分析對比。

5 結束語

常規高階交錯網格有限差分格式中出現變量對空間二次或三次導數時，采用了低階的差分格式，并沒有達到空間2N階精度。本文通過Taylor公式展開法，首次推導出了交錯網格有限差分算法中的變量任意次導數的任意偶數階精度的差分格式以及差分系數，完善了傳統時間4階、空間2N階精度的交錯網格有限差分算法，具有重要的理論意義。

[1] 牟永光，裴正林. 三維復雜介質地震數值模擬[M]. 北京：石油工業出版社，2005. Mu Yongguang,Pei Zhenglin.Seismic Numerical Simulation of 3D Complex Medium[M]. Beijing: Petroleum Industry Press, 2005.

[2] 李賓. 橫向各向同性介質有限差分法波場模擬方法研究[D]. 北京：中國石油大學，2009. Li Bin. Finite Difference Numerical Modeling of Seismic Wavefields in Tansversely Isotropic Media with a Vertical Axis[D]. Beijing: China University of Petroleum, 2009.

[3] 龍桂華，李小凡，張美根. 錯格傅里葉偽譜微分算子在波場模擬中的應用[J]. 地球物理學報，2009，52(1):193-199. Long Guihua, Li Xiaofan, Zhang Meigen. The Application of Staggered-Grid Fourier Pseudo Spectral Differentiation Operator in Wavefield Modeling[J]. Chinese Journal of Geophysics, 2009,52(1):193-199.

[4] 孟慶生，樊玉清，張珂,等.高階有限差分法管波傳播數值模擬[J].吉林大學學報：地球科學版，2011，41(1):292-298. Meng Qingsheng, Fan Yuqing, Zhang Ke, et al. Tube Wave Propagation Numerical Simulation Based on High Order Finite-Difference Method[J]. Journal of Jilin University: Earth Science Edition, 2011, 41(1):292-298.

[5] 周曉華，陳祖斌，曾曉獻，等.交錯網格有限差分法模擬微動信號[J].吉林大學學報：地球科學版，2012，42(3):852-857. Zhou Xiaohua, Chen Zubin, Zeng Xiaoxian, et al. Simulation of Microtremor Using Straggered-Grid Finite Difference Method[J]. Journal of Jilin University: Earth Science Edition, 2012, 42(3):852-857.

[6] 孫衛濤. 彈性波動方程的有限差分數值方法[M]. 北京：清華大學出版社，2009. Sun Weitao. The Finite Difference Numerical Method of Elastic Wave Equation[M]. Beijing: Tsinghua University Press,2009.

[7] Wang Y, Xu J, Gerard T S. Viscoelastic Wave Simulation in Basins by a Variable-Grid Finite-Difference Method[J]. Bull Seismol Soc Am, 2001, 91(6): 1741-1749.

[8] Narayan J P, Kumar S. A Fourth Order Accurate SH-Wave Staggered Grid Finite-Difference Algorithm with Variable Grid Size and VGR-Stress Imaging Technique[J]. Pure Appl Geophys, 2008, 165:271-294.

[9] Shin Aoi, Hiroyuki Fujiwara. 3D Finite-Difference Method Using Discontinuous Grids[J]. Bull Seismol Soc Am, 1999,89:918-930.

[10] Opr?al I, Zahradník J. Elastic Finite-Difference Method for Irregular Grids[J]. Geophysics, 1999,64:240-250.

[11] Madariaga R. Dynamics of an Expanding Circular Fault[J]. Bull Seismol Soc Am, 1976, 66: 639-666.

[12] Virieux J. SH-Wave Propagation in Heterogeneous Media: Velocity-Stress Finite-Difference Method[J]. Geophysics, 1984, 49(11): 1933-1957.

[13] Thomas B, Erik H S. Accuracy of Heterogeneous Staggered-Grid Finite-Difference Modeling of Rayleigh Waves[J]. Geophysics,2006,71(4): 109-115.

[14] Erik H S, Norbert G, Serge A S. Modeling the Propagation of Elastic Waves Using a Modified Finite-Difference Grid[J]. Wave Motion, 2000,31:77-92.

[15] Reeshidev B, Mrinal K S. Finite-Difference Modeling of S-Wave Splitting in Anisotropic Media[J]. Geophysical Prospecting, 2008, 56: 293-312.

[16] Tessmer E. Seismic Finite-Difference Modeling with Spatially Varying Time Steps[J]. Geophysics, 2000,65(4):1290-1293.

[17] Liu Y,Mrinal K S.A Practical Implicit Finite-Diffe-rence Method: Examples from Seismic Modeling[J]. Journal of Geophysics and Engineering, 2009, 6(3): 231-249.

[18] Liu Y, Mrinal K S. Finite-Difference Modeling with Adaptive Variable-Length Spatial Operators[J]. Geophysics, 2011, 76(4): 79-89.

[19] Levander A R. Fourth-Order Finite-Difference P-SV Seismograms[J]. Geophysics, 1988, 53(11):1425-1436.

[20] 董良國，馬在田，曹景中，等. 一階彈性波方程交錯網格高階差分解法[J]. 地球物理學報，2000，43(6)：856-864. Dong Liangguo, Ma Zaitian, Cao Jingzhong, et al. The Stability Study of the Staggered-Grid High-Order Difference Method of One-Order Elastic Equation[J].Chinese Journal of Geophysics, 2000,43(6):856-864.

[21] 董良國，馬在田，曹景中. 一階彈性波方程交錯網格高階差分解法穩定性研究[J]. 地球物理學報，2000，43(3):411-419. Dong Liangguo, Ma Zaitian, Cao Jingzhong. The Staggered-Grid High-Order Difference Method of One-Order Elastic Equation[J].Chinese Journal of Geophysics, 2000,43(3):411-419.

[22] Liu Y, Mrinal K S. Acoustic VTI Modeling with a Time-Space Domain Dispersion-Relation-Based Finite-Difference Scheme[J]. Geophysics,2010,75(3): 11-17.

[23] Hestholm S. Acoustic VTI Modeling Using High-Order Finite Differences[J]. Geophysics, 2009,74(5):67-73.

[24] Finkelstein B, Kastner R. Finite Difference Time Domain Dispersion Reduction Schemes[J]. Journal of Computational Physics, 2007,221: 422-438.

[25] Fei T W, Christopher L L. Hybrid Fourier Finite-Difference 3D Depth Migration for Anisotropic Media[J]. Geophysics, 2008, 73(2): 27-34.

[26] Fornberg B. Generation of Finite Difference Scheme on Arbitrarily Spaced Grids[J]. Mathematics of Computation, 1988, 51: 699-706.

[27] Collino F, Tsogka C. Application of the Perfectly Matched Absorbing Layer Model to the Linear Elasto Dynamic Problem in Anisotropic Hetergeneous Media[J]. Geophysics, 2001,66(1):294-307.

[28] Dablain M A. The Application of High-Order Diffe-rencing to the Scalar Wave Equation[J]. Geophysics,1986, 51:54-66.

《吉林大學學報(地球科學版)》征稿簡則

《吉林大學學報(地球科學版)》是以地學為特色的綜合性學術期刊，是中文核心期刊、中國科技核心期刊、中國科學引文數據庫(CSCD)重要期刊、中國學術期刊評價(2013-2014)(RCCSE)權威期刊。

本刊主要刊登地質與資源、地質工程與環境工程、地球探測與信息技術等學科領域中的最新科研成果。歡迎廣大作者踴躍投稿。

投稿注意事項如下：

1 來稿必須是未正式發表的、具有創新性成果的科技論文；必須是受到省部級及以上科研基金資助的立項課題的主要研究成果。

2 來稿要求論點明確、數據可靠、邏輯嚴密、文字精煉。每篇論文必須包含題目、作者姓名、作者單位(以上需中英文對照)及郵政編碼、中英文摘要和關鍵詞、中國圖書館分類號、正文、參考文獻及第一作者和通信作者簡介，并且在文稿首頁地腳處注明論文屬何項目(基金)資助及項目編號。作者簡介格式為：姓名(出生年-)、性別、民族(漢族可省略)、職稱、學位、研究方向、電話(可省略)、電子信箱。

3 論文摘要須寫成報道性文摘，包括目的、方法、結果、結論4部分(200～500字)，摘要應具有獨立性與自明性。摘要應由正文中重要語句組成。

4 文稿篇幅(含圖表)一般不超過1萬字，文中量和單位必須使用法定的量和計量單位，外文字母必須分清大小寫、正斜體、黑白體、上下角標字符的位置應區別明顯。

5 文中圖表應有自明性，圖表應附相應的中英文名。附圖力求簡明清晰，經緯度、比例尺、圖例等內容齊全。圖中文字(漢字用宋體，數字和字母用Times New Roman體；字號相當于Word文檔中的小5至6號字大小)、符號、縱橫坐標必須寫清，并與正文一致；坐標標值線應放在坐標軸內側，標目中的量和單位間用斜線(如：t/s，變量字母t用斜體)。凡涉及國界線的圖件必須繪制在地圖出版社公開出版的最新地理底圖上。圖版長×寬不超過23cm×16cm，照片要求清晰，層次分明。文稿中須留出插圖位置。

6 參考文獻只選最主要的列入，并要注意引用近兩年國內外的文獻。文后參考文獻表采用順序編碼制，按文中出現的先后順序編號。中文文獻要附上相應的英譯文獻，英譯時責任者的姓名首字母大寫，雙名連拼。英文文獻責任者姓名姓前名后：姓全拼，首字母大寫；名縮寫，不加縮寫點。英譯文獻及英文文獻的責任者不超過3人時全部寫出，超過時只寫前3名，后加“等”或“et al．”；題名、刊名實詞首字母均大寫。

不同文獻的著錄格式如下：

(1) 專著的著錄格式：[序號]主要責任者(作者、編者等)．題名(書名)[M]．版本(第一版不標注)．出版地：出版者，出版年．

(2) 論文集中析出文獻的著錄格式：[序號]作者．題名[C]//編者．論文集名．版本(第一版不標注)．會議地址：承辦單位(或出版地：出版者)，出版年：起止頁碼．

(3) 連續出版物文獻的著錄格式：[序號]作者．題名[J]．期刊名，年，卷(期)：起止頁碼．

(4) 專利的著錄格式：[序號]專利申請者．專利題名[P]．國別，專利文獻種類，專利號．出版日期．

(5) 學位論文的著錄格式：[序號]作者．題名[D]．保存地點：保存單位，年份．

(6) 技術標準的著錄格式：[序號]標準代號．標準順序-發布年．標準名稱[S]．出版地：出版者，出版年．

(7) 報告的著錄格式：[序號]作者. 題名[R]. 保存地：保存單位，保存年.

(8) 電子文獻的著錄格式：[序號]主要責任者(作者)．電子文獻題名[電子文獻及載體類型標識]．發表或更新日期/引用日期(任選)．電子文獻的出處或可獲得地址.

7 請勿一稿多投。投稿3個月后仍未接到反饋信息者，請與編輯部聯系。

8 本刊已加入《中國學術期刊(光盤版)》及其他數據庫，其作者文章著作權使用費與本刊稿酬一次性給付。如作者不同意文章被收錄，請在來稿時向本刊聲明，本刊將做適當處理。

請作者提供正文的電子文檔(方正、Word均可)，并提供矢量圖格式文件，圖件盡量使用CorelDraw繪制。

通信地址：吉林省長春市西民主大街938號《吉林大學學報(地球科學版)》編輯部。郵編130026。

電子信箱: xuebao1956@jlu.edu.cn 本刊網址：http://xuebao.jlu.edu.cn/dxb

電 話：0431-88502374 發行代號：12-22(國內)；BM5074(國外)

Staggered Grid Finite Difference Scheme and Coefficients Deduction of Any Number of Derivatives

Yang Qingjie, Liu Cai, Geng Meixia, Feng Xuan, Guo Zhiqi, Liu Yang

CollegeofGeoExplorationScienceandTechnology,JilinUniversity,Changchun130026,China

Staggered grid finite difference algorithm is effective, accurate and practical, so it has wide application prospect and practical significance.So far, the common high-order scheme of staggered grid finite difference algorithm is 4-order temporal and 2N-order spatial accuracy.However, the FD scheme isn’t 2N-order accuracy actually when computing second or third spatial derivative.The authors deduce the FD scheme with 2N-order accuracy and corresponding coefficients of any number of derivatives of functions which have any number of derivatives forthefirsttime.So we can consummate conventional high-order staggered grid finite difference algorithm.We make a simulation of seismic response with conventional FD scheme and the new FD scheme respectively，and compared with analytic solution respectively. As a result, the new FD scheme is more stabilized and more accurate.

staggered grid; finite difference scheme; finite difference coefficients; high-order

10.13278/j.cnki.jjuese.201401307.

2013-07-01

國家自然科學基金項目(40974054，41174080)；國家“973”計劃項目(2009CB219301)；油頁巖勘探開發利用產學研用合作創新研究項目(OSP-02，OSR-02)

楊慶節(1987-)，男，博士研究生，主要從事多波多分量地震勘探研究，E-mail:qjyang58@gmail.com

劉財(1963-)，男，教授，博士生導師，主要從事地震波場正反演理論、綜合地球物理等研究，E-mail:liucai@jlu.edu.cn。

10.13278/j.cnki.jjuese.201401307

P631.4

A

楊慶節，劉財，耿美霞，等.交錯網格任意階導數有限差分格式及差分系數推導.吉林大學學報:地球科學版,2014,44(1):375-385.

Yang Qingjie, Liu Cai, Geng Meixia,et al.Staggered Grid Finite Difference Scheme and Coefficients Deduction of Any Number of Derivatives.Journal of Jilin University:Earth Science Edition,2014,44(1):375-385.doi:10.13278/j.cnki.jjuese.201401307.