本發(fā)明屬于石油勘探領(lǐng)域，具體涉及一種劇烈起伏自由地表有限差分正演模擬系統(tǒng)及方法。

背景技術(shù)：

劇烈的起伏地表廣泛地存在于中國(guó)地區(qū)。地震勘探的前緣正逐步向復(fù)雜地表地區(qū)轉(zhuǎn)移，常規(guī)的處理和解釋方法對(duì)上述復(fù)雜起伏地表已經(jīng)不再適用。為了更好地解決勘探開(kāi)發(fā)難題，須密切研究起伏地表?xiàng)l件下介質(zhì)地震波正演模擬，對(duì)于起伏地表?xiàng)l件下的地質(zhì)構(gòu)造進(jìn)行正演模擬，可從數(shù)值模擬結(jié)果中認(rèn)識(shí)、分析地震波波場(chǎng)特征及地震波傳播規(guī)律等，更好地指導(dǎo)后續(xù)的地震勘探。

有限差分方法因其計(jì)算速度快、占用內(nèi)存低、容易實(shí)現(xiàn)等優(yōu)點(diǎn)得到廣泛應(yīng)用，但有限差分方法是基于矩形網(wǎng)格剖分的，因此在處理起伏地表時(shí)存在很大的缺陷，在模擬地震波傳播時(shí)會(huì)產(chǎn)生明顯的繞射噪音及不穩(wěn)定。一種常用的解決方法是對(duì)物理模型進(jìn)行曲網(wǎng)格剖分，并將物理域的曲網(wǎng)格變換為矩形網(wǎng)格，但該方法也存在一定的局限性，該方法對(duì)穩(wěn)定性的要求較高，無(wú)法處理劇烈的起伏地表。

技術(shù)實(shí)現(xiàn)要素：

針對(duì)現(xiàn)有技術(shù)中存在的上述技術(shù)問(wèn)題，本發(fā)明提出了一種劇烈起伏自由地表有限差分正演模擬系統(tǒng)及方法，設(shè)計(jì)合理，克服了現(xiàn)有技術(shù)的不足，具有良好的效果。

為了實(shí)現(xiàn)上述目的，本發(fā)明采用如下技術(shù)方案：

一種劇烈起伏自由地表有限差分正演模擬系統(tǒng)，包括模型輸入模塊、正演模擬模塊、速度場(chǎng)轉(zhuǎn)換模塊、波場(chǎng)計(jì)算模塊、波場(chǎng)存儲(chǔ)模塊、波場(chǎng)反變換模塊以及結(jié)果輸出模塊；

模型輸入模塊，被配置為用于輸入速度場(chǎng)和起伏地表的高程；

正演模擬模塊，被配置為用于進(jìn)行劇烈起伏自由地表有限差分正演模擬；

速度場(chǎng)轉(zhuǎn)換模塊，被配置為用于將速度場(chǎng)變換到矩形網(wǎng)格坐標(biāo)系下；

波場(chǎng)計(jì)算模塊，被配置為用于在過(guò)渡坐標(biāo)系下更新地下區(qū)域波場(chǎng)以及在過(guò)渡坐標(biāo)系下采用自由地表邊界條件更新地表波場(chǎng)；

波場(chǎng)存儲(chǔ)模塊，被配置為用于在矩形網(wǎng)格坐標(biāo)系下存儲(chǔ)波場(chǎng)；

波場(chǎng)反變換模塊，被配置為用于將波場(chǎng)變換到笛卡爾坐標(biāo)系下；

結(jié)果輸出模塊，被配置為用于輸出最終的波場(chǎng)快照和炮記錄。

此外，本發(fā)明還提到一種劇烈起伏自由地表有限差分正演模擬方法，該方法采用如上所示的一種劇烈起伏自由地表有限差分正演模擬系統(tǒng)，包括如下步驟：

步驟1：通過(guò)模型輸入模塊輸入速度場(chǎng)和起伏地表的高程，建立觀測(cè)系統(tǒng)；

步驟2：通過(guò)速度場(chǎng)轉(zhuǎn)換模塊采用如下所示的變換方程將速度場(chǎng)變換到矩形網(wǎng)格坐標(biāo)系下；

其中，x和z為笛卡爾坐標(biāo)下的空間坐標(biāo)，ξ和η為矩形網(wǎng)格坐標(biāo)系下的空間坐標(biāo)，z₀(ξ)為笛卡爾坐標(biāo)系下的高程函數(shù)，η_max為矩形網(wǎng)格坐標(biāo)系下的最大采樣點(diǎn)數(shù)；

