平原區(qū)深層隱伏巖溶的重力正演模擬

王立發(fā) 雷曉東 何祎

摘? 要：巖溶塌陷一般發(fā)生在覆蓋層較薄的地區(qū)，深層巖溶塌陷雖較少見，但若發(fā)生在城市建筑密集區(qū)其危害極大。為研究北京東部地區(qū)深層巖溶發(fā)育的地質(zhì)條件，對不同尺度的巖溶風(fēng)化殼和溶洞可能產(chǎn)生的重力異常響應(yīng)進行了正演模擬，結(jié)果表明：厚覆蓋層地區(qū)開展微重力測量識別一定規(guī)模的巖溶是可行的，覆蓋層厚度大于200m、半徑小于2m的空洞將難以通過微重力方法直接識別出來;巖溶空洞引起重力異常的幅值與其半徑、埋深及內(nèi)部充填物的密度有關(guān);實測剖面長度應(yīng)大于巖溶風(fēng)化殼寬度或空洞直徑，且需要在實測重力異常中準(zhǔn)確去除覆蓋層的影響。

關(guān)鍵詞：重力;巖溶;正演

中圖分類號：P631.1;P642.2? 文獻標(biāo)識碼：A? ? ?文章編號：1007-1903（2018）04-0019-05

0 引言

巖溶塌陷是巖溶發(fā)育區(qū)主要的地質(zhì)災(zāi)害類型之一，在碳酸鹽巖地層之上第四系覆蓋層較薄（一般小于50m）的地區(qū)尤為常見。深層巖溶塌陷發(fā)生在覆蓋層較厚的地區(qū)，一般較為少見，但若發(fā)生在城市建筑密集區(qū)，一旦發(fā)生可能會造成重大危害。國內(nèi)外巖溶發(fā)育區(qū)已發(fā)生過多起深層巖溶塌陷，北京平原區(qū)也不例外。例如通州北部的龍旺莊地區(qū)于1995年曾發(fā)生深層巖溶塌陷，塌陷深度289.29m（鐘立勛，2001），順義漢石橋濕地附近的巖溶水勘探孔也曾于2013年發(fā)生塌陷，塌陷深度283m。通州地區(qū)作為北京城市副中心建設(shè)區(qū)，同時也是深層巖溶發(fā)育的地區(qū)，該區(qū)碳酸鹽巖地層以寒武系灰?guī)r、薊縣系白云巖為主，上覆第四系松散層，松散層厚度200～700m。為避免副中心重大工程建設(shè)活動誘發(fā)深層巖溶塌陷，我們對這一帶巖溶發(fā)育條件進行了調(diào)查（李巧靈等，2018）。因該區(qū)為平原覆蓋區(qū)，調(diào)查手段以地球物理方法為主。微重力測量通過探測密度變化推斷地下異常體位置和幾何尺寸，在淺地表勘探領(lǐng)域應(yīng)用效果較好（王謙身等，1995;楊進，2011）。隨著儀器精度的提高，對異常體的識別能力越來越強，巖溶發(fā)育區(qū)深層發(fā)育大型溶洞或破碎帶時，能否通過微重力測量被發(fā)現(xiàn)？為此，我們使用正演模擬方法研究了不同類型巖溶的重力異常響應(yīng)，為后續(xù)研究巖溶發(fā)育地質(zhì)條件提供依據(jù)。

1 方法原理

