趙 頔，史浙明，王廣才，周平根，孫小龍，丁 謀

(1.中國(guó)地質(zhì)大學(xué)(北京)生物地質(zhì)與環(huán)境地質(zhì)國(guó)家重點(diǎn)實(shí)驗(yàn)室， 北京 100083；2.中國(guó)地質(zhì)大學(xué)(北京)地下水循環(huán)與環(huán)境演化教育部重點(diǎn)實(shí)驗(yàn)室， 北京 100083；3.中國(guó)地質(zhì)大學(xué)(北京)水資源與環(huán)境學(xué)院， 北京 100083；4.中國(guó)地質(zhì)環(huán)境監(jiān)測(cè)院， 北京 100081；5.中國(guó)地震局地殼應(yīng)力研究所(地殼動(dòng)力學(xué)重點(diǎn)實(shí)驗(yàn)室)， 北京 100085；6. 天津城建大學(xué)環(huán)境與市政工程學(xué)院， 天津 300384)

地震是一種極常見的自然地理災(zāi)害，多發(fā)于板塊邊緣和板塊內(nèi)部的構(gòu)造活動(dòng)帶，由于能量積累而產(chǎn)生，可以改變地質(zhì)構(gòu)造和地下水的運(yùn)動(dòng)狀態(tài)。地殼中普遍存在的地下水，形成于一定封閉條件承壓含水層時(shí)，可以起到靈敏測(cè)壓計(jì)的作用，能夠放大井-含水層系統(tǒng)壓力波動(dòng)的影響，客觀地反映地殼中應(yīng)力應(yīng)變狀態(tài)的變化情況。一旦受到外力干擾，含水層孔隙的收縮和膨脹能夠引起含水層孔隙壓力ΔP的變化，并通過滲透水流傳遞靜水壓力以井水位動(dòng)態(tài)表現(xiàn)。井水微動(dòng)態(tài)的改變相對(duì)迅速，微小的應(yīng)力應(yīng)變改變都會(huì)在封閉性較好的承壓含水層井孔位上清晰、明顯地反應(yīng)[1]。通過監(jiān)測(cè)地震地下流體動(dòng)態(tài)變化獲得水文地質(zhì)、物理應(yīng)力應(yīng)變和水化學(xué)等多方面的地球參數(shù)，對(duì)于地震的預(yù)報(bào)工作尤其對(duì)地殼活動(dòng)特征、地殼內(nèi)孔隙壓力和體應(yīng)力、偏應(yīng)力等的應(yīng)力分布、地下水徑流以及含水層性質(zhì)等方面的研究提供了非常重要的依據(jù)。對(duì)于許多地質(zhì)工程應(yīng)用來說，理解地震導(dǎo)致的地下水位變化機(jī)制對(duì)于地下水資源評(píng)價(jià)與管理，滑坡的穩(wěn)定性評(píng)估，采礦設(shè)施的穩(wěn)定性，評(píng)價(jià)水庫(kù)誘發(fā)地震等方面至關(guān)重要[2]。

