潛水位微動態固體潮影響及主要影響分波識別

史昊鑫, 許 模, 葛建宏, 郭 健

(成都理工大學 地質災害防治與地質環境保護國家重點實驗室, 成都 610059)

0 引 言

固體潮與地下水之間的聯系是水文地質領域的熱門研究課題, 20世紀以來, 不少學者對固體潮與地下水之間的關系進行了探索。Bredehoeft[1]研究了各種要素對井水位變化的影響, 并建立了地下水位與固體潮的函數關系, 成為研究地下水對氣壓效應及固體潮的理論基石; Melchior[2]、Rhoads等[3]給出了固體潮沒有滯后的解釋; Rojstaczer等[4]給出了地下水位對固體潮和氣壓效應的響應方程; 王勇[5]分析指出, 地層年代越老, 巖性越堅硬, 對固體潮響應能力越強; 王鐵城等[6]驗證了在研究水位與固體潮的關系時, 預先消除氣壓影響的必要性; Merritt[7]通過分析固體潮、海洋潮汐和氣壓效應, 估算了美國佛羅里達州含水層系統的水力特性; 呼晶磊等[8]通過對馬17井的觀測, 確定了井水位存在固體潮效應，且水位與固體潮效應之間存在正相關關系; 史浙明等[9]對不同深度承壓井的潮汐進行對比分析, 探索了微動態與地震的關系; 劉陽等[10]通過對潮汐效應的分析, 解釋了潮汐效應機理, 并得出了潮汐效應沒有滯后性; 蘇喬等[11]采用時間序列分析法和互相關分析法分析了潮汐對濱海地下水位的影響并證實了地下水位與潮位頻率相同; 賀前錢等[12]利用井水位和降雨數據模擬計算土壤含水率變化, 估計了地下水變化導致的重力效應; Kazama等[13]通過建立的水文模型, 導出了地下水引起的重力擾動, 并可用來糾正重力數據中的水文干擾; Wang等[14]通過分析潛水含水層中地下水的潮汐效應實現了對俄克拉荷馬州一地區的含水層滲漏監測。

以往對固體潮的研究主要集中在承壓含水層, 針對固體潮對潛水含水層影響的研究較少。本研究以四川省德陽市中江縣垮梁子鄉UG01監測井的地下水位高精度連續監測數據為基礎, 采用頻譜分析法和互相關分析方法, 探討固體潮對潛水位的影響, 并利用調和分析法進一步識別影響水位動態的主要固體潮分波。

1 研究區及井概況

研究區位于四川省德陽市中江縣馮店鎮垮梁子鄉, 丘陵地貌, 地勢由東北向西南逐漸降低(圖1)。區內除第四系全新統堆積層覆蓋外, 基巖主要為侏羅系蓬萊鎮組砂泥巖互層(圖2), 巖層產狀趨于水平, 裂隙較發育, 近水平的層面和兩組垂向裂隙構成了地下水運移的主要通道。

圖2 研究區地層Fig.2 Stratum in the study area

圖1 馮店鎮地貌圖Fig.1 Geomorphic map of Fengdian Town

研究區主要含水巖組為厚度為1～13 m的第四系全新統(Q4)粉質黏土夾少量塊石、碎石的耕植土覆蓋層， 以及下伏侏羅系蓬萊鎮組上段(J3p2)的厚層—巨厚層的砂巖與泥巖、粉砂巖互層, 因區內地下水位多在20 m深度以下, 故地下水類型主要為裂隙水。

監測井UG01井位于山丘斜坡中部, 位置坐標為104°53′53″N、30°38′56″E, 井口高程530 m, 初見水位506.78 m。鉆井完成后, 在34.58 m深度砂巖地層中埋入滲壓計, 并用粗砂和膨潤土進行分層回填。

根據研究區2007—2017年的降雨量數據, 降雨主要集中在5—9月(圖3), 占全年總雨量的82.47%; 10月至次年4月的降雨量僅占全年的17.53%; 10 mm以上的降雨日數年均22.5 d。為減少降雨對UG01井潛水位的影響, 此次數據采集選在秋冬季節(2018-08-20—2019-01-28), 該時段平均降雨頻率為0.1 cpd(circle per day, 日循環次數)。

圖3 2007—2017年月平均降雨量Fig.3 Average monthly rainfall from 2007 to 2017

2 數據來源和研究方法

2.1 數據來源

本次地下水位數據采自RST振弦式滲壓計(加拿大產), 振弦式壓力傳感器輸出為頻率信號, 該頻率與施加在振弦換能器膜片上的壓力成比例且不受導線阻抗和接觸電阻的影響, 這類的振弦式傳感器可以安裝在鉆孔中或打入軟土地基中。數據采集時間間隔設置為1 h, 精度為0.1 cm, 選取2018-08-20T00:00—2019-01-28T00:02時間段的監測數據進行研究, 地下水位變化曲線如圖4(圖中灰色時間段為選取的特征數據段, 具體見3.2節)。