大量研究表明，巖溶塌陷形成的關(guān)鍵因素包括地質(zhì)條件和地下水動力條件，其中地質(zhì)條件包括巖溶發(fā)育程度、覆蓋層厚度與性質(zhì)、與活動斷裂的距離等;地下水動力條件主要指水位變幅、水位降深、水位距基巖面位置等（石建省，1996;胡瑞林等，2001;雷曉東等，2010;趙博超等，2015;廖海軍等，2015）。微重力測量是采用高精度重力儀觀測重力值，探測微小目標(biāo)體，精度達到微加級的一種地球物理勘探方法，在地災(zāi)、考古、地?zé)岬阮I(lǐng)域應(yīng)用較多，它與傳統(tǒng)區(qū)域重力測量的區(qū)別是測區(qū)面積相對較小、測點分布較密（王延濤等，2012;袁炳強等，2015;Pearson-Grant， et al，2018）。應(yīng)用在巖溶探測方面其原理是利用穩(wěn)定圍巖和孔洞或破碎帶之間的密度差異，通過正反演計算確定目標(biāo)體的埋深和規(guī)模。微重力測量的主要影響因素來自于測點周圍的地形起伏，實測過程中遇到地形影響需要進行校正，校正的誤差將會影響探測效果。在淺層地質(zhì)體中發(fā)育的溶洞在實測重力異常剖面曲線上會產(chǎn)生明顯的負(fù)異常，異常極大值對應(yīng)于溶洞的中心位置，在平面等值線圖上一般形成橢圓形重力負(fù)異常區(qū)（石亞雄等，1991;陳貽祥，1995），局部發(fā)育的風(fēng)化殼產(chǎn)生較寬緩的負(fù)異常，而在深層地質(zhì)體中發(fā)育的溶洞或巖溶風(fēng)化殼產(chǎn)生的重力異常將大大減弱。

2 正演模型建立

巖溶發(fā)育程度是巖溶塌陷形成的關(guān)鍵因素，巖溶較為發(fā)育的地區(qū)廣泛分布有地下暗河、大型溶洞、管道系統(tǒng)或者巨厚風(fēng)化殼。為研究方便，正演模型分為兩類，一類是風(fēng)化殼，一類是溶洞。

（1）巖溶風(fēng)化殼模型

模型建立：假定第四系、巖溶風(fēng)化殼及基巖為密度均勻的水平層狀介質(zhì)。模型最大X坐標(biāo)范圍為±5000m，最大Y坐標(biāo)范圍為±200m，巖溶發(fā)育中心X坐標(biāo)為0m。第四系密度ρ=2.0g/cm3;巖溶風(fēng)化殼密度ρ=1.0g/cm3;基巖密度ρ=2.67g/cm3。共設(shè)定以下6種模型（圖1）：

A1：無巖溶發(fā)育，第四系厚度200m;

A2：層狀巖溶發(fā)育，第四系厚度200 m，巖溶寬度1000m，厚度50 m;

A3：局部巖溶發(fā)育，第四系厚度200 m，巖溶寬度200m，厚度50 m;

A4：層狀巖溶發(fā)育，第四系厚度200 m，巖溶寬度5000m，厚度50 m;

A5：層狀巖溶發(fā)育，第四系厚度500 m，巖溶寬度5000m，厚度100 m;

A6：層狀巖溶發(fā)育，第四系厚度500 m，巖溶寬度2000m，厚度100 m。

（2）溶洞模型

模型建立：當(dāng)?shù)貙禹敯甯浇l(fā)育有巖溶空洞，并以第四系等松散沉積物作為覆蓋層時，用水平圓柱體與水平物質(zhì)帶組合近似模擬其重力異常大小，模型如圖2所示。

3 結(jié)果分析

3.1 巖溶風(fēng)化殼

選用相同的背景場值，通過正演得到6種模型的重力場數(shù)據(jù)。假設(shè)實際測線長度1km，測線中心點位于X軸0點，截取各正演模型X軸上-500m至500m段重力場數(shù)據(jù)，用于分析重力探測隱伏巖溶風(fēng)化殼的異常大小，正演重力異常曲線見圖3。

（1）巖溶風(fēng)化殼寬度大于測線長度時，能否探測到有意義的重力異常？

