李衛(wèi)衛(wèi), 熊 鑫, 蒙愛軍
(1.浙江華東建設工程有限公司,浙江 杭州 311202;2.浙江省工程物探勘察設計院有限公司,浙江 杭州 315012)
很多建筑物的持力層都落在基巖上,若基巖面附近存在巖溶則是危害建筑物基礎的重大隱患。高密度電法[1-5]、跨孔地震 CT[6]、地質(zhì)雷達[7]以及綜合物探方法[8,9]廣泛應用于巖溶探測。當基巖上部有覆蓋層時,借助于勘探孔進行地震跨孔CT探測巖溶是一種很好的選擇。談順佳等[10]和楊永龍等[11]在對一般建筑物基礎的巖溶分布探測中采用了跨孔地震CT技術,取得了較好的效果。段成龍等[12]和郭書蘭等[13]在地鐵基礎巖溶探測中采用了跨孔地震CT技術,根據(jù)成像圖中低速帶范圍及波形分析確定了溶洞位置和充填情況。朱文仲等[14]對彈性波CT技術幾個重要問題進行了正反演研究,就彈性波 CT技術的分辨率問題、“低速色盲”問題和像元劃分問題,提出了各種觀測系統(tǒng)的最高分辨率估算公式,指出了二側透視觀測系統(tǒng)的重大缺陷、彈性波CT技術的“低速色盲”本質(zhì)以及劃分像元的一般原則。石振明等[15]基于波動方程的時域有限差分對跨孔地震巖溶探測進行正演模擬,從模擬的地震記錄中提取初至旅行時,最后以基于程函方程的初至層析方法反演得到鉆孔間的地層速度模型。潘紀順等[16]利用 Matlab 編程對堤防工程進行地震斷面CT數(shù)值模擬,通過ART和SIRT重建方法進行反演,分析得出二者在堤防工程CT應用中的優(yōu)缺點及適用條件。
為了更好、更精確地探測巖溶,本文從基巖面附近的多巖溶探測入手,依據(jù)射線理論進行正反演模擬研究,討論巖溶與基巖面的距離、圍巖與巖溶速度差異大小以及觀測系統(tǒng)對反演結果的影響。本文還對將巖土介質(zhì)一起加入反演的情況進行研究。
本文采用韓佩恩等[17]提出的改進的Moser方法[18]進行初至時間場計算,利用Vidale[19]和Asakawa和Kawanaka[20]提出的走時線性插值算法確定射線路徑。
韓佩恩等在Moser算法的基礎上,提出了加入動態(tài)節(jié)點在同樣的網(wǎng)格間距下提高計算初至時間精度的方法,如圖1所示,在加入動態(tài)節(jié)點情況下采用Moser算法計算所有節(jié)點的初至旅行時,滿足根據(jù)已知網(wǎng)格化的速度模型正演計算激發(fā)點至接收點旅行時的要求。Moser算法難以提供反演要求的準確射線路徑,走時線性插值算法則可以有效解決這個問題,該算法分為兩步:第一步是向前算法,即計算從震源出發(fā)到達網(wǎng)格各個節(jié)點的地震波初至旅行時;第二步是向后算法,即根據(jù)各個節(jié)點的旅行時和接收點位置反向計算接收點至震源的確切射線路徑。具體計算過程參考相關文獻[19,20]。
圖1 改進的Moser方法示意圖Fig.1 Schematic diagram of the improved Moser method
反演采用聯(lián)合迭代重建法(SIRT),其基本原理為:將所有射線的投影差值計算出來,然后將各個射線的投影差值反投影到單元網(wǎng)格上;將網(wǎng)格單元排序,把所有穿過網(wǎng)格單元的射線的走時投影匯總起來,利用這些匯總數(shù)據(jù)對網(wǎng)格單元地震波慢度(速度的倒數(shù))進行修正,直到獲得的結果滿足預設的精度或者達到設定的最大迭代次數(shù)。在做跨孔地震波CT過程中,網(wǎng)格化所需要勘察的區(qū)域,根據(jù)先驗信息給網(wǎng)格化后的每個網(wǎng)格單元賦予一個慢度初值。再計算每條經(jīng)過網(wǎng)格的初至波射線的走時差,然后將走時差分配到每個網(wǎng)格單元中,并對所有分配到的走時差進行統(tǒng)計,依據(jù)平均值來修正網(wǎng)格單元的地震波慢度。
設lkj為第k條射線穿過第j個網(wǎng)格(按照一維編序號)的射線長度,設Lk為第k條射線長度,即
(1)
式(1)中Mk為第k條射線所經(jīng)過的模型網(wǎng)格個數(shù)。
設Δtk為第k條射線的觀測值與理論旅行時之差,設βk為第k條射線的模型修正系數(shù),即
(2)
則第k條射線對第j個網(wǎng)格的慢度修正量為(μ為微小的量):
(3)
所以,所有射線對第j個網(wǎng)格的慢度修正量為:
(4)
式(4)中Nj表示第j個網(wǎng)格所穿過的射線條數(shù)。
模型慢度修正(反演)為:
(5)
其中q為迭代次數(shù)。
本文理論模擬的研究思路是建立多巖溶不同排列組合的正方形網(wǎng)格化孔間地震縱波速度剖面模型,根據(jù)模型設計跨孔地震波初至CT成像技術的觀測系統(tǒng),如圖2所示,采用改進的Moser方法[17]計算所有炮點(激發(fā)點)至檢波點(接收點)的初至旅行時,以此觀測系統(tǒng)和初至旅行時數(shù)據(jù)反演兩孔之間的速度分布。反演的步驟是:①根據(jù)觀測系統(tǒng)和初至旅行時建立初始反演模型;②采用改進的Moser方法[18]計算所有炮檢點的初至旅行時;③采用Vidale[19]提出的走時線性插值向后算法計算所有炮檢點的射線路徑;④修改每個網(wǎng)格的速度(慢度)值;⑤判斷是否滿足反演終止條件;⑥如果不滿足終止條件,則繼續(xù)②~⑤的步驟。
圖2 跨孔地震CT觀測系統(tǒng)示意圖Fig.2 Schematic diagram of cross-hole seismic CT observation system
本文建立初始模型的思路如下:設Ns個激發(fā)點,其縱坐標為ZSi(i=1, 2, 3,...,Ns),與其對應縱坐標的接收點為ZPj(j=1, 2, 3,...,Np),地震波初至的最小旅行時為tij,假設能夠找到ZPk=ZSi,將其簡稱為Zi(i=1, 2, 3,...,Ns),假設鉆孔間距為dx,按照直射線傳播假設,則傳播速度為V(Zi)=dx/tik(i=1, 2, 3,...,Ns),用速度序列V(Zi)進行垂向線性插值,建立橫向均勻垂向連續(xù)介質(zhì)模型。
本次理論模擬的孔間距固定為15 m,孔深為20 m,左邊為激發(fā)孔,從孔口到孔底1 m一個激發(fā)點,右邊為接收孔,從孔口到孔底1 m一個接收點,共計21個接收點。觀測系統(tǒng)采用60度張望角,排列長度(接收段長度)L=2×15 m×tan60°≈52 m>2×20 m,所以21個接收點均可接收到每一個激發(fā)點出發(fā)的地震波。模型的網(wǎng)格化間距為0.25 m×0.25 m,反演模型的網(wǎng)格化間距均為0.5 m×0.5 m,速度剖面的色標值是動態(tài)的,未統(tǒng)一在同一個值域,通常情況下,圍巖的速度設為3 000 m/s,巖溶填充物的速度為1 500 m/s。
設計一組兩個水平相鄰的巖溶模型,半徑均為R=1.5 m,巖溶中心間距均為L=4 m,巖溶中心埋深(與基巖面的距離)分別為H=4 m、5.5 m、7 m。正反演結果如圖3~圖5所示??梢钥闯觯寒攷r溶中心與基巖面的垂直距離不大于5.5 m時,反演的結果只能看到一個巖溶體,埋深較準,水平位置在兩個巖溶中間(圖3d,圖4d);當巖溶中心與基巖面的垂直距離達到7 m時(圖5d),則看到兩個巖溶體的準確中心位置,緊靠接收孔的假異常則可以通過接收鉆孔排除其為巖溶體可能。當巖溶中心與基巖面的垂直距離為7 m半徑縮為1.25 m時(圖6d),也能探測到兩個巖溶體的存在。以上研究表明,當兩個巖溶體距離基巖面太近時,無法解析多巖溶體的空間分布;當兩個巖溶體與基巖面相隔一定距離時,則可以分辨多巖溶體的空間位置。
圖3 模型1正反演結果(H=4 m, L=4 m, R=1.5 m)Fig.3 The forward and inverse result for model 1(H=4 m, L=4 m, R=1.5 m)
圖4 模型2正反演結果(H=5.5 m, L=4 m, R=1.5 m)Fig.4 The forward and inverse result for model 2 (H=5.5 m, L=4 m, R=1.5 m)
圖5 模型3正反演結果(H=7 m, L=4 m, R=1.5 m)Fig.5 The forward and inverse result for model 3 (H=7 m, L=4 m, R=1.5 m)
圖6 模型4正反演結果(H=7 m, L=4 m, R=1.25 m)Fig.6 The forward and inverse result for model 4 (H=7 m, L=4 m, R=1.25 m)
為了探討跨孔地震CT的橫向分辨能力,本文設計了兩類模型,一是兩個水平相鄰的巖溶,半徑均為R=1.5 m,巖溶中心間距L和中心埋深H可變;二是一個水平矩形條狀模型,高度h=1.5 m,其中心埋深H和寬度w可變。圖7是將圖3(a)中的巖溶中心間距從4 m增加到5 m的正反演結果,能明顯看到兩個巖溶體的存在,右邊的巖溶體位置誤差較大;如圖8所示,將圖7(a)的巖溶中心埋深從4 m減到2.5 m時,已經(jīng)探測不到兩個巖溶體的準確位置;圖9是將圖8(a)中的巖溶中心間距從5 m增加到9 m的正反演結果,能明顯地看到兩個巖溶體的存在,位置基本正確。以上研究表明:如果巖溶體距離基巖面較近,當巖溶體的中心距增加到一定程度時,跨孔地震CT也可以分辨多巖溶體的存在。
圖7 模型5正反演結果(H=4 m, L=5 m, R=1.5 m)Fig.7 The forward and inverse result for model 5 (H=4 m, L=5 m, R=1.5 m)
圖8 模型6正反演結果(H=2.5 m, L=5 m, R=1.5 m)Fig.8 The forward and inverse result for model 6 (H=2.5 m, L=5 m, R=1.5 m)
圖9 模型7正反演結果(H=2.5 m, L=9 m, R=1.5 m)Fig.9 The forward and inverse result for model 7 (H=2.5 m, L=9 m, R=1.5 m)
如圖10(a)和圖11(a)所示,水平矩形條巖溶模型的寬度均為w=8 m,中心埋深分別為H=5 m、H=10 m,反演結果如圖10(d)和圖11(d)所示,中心埋深反演地都很準確,但左右邊緣位置都相差太大。圖12是將圖10(a)的水平矩形條巖溶模型的寬度增至15 m(橫向貫穿)時的正反演結果,可以看出跨孔地震CT能夠?qū)M向貫穿性巖溶進行準確探測。以上研究表明:采用跨孔地震CT技術對條帶巖溶進行探測,當條帶巖溶水平方向貫穿到兩個鉆孔時該方法可以準確探測,而對于水平方向非貫穿巖溶則只能探測其中心埋深而不能探測其左右邊緣。
圖10 模型8正反演結果(H=5 m, w=8 m, h=1.5 m)Fig.10 The forward and inverse result for model 8 (H=5 m, w=8 m, h=1.5 m)
圖11 模型9正反演結果(H=10 m, w=8 m, h=1.5 m)Fig.11 The forward and inverse result for model 9 (H=10 m, w=8 m, h=1.5 m)
圖12 模型10正反演結果(H=5 m, w=15 m, h=1.5 m)Fig.12 The forward and inverse result for model 10 (H=5 m, w=15 m, h=1.5 m)
本文設計一組兩個水平相鄰的巖溶體,半徑均為R=1.5m,巖溶中心間距均為L=5 m,巖溶中心埋深為H=2.5 m,一共三個模型,通過調(diào)整圍巖與巖溶填充物速度大小關系,探討圍巖與巖溶充填物速度差異大小對跨孔地震CT反演的影響。圖13和圖14反演結果一定程度上恢復了兩個巖溶的空間位置,但圖15中兩個巖溶體則無法分辨出來。巖溶充填物與圍巖的速度差相對于圍巖速度分別為(2 500 m/s-1 500 m/s)/2 500 m/s=0.4、(2 000 m/s-1 500 m/s)/2 000=0.25 m/s、(2 500 m/s-2000 m/s)/2 500 m/s=0.2,依據(jù)反演結果可以看出,巖溶充填物與圍巖的速度差與圍巖速度的比值越大,跨孔地震CT技術越易于分辨兩個水平相鄰的巖溶體。
圖13 模型11正反演結果(H=2.5 m, L=5 m, R=1.5 m)Fig.13 The forward and inverse result for model 11 (H=2.5 m, L=5 m, R=1.5 m)
圖14 模型12正反演結果(H=2.5 m, L=5 m, R=1.5 m)Fig.14 The forward and inverse result for model 12 (H=2.5 m, L=5 m, R=1.5 m)
圖15 模型13正反演結果(H=2.5 m, L=5 m, R=1.5 m)Fig.15 The forward and inverse result for model 13 (H=2.5 m, L=5 m, R=1.5 m)
圖16 模型14正反演結果(H=2.5 m, L=5 m, R=1.5 m)Fig.16 The forward and inverse result for model 14 (H=2.5 m, L=5 m, R=1.5 m)
圖17 模型15正反演結果(H=4 m, L=5 m, R=1.5 m)Fig.17 The forward and inverse result for model 15 (H=4 m, L=5 m, R=1.5 m)
前文研究表明,當兩個巖溶體中心間距較小且又接近于基巖面時,在基巖中布置觀測系統(tǒng)無法分辨這兩個巖溶體?;诖耍疚脑O計含覆蓋層的模型,探討在覆蓋層中也布置激發(fā)點和接收點能否探測清楚基巖面附近的多個(大于等于2)巖溶體。參考前人研究結果[11, 13],本文將覆蓋層的速度參數(shù)亦作為反演參數(shù)。
圖18 模型16正反演結果(H=5.5 m, L=5 m, R=1.5 m)Fig.18 The forward and inverse result for model 16 (H=5.5 m, L=5 m, R=1.5 m)
圖19 模型17正反演結果Fig.19 The forward and inverse result for model 17
本文設計了4個模型,覆蓋層厚度都是5 m,前三個模型都是基巖里含兩個水平相鄰的巖溶體,半徑均為R=1.5 m,巖溶中心間距均為L=5 m,巖溶中心埋深(與基巖面的距離)分別為H=2.5 m、4 m、5.5 m,第4個模型是基巖中不含巖溶?;趫D16(a)、圖17(a)、圖18(a)模型做正演,它們的反演結果如圖16(d)、圖17(d)、圖18(d)所示,這3張速度剖面中都看不到兩個巖溶體的存在,它們都與圖19(d)相似,即都在基巖面下方約5 m處模型水平中心位置有一個低速假異常體。以上研究表明:如果將激發(fā)點和接收點也布置在覆蓋層中,那么該觀測系統(tǒng)嚴重影響了反演的結果;若將激發(fā)點和接收點布置在基巖中(圖7d)則能準確探測到兩個巖溶體。
本文對基巖面附近的多巖溶體開展跨孔地震CT正反演研究,得到如下結論:
1)當所有激發(fā)點和接收點都布置在基巖中時,水平排列的兩個巖溶體與基巖面相隔一定距離時可以分辨其空間位置;
2)當巖溶體距離基巖面較近時,兩巖溶體的中心距增加到一定程度,也可以分辨多巖溶體的存在;
3)當條帶巖溶體水平方向貫穿到兩個鉆孔時該方法可以準確探測;
4)巖溶充填物與圍巖的速度差越大,跨孔地震CT技術越易于分辨兩個水平相鄰的巖溶體。
由于控制因素較多,建議在做實際跨孔地震初至CT探測巖溶工作時,根據(jù)實際孔間距、孔深以及初步解釋的巖溶建模,進行正反演研究,以便指導實測反演速度剖面的解釋工作。