步驟3：引入過(guò)渡坐標(biāo)系，通過(guò)波場(chǎng)計(jì)算模塊在過(guò)渡坐標(biāo)系下采用如下所示的一階速度-應(yīng)力彈性波方程更新地下區(qū)域波場(chǎng)；

其中，u_x和u_z分別為水平分量和垂直分量的速度，τ_xx和τ_zz是正應(yīng)力，τ_xz是切應(yīng)力，t是時(shí)間；λ和μ是拉梅常數(shù)，ρ是密度；系數(shù)A和B由下式求得：

步驟4：通過(guò)波場(chǎng)計(jì)算模塊在過(guò)渡坐標(biāo)系下采用如下所示的自由地表邊界條件更新地表波場(chǎng)；

其中，C-J八個(gè)系數(shù)可由下式求得：

步驟5：通過(guò)波場(chǎng)存儲(chǔ)模塊在矩形網(wǎng)格坐標(biāo)系下存儲(chǔ)區(qū)域波場(chǎng)和地表波場(chǎng)；

步驟6：通過(guò)波場(chǎng)反變換模塊將區(qū)域波場(chǎng)和地表波場(chǎng)變換到笛卡爾坐標(biāo)系下；

步驟7：通過(guò)結(jié)果輸出模塊輸出最終的波場(chǎng)快照和炮記錄。

優(yōu)選地，對(duì)步驟3中過(guò)渡坐標(biāo)系下的一階速度-應(yīng)力彈性波方程的推導(dǎo)過(guò)程如下：笛卡爾坐標(biāo)系下的一階速度-應(yīng)力彈性波方程為：

應(yīng)用鏈?zhǔn)椒▌t，(2)式可變換為

對(duì)于方程(3)，我們需要計(jì)算和可由下式求得

方程(4)變換得到方程(5)

由映射方程

可得

將方程(7)代入方程(5)可得

將方程(8)帶入方程(3)，得到過(guò)渡坐標(biāo)系下的一階速度-應(yīng)力方程如方程(9)所示。

優(yōu)選地，對(duì)步驟4中過(guò)渡坐標(biāo)系下的自由地表邊界條件的推導(dǎo)過(guò)程如下：

笛卡爾坐標(biāo)系下的自由地表邊界條件為：

其中

將方程(11)代入到方程(10)中可得

對(duì)方程(12)應(yīng)用鏈?zhǔn)椒▌t可得

在起伏地表處

將方程(14)代入方程(13)可得

整理方程(15)可得

本發(fā)明所帶來(lái)的有益技術(shù)效果：

本發(fā)明能夠?qū)×移鸱乇磉M(jìn)行準(zhǔn)確模擬的有限差分正演模擬，通過(guò)引入過(guò)渡坐標(biāo)系和矩形網(wǎng)格坐標(biāo)系，在三個(gè)坐標(biāo)系之間進(jìn)行轉(zhuǎn)換，既能克服傳統(tǒng)傳統(tǒng)有限差分方法在處理起伏地表時(shí)的缺陷，又能保持有限差分方法的精度，而不明顯增加計(jì)算量，開(kāi)發(fā)劇烈起伏自由地表有限差分正演模擬方法，為研究劇烈起伏地表地區(qū)的波場(chǎng)傳播規(guī)律提供了理論基礎(chǔ)，同時(shí)，為劇烈起伏地表地區(qū)的采集數(shù)據(jù)進(jìn)行后續(xù)成像和反演工作提供了高精度的地震波正演基礎(chǔ)。

附圖說(shuō)明

圖1為本發(fā)明一種劇烈起伏自由地表有限差分正演模擬方法的流程圖。

圖2為三個(gè)坐標(biāo)系網(wǎng)格示意圖，其中，圖(a)為笛卡爾坐系網(wǎng)格示意圖；圖(b)為過(guò)渡坐標(biāo)系網(wǎng)格示意圖；圖(c)為矩形網(wǎng)格坐標(biāo)系網(wǎng)格示意圖。

圖3為加拿大逆掩斷層速度模型示意圖，其中，圖(a)為縱波速度模型示意圖；圖(b)為橫波速度模型示意圖。

圖4為起伏地表高程函數(shù)示意圖。

圖5為加拿大逆掩斷層網(wǎng)格剖分圖，其中，圖(a)為笛卡爾坐標(biāo)系網(wǎng)格剖分圖；圖(b)為過(guò)渡坐標(biāo)系網(wǎng)格剖分圖。