圖4 UG01井潛水位變化曲線Fig.4 Curve of groundwater level in Well UG01

研究中的固體潮監測數據來自中國科學院精密測量科學與技術創新研究院位于麗江的OSG-066超導重力儀, 采用的數據時間從2018-08-20—2019-01-28, 時間間隔1 h, 得到的固體潮變化曲線如圖5所示(圖中灰色時間段為選取的特征數據, 具體見3.2節)。

圖5 當地固體潮變化曲線Fig.5 Curve of local earth tide

對研究區的大氣壓強進行實時監測, 選取的監測設備為Baro-Diver氣壓計, 氣壓數據每隔1 h采集一次, 精度為0.01 cmH2O(0.009 806 hPa)。選用的2018-08-20T00:00—2019-01-28T15: 00的當地氣壓變化曲線如圖6所示。

圖6 當地氣壓變化曲線Fig.6 Curve of local barometric pressure

2.2 研究方法

1)高通濾波法。 高通濾波器是通過頻率篩選并過濾信號的裝置, 其原理是當頻率低于f0的信號輸入這一濾波器時,由于容抗很大而受阻,輸出減小,頻率越低輸出越小。 當頻率高于f0的信號輸入這一濾波器時,由于容抗很小,對信號無衰減作用。頻率f0由下式決定

f0=1/(2πR1C1),

(1)

式中:R1為電阻;C1為電容。

2)頻譜分析法。由于離散的數據信號在時域圖上特征不明顯, 可通過快速傅里葉變換(FFT)將時域信號轉換成頻域信號, 凸顯數據特征。因此, 本文在研究潛水位與固體潮相關周期性時借助快速傅里葉變換, 將離散的時域信號變換為頻域信號分析水位和固體潮的周期性。

傅里葉正變換公式

(2)

傅里葉逆變換

(3)

式中:ω—角頻率(Hz);t—時間(s)。式(2)為時間域的函數轉換為頻率域的函數的積分公式。

對潛水位和固體潮進行FFT分析可以簡單理解為: 把水位或固體潮變化看作是多個簡諧波疊加后的結果, 而FFT則是將水位或固體潮變化的一個整體疊加波拆分成多個簡諧波, 進而從頻域的角度去分析, 由此便可以看到時域上不明顯的周期性特征。

3)互相關分析法?；ハ嚓P函數(cross-correlation function, CCF)是一種表示兩個變量之間相關程度和方向的函數, 給定兩個時間序列Xt和Yt(t=0,±1,…,±n),如果不是平穩的時間序列,可以經過d階差分將其轉換為平穩的時間序列,將差分后的序列記為(x1,y1),(x2,y2),…,(xn,yn)。

樣本的互協方差函數

(4)

(5)

4)調和分析法。固體潮和水位可以看作是由一系列分波疊加的結果,在實際分析中一般選取有限個較為主要的分波。調和分析法可以將各個分波的振幅和遲角從數據中分解出來。

潮位表達式

(6)

式中:S0為余水位;fj為交點因子;μj為交點訂正角;σj是第j個分波頻率;hj、gj是分波的振幅和遲角；v0j為第j個分潮的初相位;t為時間。

2.3 數據預處理

受固體潮影響的地下水位波動為典型的微動態變化[15], 監測的水位數據中除了擁有需要的周期性變化的信號, 也具有非周期性變化的或不需要研究的周期性信號(如降雨)。一次降雨引起的地下水位變動遠大于固體潮的影響, 要找到固體潮對水位的影響, 就需要對數據進行預處理, 消除潛在因素的影響。監測時間內的降雨頻率約為0.1 cpd, 通過高通濾波器, 將頻率為0～0.1 cpd(周期大于10 d)的低頻信號消除, 去除了水位中的非周期性變化后, 進一步對水位數據中的周期性變化(氣壓效應)進行消除[16-17], 得到處理后的水位(圖7)。

圖7 預處理后的潛水位變化曲線Fig.7 Curve of groundwater level after pretreatments

3 固體潮對潛水位影響的基本特征

3.1 周期性

地球固體潮是在太陽、月亮等天體的引潮力作用下, 固體地球產生的周期性的彈性變形現象。地球周期性的膨脹和壓縮變形, 導致巖石空隙的擴大和縮小, 從而使井水位出現周期性的下降和上升現象, 其變化形態與理論固體潮形態相似, 稱為地下水位的固體潮效應[18]。用精密儀器可以觀測到地球固體表層也有和海洋潮汐相似的周期性升降現象。固體潮具有1/2日、日、半月、年、11年、18.6年的周期性變化[19]。

研究區固體潮的周期性分析采用頻譜分析法對其固體潮和地下水位數據進行快速傅里葉變換, 并繪制于半對數坐標中, 得到固體潮頻譜圖和潛水位頻譜圖。

從固體潮頻譜分布圖(圖8a)可見, 在1、2和3 cpd處有明顯的峰值, 且1、2 cpd振幅高于其余的潮峰, 表明固體潮在12、24 h的周期比8 h的周期更顯著。因此, 可以判定固體潮存在明顯的1日、1/2日、1/3日周期性。