井水位對(duì)氣壓，固體潮和地震波的響應(yīng)表明井-含水層系統(tǒng)對(duì)地殼應(yīng)力具有靈敏的響應(yīng)，同時(shí)這些響應(yīng)隱含著含水層的重要信息，提供了獲取含水層特性的方法[3～4]，因此很多學(xué)者利用地震波與水震波的關(guān)系估算巖石的水力性質(zhì)。Cooper等探究了敞口自流井水位波動(dòng)與地震波的幅頻特征，從理論推導(dǎo)和頻譜分析兩方面研究含水層特性[5]；Kano研究了封閉井中水位對(duì)地震波響應(yīng)的頻率范圍，提供了利用水震波與地震波關(guān)系估算Skempton常數(shù)的方法[1]；晏銳等分析了水位與體應(yīng)變之間的變化關(guān)系，計(jì)算出水位與體應(yīng)變之間的敏感系數(shù)[6]。在以往的研究中，許多地震工作者嘗試應(yīng)用線性孔隙彈性理論解釋井-含水層系統(tǒng)中井水位與應(yīng)力應(yīng)變的相關(guān)關(guān)系，認(rèn)為地震波通過時(shí)含水層的應(yīng)力應(yīng)變變化會(huì)導(dǎo)致孔隙壓力改變，引起井-含水層間的水流運(yùn)動(dòng)，從而導(dǎo)致井水量發(fā)生變化，表現(xiàn)為井水位的升高或降低[1,7]。其中多數(shù)研究認(rèn)為水位只響應(yīng)于體應(yīng)變[1,8]。眾所周知，由于Rayleigh波和P波影響使得含水層擴(kuò)張和收縮，然而最近研究表明由于巖石的各向異性，S和Love波對(duì)井水位也會(huì)產(chǎn)生影響[5,9～10]。前人實(shí)驗(yàn)研究成果表明，體應(yīng)力和偏應(yīng)力會(huì)對(duì)孔隙壓力造成不同程度的影響[11～12]。因此，孔隙壓力的影響因素及巖石線彈性系數(shù)BKu的計(jì)算方法有待進(jìn)一步研究。

本文以線性孔隙彈性理論為基礎(chǔ)，探究了井-含水層系統(tǒng)孔隙壓力變化與體應(yīng)變和偏應(yīng)變之間的數(shù)學(xué)表達(dá)式，通過對(duì)2011年3月11日日本Mw9.0地震引起的昌平地下水觀測(cè)井實(shí)測(cè)水位變化與曲線擬合得到的計(jì)算值對(duì)比，給出了一種求解巖石線彈性系數(shù)BKu的方法。與潮汐分析方法對(duì)比，驗(yàn)證了由地震方法推算的巖石彈性特性代表實(shí)際的含水層特征。

1 昌平井觀測(cè)概況

北京昌平地震臺(tái)位于北京市昌平縣境內(nèi)，處于北東走向的南口-山前活動(dòng)斷裂與北西走向的南口-孫河活動(dòng)斷裂交匯部位以東7.5 km處(圖1)。地表為第四系沉積層，是新構(gòu)造活動(dòng)的復(fù)雜地帶，同時(shí)也是前兆觀測(cè)的靈敏點(diǎn)[13]。

昌平地下水觀測(cè)井為中國(guó)地震局于1986年為了觀測(cè)水位而設(shè)置的鉆井，井孔柱狀圖見圖1。0～14.3 m 深處井孔直徑為170 mm，安裝168 mm套管，14.3～90.1 m深度井孔直徑為150 mm。昌平井含水層位于北京溫榆河沖洪積平原含水層系統(tǒng)的下部，主要由霧迷山組白云質(zhì)灰?guī)r組成，主要含水層厚度為30 m，在觀測(cè)井70～100 m深度鉆核為碎片狀態(tài)。觀測(cè)井由黃褐色砂質(zhì)黏土、殘坡積灰?guī)r和灰黑色致密堅(jiān)硬塊狀硅化灰?guī)r組成，硅化灰?guī)r巖石破碎，節(jié)理發(fā)育，巖石中燧石條帶明顯。按施工設(shè)計(jì)，該井水位的變化為下部含水組的壓力變化導(dǎo)致[6]。所取巖芯均為碎塊，一般約為5 cm。巖層傾斜約45°，在43.35～44.41 m處有約1.06 m厚度溶洞，固結(jié)水泥，由超聲波井下電視圖可以觀測(cè)到57.20～59.70 m處巖石完整，在70～101 m處節(jié)理裂隙及其發(fā)育。

圖1 昌平井構(gòu)造環(huán)境及其井孔結(jié)構(gòu)圖Fig.1 Structural environment and well hole structure of the Changping well

