不同高度重力數(shù)據(jù)和井中重力數(shù)據(jù)融合反演研究

高秀鶴,熊盛青,于長(zhǎng)春,孫思源

(1.中國(guó)自然資源航空物探遙感中心, 北京 100083; 2.自然資源部 航空地球物理與遙感地質(zhì)重點(diǎn)實(shí)驗(yàn)室, 北京 100083；3.中國(guó)地質(zhì)大學(xué)(北京) 地球物理與信息技術(shù)學(xué)院, 北京 100029)

0 引言

重力數(shù)據(jù)反演面臨嚴(yán)重的多解性問(wèn)題，高斯定理指出其本質(zhì)原因：有無(wú)數(shù)種不同的密度分布可以引起相同的重力異常[1-3]；此外，異常體引起的重力異常是連續(xù)分布的，而我們采集的信息只是某一平面上的部分離散數(shù)據(jù)，這種信息量的不足進(jìn)一步加重了反演結(jié)果的多解性[4-5]。這種多解性是導(dǎo)致重力數(shù)據(jù)反演結(jié)果趨膚效應(yīng)嚴(yán)重、分辨率低的本質(zhì)原因。利用已知的先驗(yàn)地質(zhì)信息對(duì)密度進(jìn)行合理約束是改善重力數(shù)據(jù)反演效果的有效手段，常用的約束方法有最小質(zhì)量約束、最小體積約束、最大光滑約束、最小慣量約束等[6-10]。另一方面，增加觀測(cè)數(shù)據(jù)的數(shù)據(jù)量也能夠改善反演效果[11-12]，增加的觀測(cè)數(shù)據(jù)可以是同種地球物理數(shù)據(jù)加密或不同平面采集，也可以是不同種類(lèi)地球物理數(shù)據(jù)[13-14]。但在實(shí)際工作中，由于經(jīng)費(fèi)或者測(cè)量條件限制，無(wú)法獲得更多觀測(cè)數(shù)據(jù)時(shí)，通過(guò)延拓技術(shù)獲取多高度面觀測(cè)數(shù)據(jù)加入到反演中是可行的[15]。

垂直分辨率差是重力數(shù)據(jù)反演中最為突出的問(wèn)題[16]，重力數(shù)據(jù)的核函數(shù)隨深度的增加逐漸衰減是密度反演結(jié)果分辨率差的直接原因。從方法技術(shù)角度出發(fā)，引入深度加權(quán)抵消核函數(shù)衰減速率，能夠有效克服反演結(jié)果的趨膚效應(yīng)，提高垂直分辨率[3,17-19]；從實(shí)測(cè)數(shù)據(jù)角度出發(fā)，引入鄰近異常源的井中觀測(cè)重力數(shù)據(jù)、井中巖石密度數(shù)據(jù)，可提高反演結(jié)果的垂直分辨率[20]。

不同高度觀測(cè)面的重力數(shù)據(jù)包含的波場(chǎng)信息不同[21]，因此，單獨(dú)反演效果存在差異，本文首先設(shè)計(jì)模型試驗(yàn)分析了不同高度面觀測(cè)重力數(shù)據(jù)的反演效果；然后，融合反演不同高度面觀測(cè)的重力數(shù)據(jù)，與單一高度面觀測(cè)重力數(shù)據(jù)反演結(jié)果對(duì)比，分析融合反演的效果，同時(shí)，多觀測(cè)面數(shù)據(jù)不可得時(shí)，利用延拓技術(shù)獲得不同高度重力數(shù)據(jù)用于反演；最后，分析井中重力數(shù)據(jù)、井中密度數(shù)據(jù)約束下，不同高度重力數(shù)據(jù)反演效果。

1 融合反演原理

離散重力數(shù)據(jù)與離散密度塊體的關(guān)系如下：

g=Aρ，

(1)

