华北平原地下水位微动态变化周期特征分析

王仕琴，宋献方，王勤学，唐常源，刘昌明

(1.中国科学院遗传与发育生物学研究所农业资源研究中心，河北 石家庄 050021;2.日本千叶大学园艺学部，松户 271-8510;3.中国科学院陆地水循环与地表过程重点实验室/中国科学院地理科学与资源研究所，北京 100101;4.日本国立环境研究所，筑波 305-8506)

地理位置、地下水埋深、含水层的封闭条件、岩性等方面的不同，地下水位体现出微观的动态特征，可以反映太阳和月亮潮汐变化、地球固体潮、蒸散发、大气压力以及其它外部荷载等影响因素对地下水的作用。其中一些影响因素如天体潮汐等具有明显周期性波动特征，导致地下水位过程线也具有周期波动的特征，可以认为是由许多不同频率、不同振幅的谐成分合成的［1～4］。由于频谱分析对周期性影响因素的显著识别，在地下水位的微动态分析中是一种有效的分析方法［5～9］，可以把一个复杂的地下水位变化过程分解为具有固定周期的一系列谐波，利用该谐波的特征确定地下水的微观影响因素。

在与海洋相接触的海岸带潜水含水层中，因潮汐作用，地下水也出现周期性正弦波动。在距海洋很远的内陆地区，受太阳和月球的引力引起深层含水层的地壳应力变化，从而地下水位随天体潮汐作用产生周期性波动也会在深层含水层中体现出来，这主要取决于含水层的封闭条件，岩石弹塑性或岩石压力传导性质或水的物理性质［3］。因此通过地下水位微动态波动对天体潮汐作用的反应程度也可以用来判断含水层性质。Godin［10］指出天体潮汐作用所反映出来的半日M2，和一日 M1，S1和 K1的频率分别为 0.08051，0.04027，0.04167和0.04188 cph，对应的周期分别为12.42082，24.83238，23.99808 和 23.87775 h。Shih［4，7］利用频谱分析对台北盆地的地下水微动态进行分析后得出，地下水的波动存在日，半日以及1/4日的显著周期，此外，还存在 8.3077，6.1714，4.8869 和 4.1065h的周期，天体的潮汐作用是引起海水和Choshuihsi冲积平原地下水位周期波动的主要原因。

研究地下水的微动态需要有较高频率的地下水位观测资料，本研究基于中国科学院地理科学与资源研究所、日本国立环境研究所和中国地调局天津地质矿产研究所合作开展黄河流域地下水资源监测和地下水资源评价数据库建设研究工作，利用日本先进的自动监测设备KADEC-MIZUⅡ型地下水位计对黄河流域包括华北平原地下水动态进行监测，是对现有监测系统的补充和提高，对掌握地下水微动态规律具有重要作用。本文在研究华北平原地下水宏观水位动态的基础上［11］，通过对地下水的微动态周期的分析，进一步认识地下水的动态变化规律。

1 研究区基本情况与地下水位监测

1.1 华北平原基本概况

华北平原在地貌上处于太行山山脉以东，燕山山脉以南，海拔不超过100m，自北、西、南西三个方向向渤海湾倾斜。从山麓至渤海海岸，分为山前冲积洪积倾斜平原、中部冲积湖积平原、东部冲积海积滨海平原三部分。多种复杂的大型地貌和小型地貌交错重叠，有大型冲积扇、扇间洼地、河道带、河间带，也有河口三角洲、岗地、浅碟状洼地和条状背河洼地等［11］(图 1)。

图1 华北平原地形地貌及地下水监测孔分布图［10］Fig.1 Location of groundwater monitoring borehores in the North China Plain

1.2 地下水位监测

2003年末～2004年7月依次在华北地区选取39眼农用或生产生活观测孔(图1)，观测孔位置和设置详细资料与王仕琴等［11］分析华北平原地下水位宏观动态的观测孔一致。选择观测孔时主要考虑了潜水含水层有关的各种因素包括监测孔所在华北平原位置的代表性，仪器管理和数据收集的便利等方面。其中在中科院崇岭实验流域山区布设两眼(29和30号观测孔)，其它观测点布设于华北平原区，其中2、6、23、25号井监测承压水，这些点主要位于北部，仅代表局部的地下水动态。

