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