張磊 邵永新 馬建英 李恩建 金大利 王靜
1)中國礦業(yè)大學(xué)(北京),北京市海淀區(qū)學(xué)院路丁11號 100083
2)天津市地震局,天津河西區(qū)友誼路19號 300201
3)山東省地震局,山東濟(jì)南 250014
伸縮儀是精密測量地殼巖體兩點間水平距離相對變化的儀器,適用于觀測地殼應(yīng)變和固體潮水平分量的連續(xù)變化,為研究地震孕育過程中水平應(yīng)變的變化規(guī)律提供數(shù)據(jù),也為地球彈性研究提供重要數(shù)據(jù)(熊仲華,2006)。抽水開采是洞體應(yīng)變觀測中的一種重要干擾源(金克儉等,1986),也是伸縮應(yīng)變觀測較為常見的干擾因素(李國斌等,2004;李國斌,2006;汪翠枝等,2009)。已有研究表明(裴曉峰,1995),薊縣地震臺定點形變觀測所受干擾程度最大、發(fā)生最為頻繁的是機(jī)井抽水的干擾。SSY型伸縮儀安裝于薊縣臺東約4.5km處的小辛莊山洞內(nèi),距洞口220~260m處,自2001年5月開始正式觀測,基線長39.7m,NS測向方位角144.87°,EW測向方位角59.94°。小辛莊山洞周邊存在多口用于日常飲水和農(nóng)田灌溉的抽水井,水井距觀測山洞最近約300m,井深一般在100~200m,多分布在山洞的南北兩側(cè),伸縮應(yīng)變受抽水干擾嚴(yán)重,表現(xiàn)為相對拉張變化,尤其以NS向更為敏感。
薊縣地震臺位于首都圈中東部,“燕山沉降帶”的“京東臺凹”的北部,是全國定點地殼形變基本臺。臺站附近斷裂構(gòu)造較為復(fù)雜,存在一系列斷裂和褶皺,臺站附近有薊縣山前斷裂、楊莊斷裂和黃崖關(guān)斷裂。小辛莊山洞山體為震旦系霧迷山組灰?guī)r(含燧石),巖層層面走向335°,傾角20°~33°,測洞南距薊縣斷裂1.5km。儀器室內(nèi)的日溫差小于0.03℃,年溫差小于0.05℃,相對濕度達(dá)98%。桑梓井位于小辛莊山洞西南約37km處的薊縣桑梓,自2004年開始正式觀測,井深300m左右,含水層屬承壓水,其北部為薊縣山前斷裂(圖1)。
圖1 天津地傾斜臺站與斷裂分布
為了研究薊縣臺小辛莊山洞伸縮應(yīng)變受周邊短時地下水抽取的影響,本文利用水位的變化獲取與地下水抽取有關(guān)的信息,但是小辛莊山洞周圍沒有水位觀測井,故文中利用天津地震前兆臺網(wǎng)中薊縣桑梓井水位數(shù)據(jù),通過分析小辛莊伸縮應(yīng)變與桑梓井水位的年動態(tài)變化特征,判斷二者對短時抽取地下水的響應(yīng)是否存在一致性,進(jìn)而建立水位變化量與應(yīng)變變化量的關(guān)系,研究伸縮應(yīng)變受周邊短時抽水干擾變化幅度,以期對在季節(jié)性抽水時段伸縮應(yīng)變出現(xiàn)短期加速拉張變化的異常判定提供幫助。
本文整理了薊縣臺小辛莊山洞伸縮應(yīng)變和薊縣臺桑梓井水位2004年1月1日~2014年6月30日以來日均值數(shù)據(jù)(圖2),由圖可知,二者均表現(xiàn)出了較好的年動態(tài)變化特征:
(1)小辛莊伸縮應(yīng)變自2004年以來NS、EW兩方向年變形態(tài)呈現(xiàn)出一定的規(guī)律性變化。NS向每年年初變化基本平穩(wěn),4月前后受季節(jié)性抽水的影響呈現(xiàn)拉張變化,并在此基礎(chǔ)上波動變化,6~8月份雨季到來受降雨影響呈現(xiàn)壓縮變化,此后逐漸恢復(fù)呈相對拉張變化;EW向每年年初整體上呈現(xiàn)壓縮變化,4月前后受地下水抽取影響呈現(xiàn)相對拉張變化,之后受降雨影響呈現(xiàn)壓縮變化,此后逐漸恢復(fù)呈現(xiàn)相對拉張變化。
(2)桑梓水位自2004年開始觀測以來具有較好的年周期變化特征。每年3月前后為水位的高值期,之后轉(zhuǎn)降;6月前后降到年度最低值,下降幅度為5m左右,之后轉(zhuǎn)升;到11月前后達(dá)到年度第2次高水位狀態(tài),上升幅度為3m左右,之后再次下降;到12月初達(dá)到年度第2次低水位期,下降幅度為1m左右,之后又開始上升,一直到轉(zhuǎn)年的高水位期。
圖2 薊縣臺桑梓井水位和小辛莊伸縮應(yīng)變NS向、EW向日均值曲線
圖3表明,薊縣臺小辛莊伸縮應(yīng)變觀測對抽水響應(yīng)敏感,尤其是NS向更為突出。受季節(jié)性抽水影響,薊縣臺小辛莊伸縮應(yīng)變在4月前后加速拉張變化顯著(圖4),水位在3月前后持續(xù)高值平穩(wěn)變化后轉(zhuǎn)降。對比2004年以來歷年年初桑梓井水位與小辛莊伸縮應(yīng)變數(shù)據(jù)(圖4)發(fā)現(xiàn),伸縮兩分向加速拉張與水位加速下降在時間上具有很好的同步性。雖然薊縣臺小辛莊山洞與桑梓井空間上還存在一定的距離,但二者對季節(jié)性抽水的響應(yīng)表現(xiàn)出很好的同步性變化。楊耀棟等(2011)發(fā)現(xiàn),天津平原區(qū)深層地下水受開采影響表現(xiàn)為徑流-越流-開采型動態(tài)特征,因受灌溉強(qiáng)開采影響,低水位多出現(xiàn)于5~6月份,有時出現(xiàn)在9~10月份;豐水期8月停采后,水位逐漸回升,有時出現(xiàn)高水位,大多至翌年1~3月水位最高,高水位期較最大降水期滯后5~7個月,年水位變幅4~6m。桑梓水位自2004年觀測以來的年變化形態(tài)與上述區(qū)域深層地下水動態(tài)變化特征一致。另外,陳鵬等(2004)在薊縣地震臺開展抽水實驗研究發(fā)現(xiàn),小辛莊伸縮應(yīng)變觀測受抽水干擾的幅度既與抽水井到形變監(jiān)測儀器距離有關(guān),也與抽水井水位的最大降深有關(guān)。
薊縣臺小辛莊伸縮應(yīng)變加速拉張變化與桑梓井水位加速下降變化同步性良好,本文利用最小二乘線性回歸分析方法研究二者之間的關(guān)系,建立伸縮應(yīng)變變化量與水位變化量之間的線性回歸模型,并利用t檢驗方法對模型進(jìn)行顯著性檢驗,最后將2014年水位實測變化量代入回歸模型計算應(yīng)變量的理論變化區(qū)間,驗證2014年薊縣臺小辛莊伸縮應(yīng)變加速拉張變化與地下水抽取的關(guān)系。
對于給定數(shù)據(jù)點 (xi,yi),i=1,2…n,滿足如下關(guān)系
圖3 2012年3月12~13日薊縣臺小辛莊伸縮應(yīng)變NS、EW向受抽水干擾數(shù)據(jù)變化
圖4 薊縣臺小辛莊伸縮應(yīng)變NS向、EW向與多年桑梓井水位日均值曲線的對比
其最小二乘回歸直線可以表示為
方差σ2的無偏估計
相關(guān)系數(shù)為
利用t檢驗方法進(jìn)行回歸模型顯著性檢驗。假設(shè):H0:β1=0,H1:β1≠ 0,當(dāng)H0成立時
對于給定顯著性水平α,記t為統(tǒng)計量(7)的觀測值,若
則拒絕H0:β1=0,即認(rèn)為自變量x對因變量y有顯著性影響。
給定顯著性水平α,在x=x0處,標(biāo)準(zhǔn)差為
2.2.1 回歸模型的建立
由于小辛莊伸縮應(yīng)變兩方向加速拉張與桑梓井水位加速下降變化在時間上具有很好的同步性,本文統(tǒng)計了2004~2013年年初桑梓水位最高時水位和伸縮應(yīng)變NS向、EW向觀測值以及伸縮應(yīng)變兩方向加速拉張變化結(jié)束時的應(yīng)變最大值和水位值,計算水位最高至伸縮兩方向拉張結(jié)束時段內(nèi)水位和應(yīng)變觀測的變化量(表1、表2)。
根據(jù)表1、表2可以得到水位變化量與伸縮應(yīng)變變化量的對應(yīng)關(guān)系(表3)。表3中xNi、yNi分別為研究水位變化量與伸縮NS向應(yīng)變變化關(guān)系的樣本點,xNi為薊縣臺桑梓井水位變化量,yNi為薊縣臺小辛莊伸縮NS向應(yīng)變變化量;xEi、yEi為研究水位變化量與伸縮 EW向應(yīng)變變化關(guān)系的樣本點,xEi為薊縣臺桑梓井水位變化量,yEi為薊縣臺小辛莊伸縮EW向應(yīng)變變化量。xNi,yNi與xEi、yEi的分布如圖5所示,兩組數(shù)據(jù)的相關(guān)系數(shù)分別為 0.8342、0.5884,可見薊縣臺桑梓井水位變化量與小辛莊伸縮應(yīng)變變化量相關(guān)性良好。
表1 薊縣臺桑梓井水位最高值和小辛莊伸縮應(yīng)變NS向拉張最大值時水位與應(yīng)變的變化
表2 薊縣臺桑梓井水位最高值和小辛莊伸縮應(yīng)變EW向拉張最大值時水位與應(yīng)變的變化
表3 桑梓井水位變化量與小辛莊伸縮應(yīng)變變化量的對應(yīng)關(guān)系
由表3數(shù)據(jù),利用2.1中的最小二乘線性回歸方法,根據(jù)式(2~4)計算得到最小二乘回歸直線相關(guān)參數(shù),如表4所示。
表4 最小二乘線性回歸結(jié)果參數(shù)表
由表4可得薊縣臺桑梓井水位變化量與小辛莊伸縮NS向、EW向應(yīng)變變化量之間的最小二乘回歸直線(圖5)表示為
相應(yīng)的式(10)回歸直線方差的無偏估計分別為
圖5 薊縣臺小辛莊伸縮應(yīng)變變化量與桑梓井水位變化量線性擬合分析結(jié)果
2.2.2 回歸模型顯著性檢驗
本文利用t檢驗方法對回歸模型(10)進(jìn)行顯著性檢驗,依據(jù)式(7~9)分別用表示回歸模型式(10)中直線xE的t統(tǒng)計量的觀測值。則有:
(1)對于薊縣臺桑梓井水位變化量與小辛莊伸縮NS向應(yīng)變變化量的回歸直線:44.6858+328.7437xN,有n=10,假設(shè):H0:β1=0;H1:β1≠0,對于給定顯著性水平α=0.05,t0.025(n-2)=t0.025(8)=2.3060,t統(tǒng)計量的觀測值,則拒絕H0,認(rèn)為自變量xN對因變量yN有顯著性影響即水位變化對伸縮NS向應(yīng)變變化影響顯著,置信水平達(dá)到了95%。
(2)對于薊縣臺桑梓井水位變化量與小辛莊伸縮EW向應(yīng)變變化量的回歸直線:有n=10,給定顯著性2.0582<t0.025(8),不滿足顯著性水平為α=0.05的假設(shè)檢驗;于是取α=0.1,有t0.05(n-2)=1.8595,),則在給定顯著性水平α=0.1條件下,自變量xE對因變量yN有顯著性影響即水位變化對伸縮EW向應(yīng)變變化影響顯著,置信水平達(dá)到了90%。
上述t檢驗結(jié)果表明,水位變化對伸縮應(yīng)變變化影響顯著,并且水位變化量與應(yīng)變變化量的回歸直線方程N(yùn)S向置信水平高于EW向,這與實際觀測中伸縮應(yīng)變NS向?qū)Τ樗憫?yīng)更為敏感這一實測現(xiàn)象也是一致的。
2.2.3 應(yīng)變量變化區(qū)間分析
薊縣臺小辛莊伸縮應(yīng)變觀測受抽水干擾影響嚴(yán)重,一般在每年4月前后開始加速拉張變化(圖6),但2014年伸縮應(yīng)變兩方向自3月中旬開始加速拉張變化,時間較往年同期有所提前,且變化較為突出(圖7)。震例研究表明(汪翠枝等,2009),薊縣伸縮應(yīng)變震前一般表現(xiàn)為短臨前兆異常,NS向在文安地震前、EW向在盧龍地震前均出現(xiàn)短期加速拉張變化。針對這種異常變化,首先對比小辛莊伸縮應(yīng)變和桑梓井水位日均值曲線,發(fā)現(xiàn)水位加速轉(zhuǎn)降與伸縮兩方向加速拉張變化開始的時間基本一致。另外,調(diào)查發(fā)現(xiàn),此期間儀器工作正常,周邊也并未發(fā)現(xiàn)新的干擾源,故初步判定此異常變化可能與季節(jié)性抽水有關(guān)。
圖6 薊縣臺小辛莊伸縮應(yīng)變NS、EW向日均值多年對比曲線
本文依據(jù)表1、表2,統(tǒng)計了2014年水位最高時水位和伸縮NS向、EW向觀測值、伸縮應(yīng)變兩方向加速拉張變化結(jié)束時的應(yīng)變最大值和水位值,分別計算出2014年該時段內(nèi)薊縣臺桑梓井水位和小辛莊伸縮 NS向、EW向的實測應(yīng)變變化量(表5),分別用(xN0,yN0)和(xE0,yE0)表示,相應(yīng)的標(biāo)準(zhǔn)差分別用δN和δE表示。如表5所示,依據(jù)式(10)和式(9)有:
(1)當(dāng)薊縣臺桑梓井水位變化量xN0=2.8841時,代入式得到2014年伸縮NS向應(yīng)變變化量的理論估計值標(biāo)準(zhǔn)差δN=408.8522,2014年薊縣臺小辛莊伸縮應(yīng)變NS向的實測變化量yN0=674.4375,落入?yún)^(qū)間即(583.9633,1401.6677)內(nèi),置信水平為 95%;
圖7 2014年薊縣臺小辛莊伸縮應(yīng)變NS向、EW向與桑梓井水位日均值對比曲線
表5 2014年桑梓井水位最高值和小辛莊伸縮應(yīng)變拉張最大值時水位與應(yīng)變變化
(2)當(dāng)薊縣臺桑梓井水位變化量xE0=2.1929時,代入式得到2014年伸縮EW向應(yīng)變變化量的理論估計值、標(biāo)準(zhǔn)差δE=472.3435,2014年薊縣臺小辛莊伸縮EW向的實測應(yīng)變變化量yE0=481.05,落入?yún)^(qū)間即(-93.7461,850.9409)內(nèi),置信水平為 90%。
對于給定任意xN0和xE0,yN0和yE0的變化區(qū)間理論上下限曲線如圖5中虛線所示。2014年薊縣臺桑梓井水位的實測變化量和小辛莊伸縮NS向、EW向的實測應(yīng)變變化量(2.8841,674.4375)、(2.1929,481.05)在理論控制線內(nèi)。因此,伸縮應(yīng)變在2014年3月中旬以來出現(xiàn)的加速拉張變化與季節(jié)性抽水有關(guān),作為地震前兆異常的信度不高。
地下水抽取過程實質(zhì)上為以井孔為中心的含水層的疏干過程,由于含水層失水產(chǎn)生地表沉降,形成以井孔為中心的局部收縮區(qū)域,受其影響井孔中心收縮區(qū)以外受到拉張作用,這就解釋了小辛莊山洞附近短時抽水對伸縮應(yīng)變的拉張影響。Terzaghi于1925年提出有效應(yīng)力原理(王大純等,2002)可以解釋抽水引起地下水位動態(tài)變化進(jìn)而引發(fā)地面沉降的問題。有效應(yīng)力即為實際作用于砂層骨架上的應(yīng)力,等于總應(yīng)力減去孔隙水壓力。抽水造成地下水水位下降,孔隙水壓力降低,有效應(yīng)力增加,砂層顆粒發(fā)生位移,顆粒接觸面積增加,排列更趨緊密,孔隙度降低,砂層受到壓密,從而產(chǎn)生地面沉降,含水層的孔隙度、給水度、滲透系數(shù)等參數(shù)均將變小。砂礫類巖土基本呈彈性形變,待砂層孔隙水壓恢復(fù)后,砂層能大體恢復(fù)原狀。
已有水文資料表明薊縣臺小辛莊山洞附近地下水類型屬于潛水,利用地下水潛水完整井穩(wěn)定流Dupuit公式可以得到抽水井影響半徑,該區(qū)域旱季和雨季井水位降深對形變干擾的影響范圍在800~1000m(陳鵬等,2004;李恩建等,2006)。楊耀棟等(2011)的研究也表明,天津地區(qū)抽水開采是影響深層地下水位動態(tài)的主要因素,地下水開采與全淡水區(qū)的淺層地下水位動態(tài)具有良好的相關(guān)性。文中回歸分析結(jié)果也表明,桑梓水位的短時下降變化與伸縮應(yīng)變的短時加速拉張變化相關(guān)性良好,二者之間的回歸模型置信水平達(dá)到了90%,同時也說明了利用桑梓水位作為區(qū)域地下水動態(tài)背景信息研究伸縮應(yīng)變受抽水影響的短時變化是具有實際意義的。
薊縣小辛莊伸縮應(yīng)變觀測受抽水干擾影響顯著,每年年初受季節(jié)性抽水影響NS、EW兩方向會表現(xiàn)出不同程度的短期加速拉張變化。桑梓水位多年來年動態(tài)特征良好且與區(qū)域深層地下水動態(tài)特征一致,在季節(jié)性抽水時段,水位加速下降與伸縮應(yīng)變加速拉張變化時間上具有很好的同步性。本文利用天津地震前兆臺網(wǎng)的觀測數(shù)據(jù)開展了地下水抽取與伸縮應(yīng)變短時變化關(guān)系的研究,得到如下結(jié)果:
(1)文中統(tǒng)計并計算了2004~2013年水位最高值至伸縮應(yīng)變拉張結(jié)束時水位變化量和應(yīng)變的變化量,利用最小二乘線性回歸方法得到了水位變化量與應(yīng)變變化量的回歸直線(式(10)),二者呈正相關(guān),即伸縮應(yīng)變受抽水影響的程度隨著水位變化量的增大而增大,這與陳鵬等(2004)抽水實驗的研究結(jié)果也是一致的。
(2)水位變化量與伸縮NS向和EW向應(yīng)變變化量的線性回歸模型分別滿足顯著性水平為α=0.05和α=0.1的t檢驗,即回歸模型的置信水平分別為95%和90%。回歸模型置信水平較高,水位變化量與NS向應(yīng)變變化量的回歸模型置信水平更高,這與實際觀測過程中伸縮應(yīng)變NS向?qū)Τ樗憫?yīng)更為敏感的客觀事實一致。
(3)伸縮應(yīng)變兩方向在2014年3月中旬開始加速拉張變化,考量年變形態(tài)發(fā)現(xiàn),時間較往年有所提前,且變化較為突出。而調(diào)查中并未發(fā)現(xiàn)新的干擾源,且儀器工作正常。利用文中的回歸分析方法計算得到2014年伸縮兩方向應(yīng)變變化量的理論變化區(qū)間,并且2014年應(yīng)變的實測變化量在變化區(qū)間內(nèi),證明了該變化為季節(jié)性抽水所致,作為地震前兆異常的信度較低。
桑梓水位反映了區(qū)域地下水的動態(tài)背景,在季節(jié)性抽水時段桑梓水位變化量與小辛莊伸縮應(yīng)變變化量的回歸模型置信水平達(dá)到了90%以上,回歸模型顯示伸縮應(yīng)變變化受水位變化影響顯著。由于小辛莊山洞附近缺乏相應(yīng)的水位觀測,因此可以利用桑梓水位作為輔助資料研究伸縮應(yīng)變受季節(jié)性抽水短時干擾的變化幅度。隨著觀測資料的積累,樣本量逐漸增加,回歸模型可以不斷進(jìn)行修正并進(jìn)行實例驗證,這有利于輔助判定在季節(jié)性抽水時段伸縮應(yīng)變出現(xiàn)短期加速拉張變化的異常信度。另外,在受抽水干擾影響顯著的形變觀測臺站,可以開展水位觀測項目,以輔助研究抽水對地形變觀測的短時影響。