一種實(shí)現(xiàn)地震波場(chǎng)數(shù)值模擬方法
【專利摘要】本發(fā)明涉及地球物理勘探【技術(shù)領(lǐng)域】，特別是關(guān)于地震波場(chǎng)數(shù)值模擬【技術(shù)領(lǐng)域】，具體的講是一種地震波場(chǎng)數(shù)值模擬的方法，包括在波數(shù)域使用低階有理分式有限差分算子逼近超高階顯式有限差分算子；計(jì)算所述低階有理分式有限差分算子的系數(shù)；將有理分式有限差分算子進(jìn)行分解，利用傅里葉變換將低階有理分式有限差分算子從波數(shù)域變換到空間域；在空間域使用有理分式有限差分法分別計(jì)算不同方向上的空間二階偏導(dǎo)數(shù)項(xiàng)；在時(shí)間方向上進(jìn)行迭代實(shí)現(xiàn)地震波場(chǎng)的數(shù)值模擬。通過本發(fā)明實(shí)施例的方法，使用低階有理分式有限差分法來實(shí)現(xiàn)地震波數(shù)值模擬，提高了地震波場(chǎng)數(shù)值模擬的精度和計(jì)算效率。
【專利說明】一種實(shí)現(xiàn)地震波場(chǎng)數(shù)值模擬方法
【技術(shù)領(lǐng)域】
[0001]本發(fā)明涉及地球物理勘探【技術(shù)領(lǐng)域】，特別是關(guān)于地震波場(chǎng)數(shù)值模擬【技術(shù)領(lǐng)域】，具體的講是一種實(shí)現(xiàn)地震波場(chǎng)數(shù)值模擬方法。
【背景技術(shù)】
[0002]地震波數(shù)值模擬或者地震波波動(dòng)方程正演是地球物理勘探領(lǐng)域的重要課題，在研究地下復(fù)雜介質(zhì)的地震響應(yīng)、復(fù)雜介質(zhì)的高精度地震成像(例如逆時(shí)偏移)、高精度的地下介質(zhì)速度建模(例如全波形反演)等領(lǐng)域有著廣泛的應(yīng)用。現(xiàn)有技術(shù)中進(jìn)行地震波數(shù)值模擬時(shí)主要采用如下幾種方法:
[0003]1、有限元法
[0004]基于變分原理和網(wǎng)格插值的有限元法具有靈活的不規(guī)則網(wǎng)格剖分方案，能從幾何形狀上適應(yīng)不規(guī)則地下介質(zhì)或起伏的地表，但是算法復(fù)雜、計(jì)算存儲(chǔ)量大、計(jì)算精度較低和計(jì)算效率低下，不適合于大規(guī)模的地震波場(chǎng)數(shù)值模擬。
[0005]2、偽譜法
[0006]利用傅里葉正變換將空間偏導(dǎo)數(shù)的計(jì)算轉(zhuǎn)換到波數(shù)域進(jìn)行，再通過傅里葉反變換得到空間域的空間偏導(dǎo)數(shù)值。理論上，偽譜法計(jì)算空間偏導(dǎo)數(shù)具有最高的數(shù)值精度。一般來說，偽譜法的計(jì)算效率要低于有限差分法。
[0007]3、譜元法
[0008]將有限元法和偽譜法結(jié)合就產(chǎn)生了譜元法，它充分利用了有限元法的靈活網(wǎng)格剖分特點(diǎn)和偽譜法數(shù)值模擬精度高的優(yōu)點(diǎn)，因此在求解波動(dòng)方程時(shí)，譜元法所用的網(wǎng)格節(jié)點(diǎn)少、數(shù)值精度高。但是，譜元法涉及到矩陣求逆，計(jì)算效率仍然非常低下，比較適合于小規(guī)模和中等規(guī)模的波動(dòng)方程數(shù)值模擬。
[0009]4、有限差分法
[0010]有限差分法是波動(dòng)方程數(shù)值模擬中實(shí)際應(yīng)用中使用最多的方法，算法相對(duì)比較簡單、與前面三種方法相比計(jì)算效率較高，適合于大規(guī)模復(fù)雜介質(zhì)的地震波場(chǎng)數(shù)值模擬。常見的有限差分法一般分為顯式有限差分法和隱式有限差分法，隱式有限差分法一般在頻率域進(jìn)行，常常用于波動(dòng)方程反演(如全波形反演)，由于只需要少數(shù)幾個(gè)頻率采樣點(diǎn)，因此計(jì)算效率相對(duì)比較高，但是對(duì)于三維問題往往帶來巨大的存儲(chǔ)量和計(jì)算量；顯式差分法一般在時(shí)間-空間域進(jìn)行，針對(duì)所有頻率采樣點(diǎn)同時(shí)計(jì)算，常常用于波動(dòng)方程正演、波動(dòng)方程成像和波動(dòng)方程反演。
[0011]對(duì)于高階顯式有限差分法，隨著差分階數(shù)的增加，計(jì)算精度也隨之增加。理論的空間二階偏導(dǎo)數(shù)對(duì)應(yīng)的差分算子響應(yīng)曲線為一拋物線，高階顯式差分算子響應(yīng)曲線對(duì)拋物線的逼近程度越高，實(shí)現(xiàn)數(shù)值模擬時(shí)的計(jì)算精度就越高。附圖3為2M階差分精度的高階顯式差分算子在波數(shù)域的響應(yīng)，隨著M值的增大，對(duì)應(yīng)的曲線越來越逼近理論的拋物線。
[0012]但是，隨著差分階數(shù)的增加，高階顯式有限差分的計(jì)算效率也隨著降低。高階顯式有限差分法的計(jì)算效率與差分階數(shù)呈反比例關(guān)系。以X坐標(biāo)軸方向的空間偏導(dǎo)數(shù)為例，2M階差分精度的高階顯式有限差分法其公式為:
[0013]
【權(quán)利要求】
1.一種實(shí)現(xiàn)地震波場(chǎng)數(shù)值模擬方法，其特征在于包括， 步驟101，在波數(shù)域使用低階有理分式有限差分算子逼近超高階顯式有限差分算子； 步驟102，計(jì)算所述低階有理分式有限差分算子的系數(shù)； 步驟103，將有理分式有限差分算子進(jìn)行分解，利用傅里葉變換將低階有理分式有限差分算子從波數(shù)域變換到空間域； 步驟104，在空間域使用有理分式有限差分法分別計(jì)算不同方向上的空間二階偏導(dǎo)數(shù)項(xiàng)； 步驟105，在時(shí)間方向上進(jìn)行迭代實(shí)現(xiàn)地震波場(chǎng)的數(shù)值模擬。
2.根據(jù)權(quán)利要求1所述的一種實(shí)現(xiàn)地震波場(chǎng)數(shù)值模擬方法，其特征在于，使用如下公式逼近超高階顯式有限差分算子:
3.根據(jù)權(quán)利要求2所述的一種實(shí)現(xiàn)地震波場(chǎng)數(shù)值模擬方法，其特征在于，M的取值為24到28之間，L的取值為I到5之間。
4.根據(jù)權(quán)利要求2所述的一種實(shí)現(xiàn)地震波場(chǎng)數(shù)值模擬方法，其特征在于，所述bm通過如下公式計(jì)算得到:
5.根據(jù)權(quán)利要求2所述的一種實(shí)現(xiàn)地震波場(chǎng)數(shù)值模擬方法，其特征在于，通過最小二乘法、模擬退火，神經(jīng)網(wǎng)絡(luò)或遺傳算法計(jì)算來計(jì)算所述低階有理分式有限差分算子的系數(shù)C1 和 d1 ο
6.根據(jù)權(quán)利要求5所述的一種實(shí)現(xiàn)地震波場(chǎng)數(shù)值模擬方法，其特征在于，使用所述最小二乘法計(jì)算所述低階有理分式有限差分算子的系數(shù)C1和Cl1，具體包括， 步驟201，將有理分式有限差分算子的分母項(xiàng)轉(zhuǎn)化為等價(jià)的L階多項(xiàng)式求和，所得到的計(jì)算公式為:
7.根據(jù)權(quán)利要求1所述的一種實(shí)現(xiàn)地震波場(chǎng)數(shù)值模擬方法，其特征在于，在將有理分式有限差分算子進(jìn)行分解，利用傅里葉變換將低階有理分式有限差分算子從波數(shù)域變換到空間域中進(jìn)一步包括: 步驟301，使用有理分式有限差分計(jì)算空間二階偏導(dǎo)數(shù)，在波數(shù)域通過下面公式來實(shí)現(xiàn):
8.根據(jù)權(quán)利要求7所述的一種實(shí)現(xiàn)地震波場(chǎng)數(shù)值模擬方法，其特征在于，所述在空間域使用有理分式有限差分法分別計(jì)算不同方向上的空間二階偏導(dǎo)數(shù)項(xiàng)進(jìn)一步包括: 步驟401，結(jié)合已知的空間域的波動(dòng)方程波場(chǎng)變量P (h)，使用顯式有限差分法計(jì)算步驟303中的第一個(gè)方程得到qQ(h)； 步驟402，結(jié)合上述步驟401中qQ(h)，依次計(jì)算步驟303中L個(gè)二階隱式方程得到Qi (h), (I = I, 2，…，L)； 步驟403，結(jié)合上述步驟402中qjh)通過步驟303中的第三個(gè)公式，得到二階空間偏導(dǎo)數(shù)結(jié)果￥;crh 步驟404，通過上述步驟401~403，計(jì)算x方向的空間偏導(dǎo)數(shù)