由潛水位頻譜分布圖(圖8b)可見, 潛水位在5個不同的周期內存在明顯的波峰, 其中1、2和3 cpd的周期性波動特征與固體引力潮重疊, 從變化周期的角度來看, 潛水位存在與固體潮相似的1日、1/2日和1/3日周期性關聯。

圖8 固體潮(a)及潛水位(b)頻譜分布圖Fig.8 Spectrogram of earth tide(a) and groundwater level(b)

3.2 負相關性

為充分研究反映潛水位與固體潮的相關性, 從時間序列中抽取2018-11-06—13和2019-01-02—09兩周的固體潮和地下水位數據進行曲線對比。水位與固體潮變化曲線如圖9、圖10。從圖形的形態特征來看, 固體潮的峰(谷)值與潛水位的谷(峰)值對應關系明確, 固體潮與潛水位之間表現出明顯的負相關性。

圖9 2018-11-06—13固體潮和潛水位周變化曲線Fig.9 Curves of earth tide and water level from 6 to13 in November of 2018

圖10 2019-01-02—09固體潮和潛水位周變化曲線Fig.10 Curves of earth tide and water level from 2 to 9 in January of 2019

利用式(4)、(5)對固體潮和地下水位作互相關函數計算, 得到分析結果如圖11所示。地下水位與固體潮在Lag=0時的相關性系數最高, 為-0.864, 反映固體潮與潛水位呈負相關且不存在明顯的時間滯后。即當固體潮增大的時候, 地下水位隨即下降, 水位埋深增大; 當固體潮減小的時候, 地下水位隨即上升, 水位埋深減小。呼晶磊等[8]在研究了馬17井承壓含水層地下水位后得出水位與固體潮呈負相關, 并解釋其原因是當潮汐應力增大時, 含水巖層發生體膨脹, 含水層的孔隙水壓降低, 井水位下降; 當潮汐應力變小時, 含水巖層發生體壓縮, 含水層的孔隙水壓升高, 井水位上升。本研究也在潛水含水層中得出了相同結果, 但是由于本研究的井結構特殊, 采用的傳感器也非常規的壓電式壓力傳感器, 故在非承壓含水層中這種負相關關系的物理原因仍需進一步探索。

圖11 固體潮與潛水位互相關系數分布圖Fig.11 Cross-relation coefficients between earth tide and groundwater level

4 主要影響分波的識別

從前文可知固體潮對潛水位存在明顯影響, 為進一步分析影響地下水位的固體潮分波, 采用調和分析法得到地下水位和固體潮的主要分波的振幅, 進而結合各波的重疊頻率域來識別影響水位波動的固體潮分波[20]。調和分析法采用MTALAB實現, 從結果中提取出日波、半日波和1/3日波的振幅, 結合各波對應的頻率, 得到各分波的頻率振幅數據(圖12)。

圖12 固體潮(a)及潛水位(b)主要分波頻譜圖Fig.12 Main frequency components of earth tide(a) and groundwater level(b)

固體潮和潛水位的頻譜分析顯示, 固體潮頻譜分析圖中的主要成分為O1(頻率 0.928 8 cpd, 振幅 31.74 μGal)、K1(頻率 1.003 2 cpd, 振幅 43.70 μGal)日波, N2(頻率1.896 0 cpd, 振幅12.65 μGal)、M2(頻率1.932 0 cpd, 振幅64.86 μGal)、S2(頻率1.999 2 cpd, 振幅30.65 μGal)半日波, M3(頻率2.899 2 cpd, 振幅1.02 μGal)1/3日波; 而潛水位頻譜分析圖中的主要成分為K1(頻率1.003 2 cpd, 振幅9.1 mm)日波, S2(頻率1.999 2 cpd, 振幅9 mm)半日波, SK3(頻率3.002 4 cpd, 振幅1.6 mm)1/3日波, 因此, 從兩者間重疊的峰值部分可知, 固體潮對水位的主要影響分波為K1波和S2波。

5 結 論

垮梁子鄉的UG01監測井為一封閉井, 其潛水位實時監測數據不受氣壓效應影響, 固體潮作用是潛水位微動態變化的主要原因, 主要結論如下:

(1)通過頻譜分析法, 得到固體潮存在明顯的8、12和24 h的周期變化, 研究區潛水位也表現出了相同的周期性波動。

(2)時間序列和互相關分析均表明, 固體潮和潛水位標高呈負相關關系, 即固體潮增大時, 水位下降, 固體潮減小時, 水位上升; 且固體潮與潛水位之間不存在明顯的滯后性, 在非承壓含水層中這種負相關關系的物理原因仍需進一步探索。

(3)利用調和分析法, 得出影響地下水位的主要固體潮分波為K1日波、S2半日波。

致謝: 中國科學院精密測量科學與技術創新研究院的陳曉東研究員提供了固體潮數據, 在此表示衷心的感謝!