2001年9月我國(guó)地震臺(tái)網(wǎng)完成了數(shù)字化改造，水位觀測(cè)采用數(shù)字化自動(dòng)觀測(cè)，實(shí)現(xiàn)了分鐘值采樣和數(shù)據(jù)自動(dòng)傳輸?shù)墓δ堋?010年投入使用LN-3類型水位傳感器觀測(cè)昌平井水位，測(cè)量精度為±0.2%，穩(wěn)定性±0.25%，儀器分辨率<1 mm，采集數(shù)據(jù)樣品率50 Hz。

2 昌平井井儲(chǔ)效應(yīng)的校正

昌平地下水觀測(cè)井成功記錄到了由日本地震造成的水位波動(dòng)，同時(shí)位于昌平井東部1.2 km的十三陵地震臺(tái)由BBVS-120類型地震儀記錄到了地震波波形的變化(圖2)，震中距為2 283 km。本研究將水位數(shù)據(jù)與3個(gè)方向地震波波形數(shù)據(jù)進(jìn)行對(duì)比分析，當(dāng)P波和S波到達(dá)時(shí)，水位波動(dòng)僅有輕微加強(qiáng)趨勢(shì)，而Love和Rayleigh波作用階段，水位波動(dòng)強(qiáng)烈，因此日本地震引起的水位改變主要受到了面波的影響。地震波傳播方向由東往西，意味著十三陵地震臺(tái)記錄到的地震波東西和垂直方向速度大于南北方向的傳播速度[14]。

圖2 2011-03-11 13:45～14:10 日本Mw9.0地震波速度與測(cè)量井水位響應(yīng)(a)徑向速度；(b)橫向速度；(c)垂向速度；(d)水位變化Fig.2 compare of groundwater levels and three component seismic velocity waveforms recorded at SSL during the period from 13:45 to 14:10 (GMT+8:00) on 11 March 2011 to the Mw9.0 Tohoku earthquake. (a) Radial; (b) Transversal; (c) Vertical component of seismogram;(d) Groundwclter level changes recorded at the Changping well

地震波到達(dá)前13:45～13:50的水位和地震波形特征判斷為非干擾階段，P、S波(13:51～13:56)和Love、Rayleigh波(13:56～14:02)時(shí)間段作為干擾階段。Love波相存在更大的橫向速度，分散性較小，而Rayleigh波存在更大的徑向和垂向速度并且有明顯的發(fā)散，面波這一傳播特點(diǎn)為此次研究提供了科學(xué)依據(jù)和重要參考[15]。本研究分析重點(diǎn)為地震波的體應(yīng)力和偏應(yīng)力對(duì)巖層孔隙壓力造成不同程度的影響。由于P、S和Love、Rayleigh波相特征不同，因此體應(yīng)力和偏應(yīng)力可以分開研究。圖2顯示地震波速度和水位變化之間具有相似性。Shalev等人曾利用位于以色列死海斷裂帶的Meizar1(敞口井)和Meizar2(封閉井)對(duì)2015年Nepal地震引起的水位震蕩進(jìn)行對(duì)比研究[16]，兩口井具有相同含水層條件，相距3 km，由水位的時(shí)間序列可以看出相對(duì)于同一條件下的封閉井，敞口井Meizar1的水位有明顯的衰減，衰減因子為12，這種衰減是由于存在井孔儲(chǔ)積效應(yīng)而形成。本研究的目的是量化水位變化與地震波速度關(guān)系。由于昌平井為敞口自流井，因此孔隙壓力受到井儲(chǔ)效應(yīng)的影響，首先校正井孔儲(chǔ)積效應(yīng)。Liu提供了校正井孔儲(chǔ)積效應(yīng)的方法[17]，Cooper也對(duì)井儲(chǔ)進(jìn)行了研究并提出了計(jì)算公式[18]：

He=Hw+3/8d

(1)

aw=rw(ωS/T)1/2

(2)