水位监测仪器采用日本国立环境研究所提供的KADEC-MIZUⅡ型地下水位计，地下水位监测频率为30min或1h(北京为1d)，连续监测地下水位的变化。

2 频谱分析方法

本文采用谐波分析方法进行周期分量的分析提取。谐波分析是利用傅立叶级数把某个资料的时间序列表示成无数个不同周期的简谐波和的形式分析序列变化规律的一种方法［12～13］。对序列yt可用L个波叠加的形式表示其估计值:

L——谐波个数，取n/2的整数部分;

k——谐波序号，k=1，2，…，L;

ak，bk——傅立叶系数。

故根据系列yt，可求得a0，ak，bk。为节省工作量，通常在L个波中选取波动比较显著的几个谐波相加估计yt，在实际应用中只需选取前6个显著谐波就能满足精度要求。

式中:s2——系列yt的方差;

α——检验的显著性水平(一般取5%)。

将显著波所对应的a0，ak，bk代入式(1)中即得周期分量的估计方程［14］。

一般取前6个显著周期为主要周期。

3 周期计算

3.1 频谱分析对周期识别的显著性验证

本文首先选取河南的33号观测孔2008年8月15日～2009年1月31日地下水位变化稳定时期的4078个地下水位监测数据进行分析，研究其频谱分析方法对周期识别的显著性和有效性。虽然识别周期显著性的数据只有半个水文年，由于频谱分析主要用于分析时间序列数据中的微观周期项成分，而且这些周期项成分长期存在，因此，对于分析微观动态来说具有一定的代表性。

频谱分析法适用于长序列平稳动态资料，即动态历时曲线的周期成分比较显著，而对于非平稳的动态历时曲线，应先对水位序列进行平稳化改造，即消除趋势项。通过初步图形直观分析，2008年8月15日～2009年1月31日地下水动态具有明显上升趋势，且具有以天为周期的动态变化(图2)。

首先，提取趋势项成分，对于趋势分量ft可用多项式逼近，即:

采用多元回归方法确定待定系数c0，c1，c2，…，ck和阶数k。通过分析以拟合度R2为0.923的线性方程拟合最好(表1)。

表1 线性回归方程拟合参数表Table 1 Linear regression equation and its parameters

对于原序列减去趋势项后得到剩余项(图3中实测值数据)，进行周期项频谱分析，利用F分布检验，在显著性水平α为0.05状态下共得到47个显著性周期(由于显著性分波参数较多，这里不列出参数)。这也说明在小频率观测情况下，可以获得地下水细微动态变化特征的数学随机表达。对趋势剩余项实测值与频谱分析计算值进行拟合，见图3，拟合程度较好，除了不能体现水位波动的随机变化外，频谱分析可以很好地体现周期变化。实测值和模拟值的平均残差为-0.0052，标准残差为0.112。由此也说明频谱分析对于地下水位微动态变化规律周期的识别也是显著的。

频谱分析结果包括12、24的显周期，这些显周期可以从实测水位中观察到(图2)，也是微动态变化的主要周期，这种周期性的变化主要是由气温的变化、抽水或天体潮汐变化引起的水位波动导致的，虽然该观测孔是潜水井，但由于位于冲洪积平原，含水层具有互层性，地下水位埋藏较深，因此也体现了由岩石弹性形变而产生的水位微动态变化［1］，12h和24h的周期波动可能由天体的潮汐变化引起。另外，频谱分析表明还包含很多隐性的观察不到的长周期，如大约5～9d、一个月的周期变化，这也可能说明固体潮的影响，进一步说明地下水对外界应力响应的灵敏性。

由于频谱分析对周期分析的显著性，可应用于整个华北平原，分析长期固定的微动态周期规律，说明华北平原不同地区不同观测点其微观动态周期变化特征。

3.2 华北平原地下水位微观周期特征

图2 地下水短周期变化图Fig.2 Short periods of groundwater levels

图3 地下水实测值与随机模型计算值拟合图Fig.3 Fitting of calculated values and observed values of groundwater levels l

