蔡小蘇 周文宗

摘 要：地震波數(shù)值模擬在地震勘探中具有重要作用。有限差分法憑借自身優(yōu)點，廣泛的應用于波動方程求解中。本文將通量校正方法和交錯網(wǎng)格高階差分法有效的結合，對均勻半空間介質二維一階速度-應力波動方程組進行數(shù)值解。得到同一時刻，時間二階，不同空間階數(shù)，不同網(wǎng)格步長的波場快照。結果表明，細網(wǎng)格步長和高階差分能較好地壓制數(shù)值頻散;通量校正法在壓制網(wǎng)格數(shù)值頻散方面有明顯的優(yōu)勢，在保持一定計算精度情況下，可采用較大的網(wǎng)格間距，提高計算效率。

關鍵詞：有限差分法;數(shù)值模擬;交錯網(wǎng)格;頻散;通量校正傳輸法

中圖分類號：P631.4 文獻標識碼：A 文章編號：1671-2064（2020）08-0218-03

0 引言

地震波數(shù)值模擬是研究各種地質條件下構造、物性和巖性等各種地質因素與地震波響應特征（運動學和動力學特征）之間關系的一門技術。地震波數(shù)值模擬是地震反演和偏移成像的基礎，對于人們理解地震波傳播規(guī)律，進行實際地震資料的地質解釋以及地球資源開發(fā)等，均具有極其重要的意義。地震波波場模擬是以彈性波理論為基礎，彈性波方程是其基本方程[1]。有限差分法憑借應用范圍廣、方便靈活、計算效率高等優(yōu)點，廣泛用于數(shù)值模擬中。然而，用有限的離散網(wǎng)格來近似連續(xù)的無限介質，以差商代替微商，會出現(xiàn)不同頻率的地震波在介質中以不同速度傳播的現(xiàn)象稱為波的頻散[2]。數(shù)值頻散嚴重的影響波動方程數(shù)值計算精度，必須給予消除。本文將通量校正方法（Flux corrected transport method，F(xiàn)CT）應用于二維一階速度一應力彈性波動方程的交錯網(wǎng)格高階有限差分求解中，進行數(shù)值模擬試驗，并與其他差分模擬結果進行對比分析。

1 二維一階應力-速度彈性波方程

由運動平衡微分方程、本構方程和幾何方程推導得到二維一階速度-應力彈性波動方程[3-4]：

表示介質的密度;t表示時間;、表示介質的拉梅常數(shù)和剪切模量;vx、vz表示質點在x，z方向的速度;和是正應力分量;是切應力分量。

2 二維交錯網(wǎng)格有限差分近似

與其他地震波數(shù)值模擬方法相比，有限差分法具有適應條件廣、方便靈活、便于實現(xiàn)等優(yōu)點，因此大量應用于數(shù)值模擬中。而其中的交錯網(wǎng)格高階差分法與常規(guī)網(wǎng)格高階差分相比，可進一步提高數(shù)值模擬的精度并有效的壓制數(shù)值頻散[5]。

2.1 時間2M階差分近似

在用交錯網(wǎng)格有限差分法解一階彈性波方程時[6]，在半節(jié)點處（）計算速度分量vx、vz，在節(jié)點處（t時刻）計算應力分量、、。利用Taylor公式將和在t時刻展開，可得到：

2.2 空間2N階差分近似

在交錯網(wǎng)格有限差分計算中，在相應的變量網(wǎng)格半節(jié)點上計算變量的導數(shù)。對于具有2N+1階導數(shù)的連續(xù)函數(shù)f（x），其一階空間導數(shù)如下：

3 通量校正傳輸方法（FCT）

由于有限差分是對時間和空間網(wǎng)格的離散，以差商代替微商，以有限近似連續(xù)介質，這樣就出現(xiàn)了模擬中不期望出現(xiàn)的波動現(xiàn)象—數(shù)值頻散。差分格式的截斷性決定了波動方程計算中的數(shù)值頻散問題是無法避免的。數(shù)值頻散嚴重干擾數(shù)值計算精度，對模擬的結果造成嚴重影響，應盡量壓制。

Boris和Book[7]等最早在求解流體動力學連續(xù)方程時提出了通量校正傳輸方法。而后將FCT法用于求解聲波方程，有效的壓制了在粗網(wǎng)格情況下的差分計算產(chǎn)生的數(shù)值頻散。

4 模型分析

設計尺寸為400m×400m均勻半空間介質模型，震源位于（200，200）處，震源為主頻30Hz的Ricker子波，時間間隔，縱波速度，橫波速度，密度為。通過交錯網(wǎng)格差分得到時間二階，不同空間階數(shù)，不同網(wǎng)格步長模型在100ms時垂直分量的波場快照圖，圖1～圖10所示[8-11]。

對比圖1、圖3、圖5，當增加空間步長，即減少一個波長內離散點數(shù)，數(shù)值頻散越嚴重。對比圖1、圖7、圖9，隨著空間差分階數(shù)的提高，產(chǎn)生的數(shù)值頻散會逐漸減小，因此可采用高階差分來減小數(shù)值頻散。對比圖1和圖2、圖7和圖8、圖9和圖10，經(jīng)過FCT校正后，數(shù)值頻散得到了很好的壓制。對比圖2、圖9，發(fā)現(xiàn)低階差分的波場模擬通過FCT校正也可以得到較高精度的模擬效果。

5 結論

本文給出了彈性波二維一階應力-速度波動方程在交錯網(wǎng)格中時間和空間差分近似，介紹了FCT方法，并將其用于波場模擬中。實際模型結果表明，采用細網(wǎng)格、高階差分能較好的壓制頻散;通量校正方法與交錯網(wǎng)格高階差分有效結合，在壓制網(wǎng)格數(shù)值頻方面有明顯的效果，并可采用較大空間步長，提高計算效率。

參考文獻

[1] 姚姚.地震波場與地震勘探[M].北京：地質出版社，2006.

[2] 董良國，李培明.地震波傳播數(shù)值模擬中的頻散問題[J].天然氣工業(yè)，2004（6）：53-65.

[3] 揚頂輝，藤吉文.各向異性介質中三分量地震記錄的FCT有限差分模擬[J].石油地球物理勘探，1997（2）：18l-188.

[4] 董良國，馬在田，曹景忠.一階彈性波方程交錯網(wǎng)格高階差分解法[J].地球物理學報，2000（3）：412-419.

[5] 蔡其新，何佩軍.有限差分數(shù)值模擬的最小頻率算法及其應用[J].石油地球物理勘探，2003（03）：247-251.

[6] 吳國忱，王華忠.波場模擬中的數(shù)值頻散分析與校正策略[J].地球物理學進展，2005（1）：58-65.

[7] Book D L， Boris J P， Hain K. Flux-corrected transport I：SHASTA， A fluid transport algorithm that works[J].Journal of Computational Physics，1973，11（1）：38-69.

[8] 鄭海山，張中杰.橫向各向同性（VTI）介質中非線性地震波場數(shù)值模擬[J].地球物理學報，2005（3）：660-671.

[9] 李勝軍，孫成禹，高建虎，等.地震波數(shù)值模擬中的頻散壓制方法分析[J].石油物探，2008（5）：444-449.

[10] 潘海濱，孫萍.交錯網(wǎng)格波場數(shù)值模擬頻散分析與校正策略[J].海洋動態(tài)地質，2009（6）：36-41.

[11] 李文杰，張改蘭.利用FCT方法壓制彈性波數(shù)值模擬中的數(shù)值頻散[J].物探化探計算技術，2011（3）：248-251.