A4、A5模型中巖溶寬度5000m，其重力異常曲線與無巖溶發(fā)育的水平地層（A1模型）相似，均為近水平直線。對比A1與A4模型，當(dāng)在200～250m深度發(fā)育巖溶時（密度ρ=1.0g/cm3，寬度5000m），引起的重力場變化為1.5001×10-5m/s2;對比A1與A5模型，當(dāng)在500～600m深度發(fā)育巖溶時（密度ρ=1.0g/cm3，寬度5000m），引起的重力場變化為2.86347×10-5m/s2;重力異常值的不同主要是巖溶埋藏深度不同的反映。可以認(rèn)為，當(dāng)巖溶的寬度遠遠超出測線長度時，難以通過重力異常曲線形態(tài)識別巖溶地層。

（2）巖溶風(fēng)化殼寬度小于等于測線長度，能否探測到有意義異常？局部發(fā)育巖溶的情況

對比A1與A2模型，當(dāng)在200～250m深度發(fā)育巖溶時（密度ρ=1.0g/cm3，寬度1000m），引起的重力場變化為1.344×10-5m/s2;對比A1與A3模型，當(dāng)在200～250m深度局部發(fā)育巖溶時（密度ρ=1.0g/cm3，寬度200m），引起的重力場變化為0.58369×10-5m/s2;由此可見當(dāng)測線長度大于風(fēng)化殼寬度時，能探測到有意義的重力異常。

（3）重力異常的影響因素分析

單一修改A1模型中第四系厚度（200m與500m）所產(chǎn)生的重力異常值變化量（△G1）為2.557×10-5m/s2;對比A2與A6模型，當(dāng)?shù)谒南岛穸取r溶風(fēng)化殼的寬度和厚度均不同時（A2：第四系厚200m、寬度1000m、厚度50m;A6：第四系厚500m、寬度2000m、厚度100m），A2模型的重力異常為-5.896×10-5m/s2、A6模型的重力異常為-8.038×10-5m/s2，兩者相差2.143×10-5m/s2（△G2）。△G1和△G2差值較小，由此可知，在探測深層隱伏巖溶目標(biāo)時，重力異常值大小主要取決于第四系覆蓋層的厚度;而巖溶風(fēng)化殼的寬度和厚度變化對重力異常的貢獻相對較小。因此要從微重力觀測異常中提取出巖溶目標(biāo)體有效異常，最重要的處理步驟是分離第四系覆蓋層的影響。

3.2 巖溶空洞

第四系厚度（空洞中心埋深）相同，空洞半徑不同時的重力異常大小對比。假定第四系厚度H=200m，空洞密度1.0 g/cm3，空洞半徑分別為5m、10m、50m、100m，重力異常值隨著空洞半徑R的增加而逐漸增加，最大重力異常分別為-5.588572×10-5m/s2、

-5.613364×10-5m/s2、-6.279914×10-5m/s2、

-7.912962×10-5m/s2（沒有去除覆蓋層的影響）;設(shè)定的最大與最小半徑溶洞（R=100m與R=5m）重力異常值之差為2.324390×10-5m/s2。

不同第四系厚度（200m、300m、400m、500m），不同空洞半徑（R=5 m、10 m、50 m、100m）、固定空洞密度（1.0g/cm3）時重力正演異常（已去除區(qū)域場），由表1可見，空洞半徑5m、第四系厚度500m時，空洞引起的重力異常為-0.003465mGal，理論上是可以被識別的。不同半徑空洞引起的重力異常隨空洞半徑增大而增大，隨空洞埋深（上覆第四系厚度）的增加而減小。

改變空洞密度。空洞密度為1.5 g/cm3時的重力異常值在去除第四系覆蓋層影響后，得到純由巖溶空洞引起的重力異常值（表2）。對比表2和表1，空洞密度變大后，不同埋深空洞產(chǎn)生的重力異常相對變小，但理論上仍然可被識別。空洞中極少或無充填物的情況下（密度0.1 g/cm3），一般空洞的半徑較小，本次模擬了0.5～4m不同半徑的、不同埋深的空洞引起的重力異常極值（沒有去除覆蓋層影響）。由表3可知，這些重力異常在-0.000054 ～ -0.008448mGal之間，當(dāng)埋深大于200m或空洞半徑小于2m時，空洞產(chǎn)生的重力異常小于0.001mGal，現(xiàn)有儀器將難以分辨出其埋深和形態(tài)。同一埋深不同半徑空洞（0.5m、4m）之間重力異常存在0.003365～0.008314mGal的差異;同一半徑不同埋深（200m、500m）的空洞之間的重力異常存在0.00008～0.005029mGal的差異。由此可知，低密度空洞重力異常的差異受半徑的影響比埋深的影響大。