式中：ρ是密度向量，g是某一觀測(cè)面的觀測(cè)重力數(shù)據(jù)向量，A是正演算子矩陣。反演過(guò)程中，已知觀測(cè)重力數(shù)據(jù)向量g和正演算子矩陣A，求解密度向量ρ。由于觀測(cè)數(shù)據(jù)中不可避免含有噪聲干擾，使得反演問(wèn)題病態(tài)，導(dǎo)致解不穩(wěn)定，吉洪諾夫正則化方法作為一種常規(guī)的處理手段，可有效解決這一問(wèn)題，求式(1)中ρ的反演問(wèn)題可轉(zhuǎn)化為求解目標(biāo)函數(shù)最優(yōu)解的問(wèn)題，目標(biāo)函數(shù)如下：

φ=φg+αφρ，

(2)

式中：φg是數(shù)據(jù)擬合函數(shù)，φρ是模型擬合函數(shù)，α代表正則化參數(shù)。數(shù)據(jù)擬合函數(shù)為

(3)

模型擬合函數(shù)為最小模型約束：

(4)

在模型擬合項(xiàng)中加入深度加權(quán)函數(shù)W，可有效克服重力反演結(jié)果中的趨膚問(wèn)題。反演結(jié)果中異常體的埋深很大程度上受深度加權(quán)函數(shù)影響，弱的深度加權(quán)函數(shù)容易獲得較淺的異常體埋深，強(qiáng)的深度加權(quán)則容易獲得較深的異常體埋深。

引入深度加權(quán)函數(shù)后，模型擬合函數(shù)可以寫(xiě)成

(5)

不同平面(高度或深度為z)重力數(shù)據(jù)融合反演時(shí)，數(shù)據(jù)擬合函數(shù)不同。觀測(cè)重力數(shù)據(jù)gz表示觀測(cè)面高度或深度為z的重力數(shù)據(jù)。將gz作為觀測(cè)數(shù)據(jù)加入到目標(biāo)函數(shù)的數(shù)據(jù)擬合函數(shù)中，則數(shù)據(jù)擬合函數(shù)為

(6)

式中：Az是與gz對(duì)應(yīng)的正演算子矩陣。此外，還可以在數(shù)據(jù)擬合函數(shù)中加入井中重力數(shù)據(jù)gb、井中巖石密度數(shù)據(jù)ρd，則有

(7)

其中：gb代表井中重力數(shù)據(jù)，Ab是與gb對(duì)應(yīng)的正演算子，ρd代表井中巖石密度數(shù)據(jù)，Ad可以看成是與ρd對(duì)應(yīng)的正演算子，它是一個(gè)稀疏矩陣，矩陣每行僅有一個(gè)非0元素，與取樣巖石位置有關(guān)，其值為1。式(7)所示的數(shù)據(jù)擬合函數(shù)加入到目標(biāo)函數(shù)中，可實(shí)現(xiàn)不同高度觀測(cè)面、井中數(shù)據(jù)的融合反演，通過(guò)調(diào)整系數(shù)λ1、λ2、λb和λd的大小可調(diào)整不同高度面觀測(cè)重力數(shù)據(jù)、井中重力數(shù)據(jù)、井中巖石密度數(shù)據(jù)對(duì)目標(biāo)函數(shù)的作用大?。合禂?shù)越小，相應(yīng)數(shù)據(jù)對(duì)目標(biāo)函數(shù)影響越??；系數(shù)越大，相應(yīng)數(shù)據(jù)對(duì)目標(biāo)函數(shù)影響越大；當(dāng)系數(shù)為0時(shí)，相應(yīng)數(shù)據(jù)對(duì)目標(biāo)函數(shù)不起作用，即相應(yīng)數(shù)據(jù)不參與反演過(guò)程。

2 模型試驗(yàn)

2.1 模型一：?jiǎn)瘟⒎襟w模型

為了驗(yàn)證不同平面觀測(cè)重力數(shù)據(jù)三維反演效果，建立一個(gè)立方體模型，如圖1a所示，立方體模型的大小為400 m×400 m×400 m，頂面埋深為300 m，中心位置為(2.5 km, 2.5 m, 0.5 km)，密度差為1 g·cm-3，背景密度為0。將地下三維空間剖分為50×50×15個(gè)大小為100 m×100 m×100 m的立方體單元，50×50個(gè)觀測(cè)點(diǎn)位于每個(gè)單元的中心。正演計(jì)算得到的不同高度觀測(cè)面的重力數(shù)據(jù)如圖1b、c、d、e、f所示。