由于本文所要揭示的是时间序列波动所体现的微观动态变化，如果微动态的周期波动受天体的潮汐作用影响其周期常年存在，而浅层地下水由于蒸散发造成的微动态变化则主要发生在夏季，因此在计算中选择夏季7月1日～7月10日连续10天的实际观测序列来分析。从实测的水位波动以及相对应的去除趋势项的平稳序列来看，有些观测点体现出明显的12h或24 h 的周期变化(图4 中 10、11、12、13、16、18、19、20、21、22、23、25、27、29、31、33、34、37、39)。通过频谱分析计算结果可以得到其真正的显著周期变化。分析过程中认为得到的太小的周期(小于4h)以及时间太长的周期(大于48 h)不列入本文研究的周期。

图4 华北平原微动态各观测孔原始观测序列和平稳序列曲线图Fig.4 Observed series and smoothed series of micro-dynamic curves of groundwater levels in the North China Plain

通过周期分析结果把华北平原(北京除外)观测 点分为以下几类(表2):

表2 华北平原地下水微动态频谱分析结果Table 2 Results of spectrum analysis of themicro-dynamics of groundwater in the North China Plain

(1)沿海地区:沿海地区的地下水观测孔埋深较浅，主要为潜水含水层，一类与海潮潮汐作用有关，另一类与蒸散发作用有关。

(a)受海潮潮汐影响显著。如11、12、26、27和38号观测点埋深较浅(小于10m)，具有与海潮潮汐变化一致的周期(24、12、8、6 或4 h)，其中11 号观测孔，虽然距离海岸较远，但是其位于河道带，也体现出一定的潮汐作用。

(b)受蒸散发作用显著。如13和37号观测点埋深很浅(小于3 m)，主要存在24h的周期变化，受蒸发的影响较大，而13号距离海岸较近也会受到潮汐作用的影响，体现出24 h和12 h的周期。

(2)内陆地区:一般情况下，在内陆的深层承压含水层由于岩石的应力状态受天体潮汐的作用发生变化从而导致地下水位也会呈现出周期性的波动，实际情况下，浅层含水层也有可能出现由于应力变化而产生的周期波动，而深层地下水位也可能不会有周期性波动，这反映了岩石或含水层的特征。

(a)地下水埋藏较浅区(小于5 m)但具有与潮汐变化一致的周期。如22、29和31号观测点，存在与潮汐变化一致的周期，推断含水层为刚性体所造成。其中29号观测孔，从日水位波动来看，可能受到抽水的作用影响较大，体现一些类似天体潮汐波动的周期变化，需要进一步区分判别。

(b)含水层埋藏浅(小于10 m)对外界反映不灵敏，如24、28、30、31、32、36 和 39 号孔地下水埋深，没有相对于潮汐周期变化的波动，每日水位波动没有明显特征。

(c)地下水埋藏较深区(大于10 m)潮汐作用显著的观测点，如 10、16、17、18、19、20、25、33、34 号观测点，体现为具有 24、12、8、6、4 h 的周期，说明深层含水层岩体的弹塑性较大，受太阳和月球引力作用造成的含水层形变而导致地下水位的波动显著。其中17和34号观测孔，从日水位波动来看，可能受到抽水的作用影响较大，体现一些类似天体潮汐波动的周期变化，需要进一步区分判别。

(d)地下水埋藏较深(大于10 m)，但受潮汐作用周期波动不显著的点，如14和21号观测点，说明含水层不受外界影响，其含水层岩性对应力-应变的反应不灵敏，或者含水层封闭性差等原因造成。

4 结论

(1)近海地区:潜水含水层地下水微动态体现出24、12、8、6或4 h的周期性变化，说明地下水位受到海潮潮汐作用明显，而地下水埋深非常浅的观测孔(小于5 m)其受到蒸发作用的影响体现出24 h的微动态。

(2)内陆地区:地下水埋藏较深的观测孔(大于10 m)，大多数观测孔具有与天体潮汐一致的周期变化，如12，24 h的周期性。只有个别观测孔未体现出类似的周期波动，说明观测孔所在含水层的岩性弹性较差，对天体等外界应力的反映不灵敏。内陆地区地下水埋藏较浅的观测孔，大多数不存在潮汐周期作用下的波动，符合一般规律，只有个别观测孔在地下水浅埋条件下能体现出与天体潮汐作用一致的周期，说明含水层岩性的弹性条件或含水层封闭性较好。