ωw=(g/He)1/2

(3)

式中：x0/h0——地震波的周期作用；

Ker、Kei——0階Kelvin作用；

g——重力加速度；

He——井水有效高度；

Hw——上層含水層的水體高度；

d——含水層厚度；

rw——井孔半徑;

aw——與井-含水層系統(tǒng)和地震波有關(guān)的參數(shù)；

S——儲(chǔ)水系數(shù)；

T——導(dǎo)水系數(shù)；

ω——地震波周期。

昌平井的衰減因子為38.99，對(duì)水位數(shù)據(jù)校正，排除井儲(chǔ)影響。校正井孔儲(chǔ)積效應(yīng)后可得到敞口井真實(shí)的水位變化與地震波速度關(guān)系的耦合。

3 孔隙壓力與體應(yīng)力、偏應(yīng)力的耦合

井中水位響應(yīng)于地震波的振幅取決于地震波特征，含水層特性和井孔的幾何特性。本文從井-含水層系統(tǒng)對(duì)地震波的多孔性響應(yīng)特征角度重點(diǎn)研究。Skempton認(rèn)為偏應(yīng)力也能影響孔隙壓力[18]， Holtz和Kovacs闡釋了三維空間內(nèi)偏應(yīng)力和孔隙壓力的相關(guān)關(guān)系[19]：

(4)

式中：ΔPf——孔隙壓力變化；

A和B——Skempton系數(shù)；

σm——平均應(yīng)力；

σd——單軸偏應(yīng)力且σ1=σ2[20]，對(duì)于地震波此應(yīng)力應(yīng)變關(guān)系應(yīng)簡(jiǎn)化為二維空間關(guān)系：

(5)

式中：Ku——不排水條件下的體積模數(shù)；

εν——體應(yīng)變；

N——剪應(yīng)力耦合系數(shù)；

μ——剪力模數(shù)；

εd——偏應(yīng)變。

3.1 潮汐分析估算BKu

封閉性良好的井-含水層系統(tǒng)中，外力導(dǎo)致含水層孔隙介質(zhì)發(fā)生形變，井水位隨之改變。在多孔介質(zhì)中孔隙壓力和體應(yīng)變的關(guān)系可表達(dá)為：

p=BKuε

(6)

式中：p——孔隙壓力；

BKu——潮汐系數(shù)；

ε——體應(yīng)變。

井-含水層系統(tǒng)在潮汐應(yīng)變周期性作用下，含水層中的孔隙壓力會(huì)相應(yīng)產(chǎn)生周期性的變化，由于含水層與井孔之間發(fā)生水量交換，使得井中水位產(chǎn)生與潮汐變化類似的狀態(tài)：

p=ρgh

(7)

聯(lián)立式(6)、(7)得到：

(8)

ρ——水的密度ρ=1 g/cm3；

BKu——潮汐系數(shù)，是表征巖石彈性大小的參數(shù)量。

B是基于孔彈性理論的重要參數(shù)，在0～1范圍內(nèi)，無量綱，與巖石巖性相關(guān)，由水位變化與張力關(guān)系計(jì)算，Ku為正值。

3.2 地震方法估算BKu

3.2.1面波應(yīng)力組成

Love波和Rayleigh波在地球的表面?zhèn)鞑?，在徑向和垂直方向上均有速度，橫向應(yīng)力εTT為0。

σzz=λeεv+2μeεzz+BKuαBεv=0

(9)

σRz=2μeεRZ=0

(10)

σTz=2μeεTZ=0

(11)

εv=εRR+εTT+εZZ

(12)

其中，λe、μe為有效Lame模數(shù)，αB為Biot系數(shù)。利用εTT=0，聯(lián)立式(9)和(12)得：

(13)

(14)

偏應(yīng)力可由地震波速度表示：

(15)

(16)

εRT=εd

(17)

式中：uR和vR——地震波徑向速度；

uT和vT——橫向速度；

