馬國慶,李宗睿,李麗麗
吉林大學(xué) 地球探測科學(xué)與技術(shù)學(xué)院,長春 130026
傾斜臺階可用于模擬實際地層的斷裂、接觸帶等,利用其重力正演公式可獲得斷裂等相似地質(zhì)體的重力異常特征,從而為數(shù)據(jù)解釋提供依據(jù)。對于臺階正演結(jié)果的反演可模擬對于斷裂等相關(guān)地質(zhì)體的反演,從而對反演的可靠性進(jìn)行評估;在一些需要迭代計算的反演中還需要正演公式參與運算。因此準(zhǔn)確的臺階公式至關(guān)重要。
傾斜臺階的重力異常表達(dá)式在國內(nèi)外眾多教材及文章[1-8]中出現(xiàn)并應(yīng)用。很多斷裂特征反演的方法[9-11]大多是通過臺階正演計算公式來進(jìn)行推導(dǎo),因此重力異常正演公式的正確性至關(guān)重要。
對現(xiàn)有教材及論文中正演公式進(jìn)行模型試驗表明,在一定條件下傾斜臺階的重力異常曲線會出現(xiàn)不規(guī)則的畸變。馮蘭天[12]和蘇和朋等[13]也通過模型正演證明了臺階重力異常計算公式所存在的畸值問題,但其未針對產(chǎn)生畸值的條件進(jìn)行分析,且在公式推導(dǎo)中存在一定的錯誤。
筆者對現(xiàn)有的正演計算公式及結(jié)果展開分析,明確了現(xiàn)有公式產(chǎn)生畸值的原因及條件,并基于臺階的正演計算積分展開式,推導(dǎo)獲得了新的臺階重力異常正演計算公式。模型試驗表明重新推導(dǎo)公式所獲得數(shù)據(jù)不存在畸點,且發(fā)現(xiàn)原有公式是由于在合并中未考慮反三角函數(shù)值域而出現(xiàn)的誤差。
以水平方向為x軸,垂直方向為y軸,建立坐標(biāo)系,繪制傾斜臺階模型(圖1)。
圖1 傾斜臺階模型Fig.1 Tilt step model
在圖1所示情況下,傾斜臺階常規(guī)計算公式為:
(1)
式(1)中,G為萬有引力常數(shù);ρ為場源與圍巖之間的密度差;x為地面上觀測點的坐標(biāo);h為傾斜臺階上頂面的埋深;H為臺階下底面的埋深;α為臺階的傾斜角。
根據(jù)式(1),設(shè)傾斜臺階模型參數(shù)為:h=10 m,H=100 m,傾角分別取45°、60°、90°、120°、135°,場源與圍巖之間的密度差異設(shè)置為1 000 kg/m3,萬有引力常數(shù)設(shè)為6.670 8×10-11m3kg-1s-2,根據(jù)式(1)得到不同傾角的傾斜臺階重力異常曲線(圖2)。
圖2 傾斜臺階模型重力異常曲線(α=45°、60°、90°、120°、135°)Fig.2 Gravity anomaly curve of tilt step model
從圖2中可以看到在傾角為60°、90°、120°時使用該公式,重力異常曲線并未產(chǎn)生畸變。因此并非所有情況下使用該公式都會產(chǎn)生畸變。為進(jìn)一步進(jìn)行分析,保持其余參數(shù)不變傾角分別取30°、45°、135°、150°,得到相應(yīng)的重力異常曲線(圖3)。
圖3 傾斜臺階模型重力異常曲線(α=30°、45°、135°、150°)Fig.3 Gravity anomaly curve of tilt step model
為研究畸點產(chǎn)生的原因,針對式(1)進(jìn)行分析。
由于該公式項數(shù)過多,不便于整體分析,所以考慮逐項討論其連續(xù)性。又因為常數(shù)項不影響函數(shù)的連續(xù)性,因此僅對含自變量的項進(jìn)行分析。
令:
(2)
(3)
(4)
(5)
則:
Δg=GΔρ[π(H-h)+f1-f2+f3-f4]
(6)
顯然分析Δg的連續(xù)性可轉(zhuǎn)化為分析f1,f2,f3,f4的連續(xù)性。
(7)
為解決畸變問題,對傾斜臺階重力異常公式重新進(jìn)行推導(dǎo)。
引力場F是保守場(沿閉合路線做功為0)或無旋場,考慮到標(biāo)量函數(shù)梯度的旋度恒等于零,可引入引力位V(標(biāo)量函數(shù)):
(8)
式(8)說明引力的方向始終指向引力位增加最快的方向。
對于質(zhì)量為m的質(zhì)點引力場外,某點的引力位的定義為:將單位質(zhì)量的質(zhì)點從無窮遠(yuǎn)移至該點時引力場所做的功,即
(9)
質(zhì)體外的引力位為:
(10)
對于二度體來說(所謂的二度體就是橫截面的形狀和深度沿某一水平方向不變且該方向無限延伸的物體),以地面上某一點O為坐標(biāo)原點,z軸鉛直向下即沿重力方向建立笛卡爾坐標(biāo)系。若物體的剩余密度均勻,則二度體的重力異??杀硎緸椋?/p>
(11)
在圖1所示模型下,臺階在x軸上任意一點P(x,0)的重力異常應(yīng)為:
(12)
得出重力異常的表達(dá)式,需要解此定積分,詳細(xì)過程為:
(13)
帶入到重力異常表達(dá)式(12)中有:
(14)
對于式(14)中積分的前一項:
(15)
(16)
(17)
運用分部積分法,則有:
(18)
(19)
則有:At2+Bt-Atcotα-Bcotα+ct2+c=x,得到:
(20)
解得:
(21)
因此有:
(22)
再帶入到原來的重力異常表達(dá)式(14)中得:
(23)
對比傳統(tǒng)的傾斜臺階重力正演公式,由于:
(24)
其中:a,b∈R,c為整數(shù)。
傳統(tǒng)公式實際上是對新公式后兩項進(jìn)行合并后的結(jié)果,而合并時未考慮c,因此產(chǎn)生了間斷點。傳統(tǒng)的傾斜臺階重力異常公式本身并沒有錯誤,只是需要在合并時對c進(jìn)行分析。
由反正切函數(shù)性質(zhì)可知:
(25)
則:
arctan(a)+arctan(b)∈(-π,π)
(26)
(27)
又由于:
(28)
因此:
(29)
即:
(30)
同理可得:
(31)
(32)
且證:
(33)
(34)
綜上所述:
arctan(a)+arctan(b)
(35)
因此,若要傳統(tǒng)公式(1)與新公式(23)等價可將其寫為:
(36)
式(36)中u=x(H-h)sin2α,v=x2sin2α+(H+h)xsinαcosα+Hh。
運用新推導(dǎo)的傾斜臺階重力異常公式(23)以及加入判定條件后的傳統(tǒng)公式(36)分別用Mathematics進(jìn)行正演計算。使用相同的參數(shù),即設(shè)傾斜臺階模型參數(shù)為:h=10 m,H=100 m,傾角分別取15°、45°、90°、135°、165°,場源與圍巖之間的密度差異設(shè)置為1 000 kg/m3,萬有引力常數(shù)設(shè)為6.670 8×10-11m3kg-1s-2,分別得到圖4、圖5。
圖4 新傾斜臺階重力異常公式得到的重力異常曲線Fig.4 Gravity anomaly curve obtained by new formula of gravity anomaly of inclined steps
圖5 加入限定條件后傳統(tǒng)臺階重力異常公式得到的重力異常曲線Fig.5 Gravity anomaly curve obtained by traditional step gravity anomaly formula after adding qualified conditions
從圖4可以看出新推導(dǎo)出的傾斜臺階重力異常公式是正確的。而圖5得到了與圖4完全相同的結(jié)果,因此傳統(tǒng)傾斜臺階重力異常公式在加入適當(dāng)?shù)南薅l件后同樣正確。且當(dāng)傾角為任意值時,曲線均是連續(xù)的,有效地解決了以往公式中出現(xiàn)的畸變問題。
(1)現(xiàn)有傾斜臺階的重力正演公式由于在推導(dǎo)過程中未考慮反三角函數(shù)的值域變化造成了畸點的出現(xiàn)。
(2)從臺階重力異常積分公式出發(fā),對傾斜臺階重力異常正演公式重新進(jìn)行了正確的推導(dǎo),模型試驗表明新推導(dǎo)公式有效地避免了畸點的產(chǎn)生。