利用正則化方法反演不同高度面觀測(cè)重力數(shù)據(jù)。為了對(duì)比反演效果，對(duì)于不同高度面的觀測(cè)重力數(shù)據(jù)反演，引入相同的密度上下限 (0, 1 g·cm-3)和相同的模型加權(quán)函數(shù)W=(ATA)1/2，其中，A表示地面重力數(shù)據(jù)的正演算子。利用共軛梯度算法迭代求解目標(biāo)函數(shù)最優(yōu)解，得到的反演結(jié)果切片如圖2所示。觀測(cè)面高度越高，相應(yīng)的反演結(jié)果越發(fā)散；觀測(cè)面高度越低，相應(yīng)的反演結(jié)果越收斂，且反演結(jié)果與理論模型的分布范圍、密度值大小越接近。

不同高度面的觀測(cè)重力數(shù)據(jù)融合反演結(jié)果如圖3所示。圖3a與圖2c相比，由于地面以下300 m觀測(cè)數(shù)據(jù)的加入，密度分布更收斂，密度值更接近理論密度值；圖3a與圖2f相比，由于地面觀測(cè)數(shù)據(jù)的加入，壓制了異常體兩側(cè)的假異常。這初步驗(yàn)證了不同高度面觀測(cè)重力數(shù)據(jù)融合反演的良好效果。相比于圖2c，圖3a顯然與圖2f更接近，說(shuō)明近源平面數(shù)據(jù)對(duì)反演結(jié)果影響更大。圖3b與圖2b、d相比，密度分布更收斂，密度值更接近理論密度值，與圖2f相比，壓制了異常體兩側(cè)的假異常。這再次驗(yàn)證了不同高度面觀測(cè)重力數(shù)據(jù)融合反演的良好效果。圖3c與圖3b相比，增加了兩個(gè)面的觀測(cè)數(shù)據(jù)，反演效果并未見(jiàn)明顯提高。這說(shuō)明為保證計(jì)算效率，需選擇適當(dāng)、適量平面的觀測(cè)重力數(shù)據(jù)進(jìn)行融合反演，過(guò)多的數(shù)據(jù)將導(dǎo)致計(jì)算量增大，但反演效果并沒(méi)有明顯的提高。

僅有一個(gè)觀測(cè)面的重力數(shù)據(jù)時(shí)，可通過(guò)延拓技術(shù)得到多個(gè)不同高度面的延拓重力數(shù)據(jù)。使用延拓?cái)?shù)據(jù)反演需注意幾個(gè)問(wèn)題:一是為提高反演效果，使用向下延拓重力數(shù)據(jù)進(jìn)行反演，原因是：?jiǎn)蝹€(gè)平面重力數(shù)據(jù)反演驗(yàn)證，觀測(cè)面高度越低，反演效果越好；二是為保證反演可靠性，下延平面應(yīng)高于異常體頂面，原因是：穿過(guò)或低于異常體的延拓?cái)?shù)據(jù)不準(zhǔn)確，不適合用于反演；三是為了擬合原始觀測(cè)重力數(shù)據(jù)，使用延拓?cái)?shù)據(jù)與原始觀測(cè)數(shù)據(jù)融合反演，原因是：向下延拓得到的重力數(shù)據(jù)存在誤差，單獨(dú)反演延拓?cái)?shù)據(jù)將導(dǎo)致反演結(jié)果可能無(wú)法擬合原始觀測(cè)數(shù)據(jù)。

圖1 立方體模型及正演重力數(shù)據(jù)Fig.1 Cube model and its forward gravity data

圖2 不同高度面觀測(cè)重力數(shù)據(jù)單獨(dú)反演結(jié)果的豎直切片(y=2.5 km) Fig.2 Vertical slices of inversion results of gravity data at different heights