圖6為采用本發(fā)明方法得到的波場(chǎng)快照?qǐng)D，其中，圖(a)為0.625s水平分量的波場(chǎng)快照?qǐng)D；圖(b)為0.625s垂直分量的波場(chǎng)快照?qǐng)D；圖(c)為0.7s水平分量的波場(chǎng)快照?qǐng)D；圖(d)為0.7s垂直分量的波場(chǎng)快照?qǐng)D。

圖7為采用傳統(tǒng)方法得到的波場(chǎng)快照?qǐng)D，其中，圖(a)為0.625s水平分量的波場(chǎng)快照?qǐng)D；圖(b)為0.625s垂直分量的波場(chǎng)快照?qǐng)D；圖(c)為0.7s水平分的波場(chǎng)快照?qǐng)D量；圖(d)為0.7s垂直分量的波場(chǎng)快照?qǐng)D。

圖8為采用本發(fā)明方法得到的炮記錄。

圖9為本發(fā)明一種劇烈起伏自由地表有限差分正演模擬系統(tǒng)的結(jié)構(gòu)示意圖。

具體實(shí)施方式

下面結(jié)合附圖以及具體實(shí)施方式對(duì)本發(fā)明作進(jìn)一步詳細(xì)說(shuō)明：

實(shí)施例1

如圖9所示，一種劇烈起伏自由地表有限差分正演模擬系統(tǒng)，包括模型輸入模塊、正演模擬模塊、速度場(chǎng)轉(zhuǎn)換模塊、波場(chǎng)計(jì)算模塊、波場(chǎng)存儲(chǔ)模塊、波場(chǎng)反變換模塊以及結(jié)果輸出模塊；

模型輸入模塊，被配置為用于輸入速度場(chǎng)和起伏地表的高程；

正演模擬模塊，被配置為用于進(jìn)行劇烈起伏自由地表有限差分正演模擬；

速度場(chǎng)轉(zhuǎn)換模塊，被配置為用于將速度場(chǎng)變換到矩形網(wǎng)格坐標(biāo)系下；

波場(chǎng)計(jì)算模塊，被配置為用于在過(guò)渡坐標(biāo)系下更新地下區(qū)域波場(chǎng)以及在過(guò)渡坐標(biāo)系下采用自由地表邊界條件更新地表波場(chǎng)；

波場(chǎng)存儲(chǔ)模塊，被配置為用于在矩形網(wǎng)格坐標(biāo)系下存儲(chǔ)波場(chǎng)；

波場(chǎng)反變換模塊，被配置為用于將波場(chǎng)變換到笛卡爾坐標(biāo)系下；

結(jié)果輸出模塊，被配置為用于輸出最終的波場(chǎng)快照和炮記錄。

實(shí)施例2

在上述實(shí)施例的基礎(chǔ)上，本發(fā)明還提到一種劇烈起伏自由地表有限差分正演模擬方法，其流程圖如圖1所示，包括如下步驟：

步驟1：通過(guò)模型輸入模塊輸入速度場(chǎng)和起伏地表的高程，建立觀測(cè)系統(tǒng)；

步驟2：通過(guò)速度場(chǎng)轉(zhuǎn)換模塊采用如下所示的變換方程將速度場(chǎng)變換到矩形網(wǎng)格坐標(biāo)系下；

其中，x和z為笛卡爾坐標(biāo)下的空間坐標(biāo)，ξ和η為矩形網(wǎng)格坐標(biāo)系下的空間坐標(biāo)，z₀(ξ)為笛卡爾坐標(biāo)系下的高程函數(shù)，η_max為矩形網(wǎng)格坐標(biāo)系下的最大采樣點(diǎn)數(shù)；

步驟3：引入過(guò)渡坐標(biāo)系，通過(guò)波場(chǎng)計(jì)算模塊在過(guò)渡坐標(biāo)系下采用如下所示的一階速度-應(yīng)力彈性波方程更新地下區(qū)域波場(chǎng)；

其中，u_x和u_z分別為水平分量和垂直分量的速度，τ_xx和τ_zz是正應(yīng)力，τ_xz是切應(yīng)力，t是時(shí)間；λ和μ是拉梅常數(shù)，ρ是密度；系數(shù)A和B由下式求得：

步驟4：通過(guò)波場(chǎng)計(jì)算模塊在過(guò)渡坐標(biāo)系下采用如下所示的自由地表邊界條件更新地表波場(chǎng)；

其中，C-J八個(gè)系數(shù)可由下式求得：

步驟5：通過(guò)波場(chǎng)存儲(chǔ)模塊在矩形網(wǎng)格坐標(biāo)系下存儲(chǔ)區(qū)域波場(chǎng)和地表波場(chǎng)；