对于地下水位出现与天体潮汐影响一致的周期性特征的分析，在内陆地区需要结合固体潮汐的观测，而在沿海地区需要结合海潮潮汐的观测进一步分析和明确。此外，在收集气压气温等条件的基础上，进一步综合研究不同地区的微动态变化规律。

［1］ 周鹏鹏，李国敏，卢耀东，等.大陆岛地下水动力学特征——以湛江东海岛为例［J］.水文地质工程地质，2013，40(1):12-18.［ZHOU P P，LI G M，LU Y D，et al.The hydrodynamic characteristics of groundwater in continental island—an example from Donghai Island，Zhanjiang［J］.Hydrogeology ＆ Engineering Geology，2013，40(1):12-18.(in Chinese)］

［2］ 高茂生，叶思源，史贵军，等.潮汐作用下的滨海湿地浅层地下水动态变化［J］.水文地质工程地质，2010，37(4):24-27.［GAO M S，YE S Y，SHI G J，et al.Oceanic tide-induced shallow groundwater regime fluctuations in coastal wetland ［J］.Hydrogeology ＆ Engineering Geology，2010，37(4):24-27.(in Chinese)］

［3］ 陈葆仁，洪再吉.地下水动态及其预测［M］.北京:科学出版社，1988:1-382.［CHEN B R，HONG Z J. GroundwaterDynamicsand Prediction［M］.Beijing:Science Press，1998:1-382.(in Chinese)］

［4］ 于福荣，卢文喜，李平，等.地下水管理模型求解方法综述［J］.水文地质工程地质，2010，37(2):27-31.［YU F R，LU W X，LI P，et al.Review on optimization methodsforgroundwatermanagement models［J］.Hydrogeology ＆ Engineering Geology，2010，37(2):27-31.(in Chinese)］

［5］ Shih DC F，Chiou K F，Lee C D，et al.Spectral responses of water level in groundwater and tidal river［J］.Hydrol Processes，1999，13(6):889-911.

［6］ Shih DC-F.Determination of hydraulic diffusivity of aquifers by spectral analysis［J］.Stochastic Environ Res Risk Assessment，1999a，13(1/2):85-99.

［7］ Shih DC-F.Inverse solution of hydraulic diffusivity determined by water level fluctuation［J］.J Amer Water Resour Assoc，1999b，35(1):37-48.

［8］ Shih DC F，Lee C D，Chiou K F，et al.Spectral analysis of tidal fluctuations in groundwater level［J］.J Amer Water Resour Assoc，2003，6(5):1087-1100.

［9］ Stephen McLin.Spectral analysis of observed aquifer water level fluctuations［J］. Computers and Geosciences，2012，46:264-271.

［10］ Godin G.The Analysis of Tides［M］.University of Toronto Press，1972:232-233.

［11］ Wang Shiqin，Song Xianfang，Wang Qinxue，et al.Shallow groundwater dynamics in North China Plain，China［J］.Acta Geographica Sinica，2009，19(4):175-188.

［12］ 杨忠平，卢文喜，李平.时间序列模型在吉林西部地下水动态变化预测中的应用［J］.水利学报，2005，36(12):1475-1479.［YANG Z P，LU W X，LI P.Application of time-series model to predict groundwater regime［J］.Journal of Hydraulic Engineering，2005，36(12):1475-1479.(in Chinese)］

［13］ 李平，卢文喜，杨忠平.频谱分析法在吉林西部地下水动态预报中的应用［J］.水文地质工程地质，2005，32(4):70-73.［LI P，LU W X，YANG Z P.Application ofspectrum analysismethod to the prediction of groundwater regime in west Jilin province［J］.Hydrogeology＆ Engineering Geology，2005，32(4):70-73.(in Chinese)］

［14］ 杨志霞.用时间序列模型预测深层地下水水位［J］.水资源研究，2010(4):1-4.［YANG Z X.With a time seriesmodelto predictthe deep groundwater level［J］.Water Resources Research，2010(4):1-4.(in Chinese)］