圖3 不同高度面觀測(cè)重力數(shù)據(jù)融合反演結(jié)果的豎直切片(y=2.5 km)Fig.3 Vertical slices of inversion results of joint gravity data at different heights(y=2.5 km)

觀測(cè)重力數(shù)據(jù)及其延拓重力數(shù)據(jù)融合反演結(jié)果如圖4所示，圖4a與圖3a相比，密度分布較發(fā)散，密度值較小，這是由于向下延拓?cái)?shù)據(jù)中存在誤差，說(shuō)明延拓?cái)?shù)據(jù)不能取代觀測(cè)數(shù)據(jù)；與圖2c相比，密度分布更收斂，密度值更接近真實(shí)密度值，反演效果明顯提高，這說(shuō)明，在多平面數(shù)據(jù)不可得時(shí)，延拓?cái)?shù)據(jù)的加入可提高反演效果。圖4b與圖3b相比，密度分布較發(fā)散；但與圖2b相比，反演效果明顯提高，再次驗(yàn)證了加入延拓?cái)?shù)據(jù)的意義。圖4c與圖4b對(duì)比，反演效果提高并不明顯。因此，為保證計(jì)算效率，選擇適當(dāng)、適量的延拓?cái)?shù)據(jù)加入反演是至關(guān)重要的，過(guò)多的延拓?cái)?shù)據(jù)導(dǎo)致計(jì)算量增大，但反演效果并沒(méi)有明顯提高。

圖4 觀測(cè)重力數(shù)據(jù)及其延拓重力數(shù)據(jù)融合反演結(jié)果的豎直切片(y=2.5 km) Fig.4 Inverted density models using the joint observed gravity and continuation gravity data(y=2.5 km)

利用單立方體模型驗(yàn)證加入井?dāng)?shù)據(jù)約束對(duì)反演效果的影響。模擬設(shè)計(jì)3口豎直井，3口豎直井分別位于立方體中心、立方體外100 m、立方體外300 m，其坐標(biāo)分別為(2.5 km,2.5 km)、(2.8 km,2.5 km)、(3.0 km,2.5 km)。井中觀測(cè)重力數(shù)據(jù)、巖石密度數(shù)據(jù)如圖5所示。在地面、地面以下300 m重力數(shù)據(jù)融合反演中，分別加入3口井?dāng)?shù)據(jù)做約束，得到的反演結(jié)果如圖6所示。與圖3a相比，反演效果明顯提高，這說(shuō)明加入位于立方體中心的井?dāng)?shù)據(jù)做約束可有效改善反演效果；加入位于立方體外100 m的井?dāng)?shù)據(jù)做約束對(duì)鄰近井位的密度分布有一定改善作用，對(duì)整體反演結(jié)果影響較?。患尤胛挥诹⒎襟w外300 m的井?dāng)?shù)據(jù)做約束對(duì)反演結(jié)果的影響很小。這說(shuō)明，井口穿過(guò)異常體，能有效改善反演效果；井口距離異常體越遠(yuǎn)，對(duì)反演效果的有益影響越小。

圖5 不同位置井的井中重力和巖石密度Fig.5 Well gravity and well density at different locations

圖6 地面、地面以下300 m重力數(shù)據(jù)和井中重力、巖石密度數(shù)據(jù)融合反演結(jié)果(白色實(shí)線代表鉆井位置)Fig.6 Recovered density models in the constraint of well data(the solid white line represents the drilling location)

2.2 模型二：雙長(zhǎng)方體模型