步驟6：通過(guò)波場(chǎng)反變換模塊將區(qū)域波場(chǎng)和地表波場(chǎng)變換到笛卡爾坐標(biāo)系下；

步驟7：通過(guò)結(jié)果輸出模塊輸出最終的波場(chǎng)快照和炮記錄。

優(yōu)選地，對(duì)步驟3中過(guò)渡坐標(biāo)系下的一階速度-應(yīng)力彈性波方程的推導(dǎo)過(guò)程如下：

笛卡爾坐標(biāo)系下的一階速度-應(yīng)力彈性波方程為：

應(yīng)用鏈?zhǔn)椒▌t，(2)式可變換為

對(duì)于方程(3)，我們需要計(jì)算和可由下式求得

方程(4)變換得到方程(5)

由映射方程

可得

將方程(7)代入方程(5)可得

將方程(8)帶入方程(3)，得到過(guò)渡坐標(biāo)系下的一階速度-應(yīng)力方程如方程(9)所示。

優(yōu)選地，對(duì)步驟4中過(guò)渡坐標(biāo)系下的自由地表邊界條件的推導(dǎo)過(guò)程如下：

笛卡爾坐標(biāo)系下的自由地表邊界條件為：

其中

將方程(11)代入到方程(10)中可得

對(duì)方程(12)應(yīng)用鏈?zhǔn)椒▌t可得

在起伏地表處

將方程(14)代入方程(13)可得

整理方程(15)可得

圖2為三個(gè)坐標(biāo)系網(wǎng)格示意圖，圖(a)笛卡爾坐系網(wǎng)格示意圖；圖(b)過(guò)渡坐標(biāo)系網(wǎng)格示意圖；圖(c)矩形網(wǎng)格坐標(biāo)系網(wǎng)格示意圖；

本發(fā)明一種劇烈起伏自由地表有限差分正演模擬方法，應(yīng)用于國(guó)際標(biāo)準(zhǔn)的起伏地表加拿大逆掩斷層模型(如圖3所示)，取得了理想的計(jì)算效果。輸入速度場(chǎng)和起伏地表的高程(如圖4所示)，建立觀測(cè)系統(tǒng)；將速度場(chǎng)變換到矩形網(wǎng)格坐標(biāo)系下(如圖5b所示)；引入過(guò)渡坐標(biāo)系，在過(guò)渡坐標(biāo)系下(如圖5a所示)更新地下區(qū)域波場(chǎng)；在過(guò)渡坐標(biāo)系下采用自由地表邊界條件更新地表波場(chǎng)(如圖6所示)；在矩形網(wǎng)格坐標(biāo)系下存儲(chǔ)波場(chǎng)；將波場(chǎng)變換到笛卡爾坐標(biāo)系下；輸出最終的波場(chǎng)快照和炮記錄(如圖8所示)。本發(fā)明的劇烈起伏自由地表有限差分正演模擬方法所得的波場(chǎng)快照(如圖6所示)，波形清晰，沒(méi)有出現(xiàn)繞射噪音和不穩(wěn)定現(xiàn)象。而傳統(tǒng)有限差分方法得到的波場(chǎng)快照(如圖7所示)出現(xiàn)了明顯的不穩(wěn)定。

本發(fā)明能夠?qū)×移鸱乇磉M(jìn)行準(zhǔn)確模擬的有限差分正演模擬，通過(guò)引入過(guò)渡坐標(biāo)系和矩形網(wǎng)格坐標(biāo)系，在三個(gè)坐標(biāo)系之間進(jìn)行轉(zhuǎn)換，既能克服傳統(tǒng)傳統(tǒng)有限差分方法在處理起伏地表時(shí)的缺陷，又能保持有限差分方法的精度，而不明顯增加計(jì)算量，開(kāi)發(fā)劇烈起伏自由地表有限差分正演模擬方法，為研究劇烈起伏地表地區(qū)的波場(chǎng)傳播規(guī)律提供了理論基礎(chǔ)，同時(shí)，為劇烈起伏地表地區(qū)的采集數(shù)據(jù)進(jìn)行后續(xù)成像和反演工作提供了高精度的地震波正演基礎(chǔ)。

當(dāng)然，上述說(shuō)明并非是對(duì)本發(fā)明的限制，本發(fā)明也并不僅限于上述舉例，本技術(shù)領(lǐng)域的技術(shù)人員在本發(fā)明的實(shí)質(zhì)范圍內(nèi)所做出的變化、改型、添加或替換，也應(yīng)屬于本發(fā)明的保護(hù)范圍。