步驟405，通過步驟401~步驟403，計(jì)算y方向的空間偏導(dǎo)數(shù);


步驟406,通過步驟401~步驟403,計(jì)算ζ方向的空間偏導(dǎo)數(shù)。





ζ
9.根據(jù)權(quán)利要求8所述的一種實(shí)現(xiàn)地震波場(chǎng)數(shù)值模擬方法，其特征在于，在所述步驟402中通過追趕法、高斯消元法或矩陣求逆法來求解L個(gè)二階隱式方程。
10.根據(jù)權(quán)利要求8所述的一種實(shí)現(xiàn)地震波場(chǎng)數(shù)值模擬方法，其特征在于，在所述步驟402中通過改進(jìn)追趕法求解L個(gè)二階隱式方程: 步驟601，預(yù)先計(jì)算求解L個(gè)三對(duì)角線性方程組所需的L個(gè)常向量Ci1，通過以下遞推公式來實(shí)現(xiàn): αι ⑴== I) 其中Cl1為前面所求的隱式差分系數(shù)；α I表示采樣長度為Nh的向量，Nh取空間采樣網(wǎng)格數(shù)K、Ny和Nz中的最大值，Nx, Ny和Nz分別為X、y和ζ坐標(biāo)軸方向的空間網(wǎng)格數(shù)； 步驟602，實(shí)現(xiàn)改進(jìn)追趕法的追過程，通過以下遞推公式進(jìn)行:
11.根據(jù)權(quán)利要求8所述的一種實(shí)現(xiàn)地震波場(chǎng)數(shù)值模擬方法，其特征在于，在時(shí)間方向上進(jìn)行迭代實(shí)現(xiàn)地震波場(chǎng)的數(shù)值模擬進(jìn)一步包括:
【文檔編號(hào)】G06F19/00GK103699798SQ201310726137
【公開日】2014年4月2日 申請(qǐng)日期:2013年12月25日 優(yōu)先權(quán)日:2013年12月25日
【發(fā)明者】秦臻, 宋建勇, 馬曉宇, 宋雪娟 申請(qǐng)人:中國石油天然氣股份有限公司