為進(jìn)一步驗(yàn)證不同高度重力數(shù)據(jù)和井?dāng)?shù)據(jù)融合反演的效果，設(shè)置不同埋深的雙長(zhǎng)方體模型，左側(cè)長(zhǎng)方體頂面埋深300 m，中心位置坐標(biāo)為(550 m,1 050 m,500 m)，大小為300 m×300 m×400 m，密度為1 g·cm-3，右側(cè)長(zhǎng)方體頂面埋深為500 m，中心位置坐標(biāo)為(1 350 m,1 050 m,700 m)，其大小、密度與左側(cè)長(zhǎng)方體相同。地下空間剖分為38×42×15個(gè)大小為50 m×50 m×100 m的立方體單元，38×42個(gè)地面觀測(cè)數(shù)據(jù)采樣間距為50 m，井中觀測(cè)數(shù)據(jù)采樣間距為50 m。加入3%的高斯噪聲的正演理論地面重力數(shù)據(jù)如圖7a所示，將地面重力數(shù)據(jù)向下延拓300 m 得到的延拓重力數(shù)據(jù)如圖7b所示。模擬兩口豎直井：井A、井B，其水平位置分別為(550 m,1 050 m)、(1 350 m,1 050 m)，即位于兩個(gè)長(zhǎng)方體中心，填加3%高斯噪聲的井中重力、巖石密度數(shù)據(jù)如圖8所示。

地面重力數(shù)據(jù)單獨(dú)反演結(jié)果如圖9b所示，左側(cè)埋深較淺的異常體中心位置識(shí)別準(zhǔn)確，右側(cè)埋深較深的異常體未被有效識(shí)別，且無(wú)法區(qū)分兩個(gè)異常體，反演結(jié)果整體較為發(fā)散，分辨率較低。加入向下延拓300 m重力數(shù)據(jù)與地面重力數(shù)據(jù)融合反演結(jié)果如圖9c，與圖9b相比，左側(cè)異常體中心位置定位準(zhǔn)確，且異常體邊界清晰，密度值的最大值達(dá)到真實(shí)值，右側(cè)異常體反演效果仍不理想，且由于下延數(shù)據(jù)中噪聲很大，導(dǎo)致反演結(jié)果中也出現(xiàn)了假異常。再加入兩口井的井中重力、巖石密度數(shù)據(jù)約束，反演結(jié)果如圖9d，準(zhǔn)確定位了兩個(gè)不同埋深的異常體的中心位置，邊界清晰，分辨率高，兩個(gè)異常體的最大異常值都達(dá)到了理論值，與圖9a所示的理論模型剖面匹配度較高。

圖7 雙長(zhǎng)方體模擬重力數(shù)據(jù)Fig.7 Gravity data based on double prisms

圖8 含3%高斯噪聲的模擬井中觀測(cè)數(shù)據(jù)Fig.8 The simulated data observed in the wells

圖9 雙立方體模型反演結(jié)果 (白色實(shí)線代表鉆井位置)Fig.9 Recovered double-prism models(the solid white line represents the drilling location)

3 結(jié)論

1) 單一觀測(cè)面重力數(shù)據(jù)三維反演效果與測(cè)量面高度有關(guān)，觀測(cè)面高度越低(觀測(cè)面距離異常源越近)，反演效果越好。

2) 與單一觀測(cè)面重力數(shù)據(jù)反演相比，多個(gè)高度面觀測(cè)重力數(shù)據(jù)融合反演的效果可得到改善，其中，近源高度面數(shù)據(jù)比遠(yuǎn)源高度面數(shù)據(jù)對(duì)反演結(jié)果影響大。

3) 過(guò)多個(gè)觀測(cè)面的重力數(shù)據(jù)對(duì)反演效果改善不明顯，因此，為保證計(jì)算效率，應(yīng)選擇適當(dāng)、適量的觀測(cè)面的重力數(shù)據(jù)進(jìn)行融合反演。

4) 當(dāng)多個(gè)高度面的觀測(cè)重力數(shù)據(jù)不可得時(shí)，通過(guò)向下延拓技術(shù)可得到多個(gè)平面的延拓重力數(shù)據(jù)(延拓平面應(yīng)高于異常體頂面)，觀測(cè)重力數(shù)據(jù)與向下延拓重力數(shù)據(jù)融合反演，也可改善反演效果。

5) 井中重力、巖石密度數(shù)據(jù)約束，可一定程度改善反演效果，改善程度受鉆井位置影響明顯，穿過(guò)異常體的井?dāng)?shù)據(jù)可明顯改善反演效果，鉆井位置距離異常體越遠(yuǎn)，對(duì)反演效果影響越小。