4 結(jié)論

（1）重力正演計算表明，通過微重力測量可識別出厚覆蓋層之下發(fā)育的規(guī)模較大的巖溶破碎帶或空洞，實測剖面長度應(yīng)大于巖溶破碎帶寬度或空洞直徑。

（2）巖溶風(fēng)化殼或空洞引起重力負(fù)異常的幅值與風(fēng)化殼埋深和厚度，空洞的半徑、埋深及內(nèi)部充填物的密度有關(guān)。

（3）覆蓋層厚度大于200m，半徑小于2m的空洞將難以通過微重力方法直接識別出來。深層隱伏巖溶空洞探測需要在實測重力異常中準(zhǔn)確去除覆蓋層的影響，因覆蓋層厚度大，巖性和結(jié)構(gòu)變化復(fù)雜，導(dǎo)致實現(xiàn)這一點比較困難。

參考文獻

陳貽祥，1995. 地面微重力方法在地質(zhì)災(zāi)害調(diào)查中的應(yīng)用效果[J]. 中國巖溶，（2）： 176-185.

胡瑞林，王思敬，李焯芬，等，2001. 唐山市巖溶塌陷區(qū)域風(fēng)險評價[J]. 巖石力學(xué)與工程學(xué)報，20（2）： 180-180.

雷曉東，姚春梅，邵景力，等，2010. 臨沂市城區(qū)巖溶塌陷危險性控制[J]. 地學(xué)前緣，17（6）： 271-277.

李巧靈，雷曉東，李晨，2018. 抗干擾編碼電法在通州深部巖溶發(fā)育區(qū)地質(zhì)構(gòu)造探測中的應(yīng)用[J]. 工程勘察，（2）： 71-78.

廖海軍，馬志飛，2015. 北京市海淀區(qū)普安店地區(qū)巖溶塌陷災(zāi)害成因分析及防治對策[J]. 城市地質(zhì)，10（S1）： 218-221.

石建省，1996. 唐山市活動性斷裂與巖溶塌陷關(guān)系的空間分析[J]. 中國地質(zhì)災(zāi)害與防治學(xué)報，（3）： 87-91.

石亞雄，孫少安，吳麗華，等，1991. 微重力測量在溶洞探測中的應(yīng)用[J]. 物探與化探，15（6）： 468-470.

王謙身，周文虎，武傳真，等，1995. 微重力方法在考古工程中的應(yīng)用──明茂陵地下陵殿探查[J]. 地球物理學(xué)進展，10（2）： 85-94.

王延濤，潘瑞林，2012. 微重力法在采空區(qū)勘查中的應(yīng)用[J]. 物探與化探，36（S1）： 61-64.

楊進，2011. 環(huán)境與工程地球物理[M]. 北京： 地質(zhì)出版社.

袁炳強，劉士毅，于國明，等，2015. 秦始皇帝陵遺址考古重磁力場特征[J]. 地球科學(xué)（中國地質(zhì)大學(xué)學(xué)報），（10）： 1616-1620.

趙博超，朱蓓，王弘元，等，2015. 淺談巖溶塌陷的影響因素與模型研究[J]. 中國巖溶，34（5）： 515-521.

鐘立勛，2001. 謹(jǐn)防京郊平原發(fā)生深層隱伏巖溶塌陷[J]. 中國地質(zhì)災(zāi)害與防治學(xué)報，12（4）： 26-29.

Pearson-Grant S C， Franz P， Clearwater J， 2018. Gravity measurements as a calibration tool for geothermal reservoir modelling[J]. Geothermics. 73： 146-157.