VR——Rayleigh波速度為3.5 km/s；

VL——Love波速度3 km/s。

將式(15)帶入式(14)得：

(18)

聯(lián)立式(5)、(16)、(17)、(18)可獲得孔隙壓力與體應(yīng)力和偏應(yīng)力的相關(guān)關(guān)系計(jì)算式：

(19)

當(dāng)只考慮體應(yīng)力對(duì)孔隙壓力的影響時(shí)：

(20)

當(dāng)只考慮偏應(yīng)力對(duì)孔隙壓力的影響時(shí)：

(21)

3.2.2體波應(yīng)力組成

體應(yīng)力可由地震波速度表示：

(22)

(23)

其中，VP,S為P波速2.5 km/s，S波速1.5 km/s。P、S波孔隙壓力與體應(yīng)力和偏應(yīng)力的相關(guān)關(guān)系為：

(24)

4 結(jié)果討論

4.1 潮汐分析

4.2 地震分析

4.2.1Love和Rayleigh波

選擇Love波和Rayleigh波出現(xiàn)的時(shí)間段13:56～14:02作為研究對(duì)象，對(duì)日本地震引起的昌平地下水觀測(cè)井校正井儲(chǔ)效應(yīng)后的實(shí)際水位變化與計(jì)算值曲線擬合，擬合結(jié)果見圖3。

圖3 體應(yīng)力、偏應(yīng)力、體應(yīng)力和偏應(yīng)力對(duì)孔隙壓力的影響Fig.3 Assuming the volumetric strain(a), assuming only deviatoric stress(b) and assuming both volumetric strain and deviatoric stress(c)

4.2.2P和S波

同時(shí)考慮P和S波體應(yīng)力和偏應(yīng)力對(duì)孔隙壓力的影響結(jié)果見圖4。

圖4 體波體應(yīng)力和偏應(yīng)力對(duì)孔隙壓力的影響Fig.4 Assuming both volumetric strain and the deviatoric stress

4.3 潮汐分析與地震分析結(jié)果對(duì)比

潮汐分析得到BKu=6.42 GPa(地震前)及BKu=7.42 GPa(地震后)。地震分析中，綜合考慮體應(yīng)力和偏應(yīng)力對(duì)孔隙壓力的共同影響，利用Love和Rayleigh波為研究對(duì)象，得到BKu=3.32 GPa。潮汐分析與地震分析結(jié)果對(duì)比，肯定了由井-含水層系統(tǒng)對(duì)地震波的響應(yīng)估計(jì)BKu系數(shù)方法的可行性。昌平井以灰?guī)r為主，灰?guī)r實(shí)驗(yàn)室參數(shù)測(cè)量結(jié)果進(jìn)一步驗(yàn)證了地震方法的可行性[21]。

5 結(jié)論

本研究以線性孔隙彈性理論為基礎(chǔ)，探究了井-含水層系統(tǒng)孔隙壓力變化與體應(yīng)變和偏應(yīng)變之間的數(shù)學(xué)關(guān)系，通過對(duì)2011年3月11日日本Mw9.0地震引起的昌平地下水觀測(cè)井實(shí)測(cè)水位變化與曲線擬合得到的計(jì)算值對(duì)比，得到以下幾點(diǎn)認(rèn)識(shí)：

(1)對(duì)于敞口井，井孔儲(chǔ)積效應(yīng)的存在會(huì)削弱水位和水壓，導(dǎo)致測(cè)量結(jié)果偏小。在研究前校正井孔儲(chǔ)積效應(yīng)時(shí)，本文提供了一種較為科學(xué)的方法進(jìn)行校正。

(2)當(dāng)同時(shí)考慮體應(yīng)力和偏應(yīng)力對(duì)孔隙壓力的影響時(shí)，相關(guān)系數(shù)最高，表明孔隙壓力受到體應(yīng)力和偏應(yīng)力的共同影響。