0 引言

地下水是反映地壳应力与固体变形最敏感的物质之一，当井-含水层系统处于封闭性良好的承压体系中时，井-含水层系统即为一个天然体应变仪(Hsieh et al，1987)，其对地壳体应变响应灵敏度可达10^-10量级(晏锐，2008)。外力作用引起含水层介质体应变的微小变化，井水位都能灵敏地作出响应(刘序俨等，2009)。研究表明，地震引起的各种水文响应中，井水位的变化最为普遍，地震孕育过程中的应力积累、地震发生后断层位错引起的静态应力以及地震波传播引起的动态应力，都会引起局部或区域尺度地壳介质结构的改变，这就不可避免地导致岩体中孔隙压力和含水层介质特性(如渗透率)的改变，从而引起井孔水位的变化(Cooper et al，1965；King et al，1999；Manga et al，2007；Shi et al，2013)。因此，研究地震引起的井水位同震变化能够了解地下介质对应力-应变的响应过程、孔隙压变化特征及介质参数变化的信息等，对于理解流体在地震孕育、发生过程中的作用，认识与地震相关的水文响应机理具有重要意义(Wang et al，2014；Montgomery et al，2003)。

水位的同震变化是地震波作用于井-含水层系统最直接的体现形式之一。国内外学者在其变化机理研究方面取得了很多进展，例如，Cooper等(1965)研究了水位对地震波的影响因素，认为影响水位对地震波响应的因素有井孔的尺寸、含水层的导水系数、贮水系数、孔隙度以及波的类型等，并得到了井水位对地震波作用的理论解析解；Manga等(2007)研究指出，远场井水位的响应幅度取决于水位埋深、井半径、含水层厚度和水文地质参数等，通常含水层导水系数越大，引起的水位变化幅度也越大，地震引起水位持续变化的原因是含水层系统渗透性的改变(Manga et al，2016；Shi et al，2016；Sun et al，2015；Yan et al，2016)；国内学者车用太等(1989)、鱼金子等(1993)研究认为，井孔水震波记震能力的强弱与渗透性的强弱有关；汪成民等(1983)提出，井孔反映地震波的能力受井孔自身特性和含水层岩体中孔隙裂隙水的流动条件以及水进出井孔的速度有关；张昭栋等(1999)应用井水位的振荡曲线计算了含水层的导水系数，发现井含水层的导水系数越大，井水位对地震波的响应幅度越大，并通过进一步的SLUG实验(张昭栋等，2000)证明了当固有振动周期接近瑞利波的振动周期(20s)时，振动阻尼系数越小、渗透系数越大时，井水位对地震波的响应越强烈。综上可知，地震波不仅能引起井水位出现水震波，还可以改变井-含水层系统的水文地质参数，如导水系数、储水系数、渗透系数等(Brodsky et al，2003；Elkhoury et al，2006；Manga et al，2012)。因此，分析地震引起的井水位水震波响应特征，是了解井区地质构造和水文地质条件，认识地震波作用下含水层介质的应力应变，探究井水位同震变化机理的有效手段之一。

除了机理研究的进步之外，水位数据采集技术也有所进步，从最早的模拟记录，发展到了当前的数字化分钟值及秒值记录。研究表明，水震波的周期一般在十几秒至二十几秒之间，而地震波也是频率较高的信号，相比传统的模拟记录或数字化分钟值采样，秒采样率的高频水位记录对P波、S波和面波等都有直接反映，能更清晰完整地记录到地震引起的水震波(He et al，2016；舒优良等，2006)，为研究水震波机理及影响因素、确定水文地质参数，提供了更为丰富的资料(Shih et al，2013；Sun et al，2015；舒优良等，2014)。

研究表明，新10井对地震波的响应较灵敏，记震能力较强(向阳等，2017)。配备数字化高频采样水位仪后，新10井同时记录到了2017年8月8日21时19分四川九寨沟县(33.20°N，103.82°E)M_S7.0地震和8月9日7时27分新疆博尔塔拉州精河县(44.27°N，82.89°E)M_S6.6地震所引起的水震波。但二者的同震响应形态不同，前者表现为井水位的振荡—平稳，后者表现为井水位的振荡—上升。本文对比分析了新10井秒采样记录到的水震波与地震波(地表垂向运动)之间的相关性特征，以及其与井-含水层系统介质参数间的关系，定量地研究了井水位同震响应特征，并讨论了其成因机理。

1 地质背景及井孔条件

新10井位于乌鲁木齐市南部，其地理坐标为43.70°N、87.62°E，海拔高度为1056m。构造上位于柳树沟-红雁池逆冲断裂及其派生断裂的交汇部位(图 1(a))，沿断裂挤压破碎带广泛出露泉水，多数泉水流量常年稳定，涌水量达20~30L/s。新10井于1980年9月成井，其井孔结构如图 1(b)所示。井深28m，开孔孔径146mm至9.2m处，130mm孔径从9.2m至13.04m，13.04m以下孔径为110mm。井孔16~24m岩芯破碎，为主要出水段，在该处设置滤水管。含水层为二叠系红雁池组硅质砂岩和砾岩，厚度约8m，井水位埋深1m左右。2016年10月底对其井水位观测仪器系统进行了改造，安置SWY-Ⅱ型数字水位仪，分辨率为1mm，采样率1s，信息记录能力明显提高。

2 同震响应特征对比分析

图 2(a)为新10井水位对2017年8月8日九寨沟M_S7.0、2017年8月9日精河M_S6.6地震的同震响应分钟值曲线。由图 2(a)可见，该井水位对九寨沟地震的响应形态为振荡型，而对精河地震的响应形态为振荡-阶升型，2次地震井水位同震响应形态存在差异。而与新10井水位分钟值表现形态有所不同的是，同观测井的秒采样水位仪清晰记录到了2次地震引起的井水位水震波波形(图 2(b)、2(c))，所以，九寨沟地震和精河地震都引起了新10井水位的振荡变化，前者表现为井水位的振荡，振荡后水位很快恢复平静；后者亦表现为井水位的振荡，但振荡后水位出现升高。以上说明井水位分钟值采样记录到的水震波波形不完整，存在信息缺失，而高采样率的秒采样则能更完整地记录到水震波信息，提高了同震响应分析结果的准确性。

图 2(b)为新10井水位记录到的九寨沟M_S7.0地震的水震波波形，震中距为1830km，水震波振荡持续时间约为3min；图 2(c)为新10井水位记录到的精河M_S6.6地震的水震波波形，震中距为390km，水震波振荡持续时间小于1min，持续时间较短。此外，由图 2还可见，九寨沟地震引起的水震波振荡周期更长。说明新10井水位经历九寨沟M_S7.0地震时，地震波作用时间较长，水震波波形发育，水位振荡持续时间较长；经历精河M_S6.6地震时，由于震中距近，地震波作用时间较短，水震波波形不发育，水位振荡持续时间很短。前人研究表明，水震波的持续时间和振荡周期与地震波密切相关(Cooper et al，1965；Sun et al，2015)。

为了进一步对比分析该井水震波对地震波的响应特征，本文收集了新10井东北方向相距约10km处的水磨沟台SLJ-100型强震计记录到的垂向速度波形数据。图 3(a)、3(e)分别为九寨沟M_S7.0地震引起的新10井水位同震响应和相应的水磨沟台垂向分量位移波形；图 3(b)、3(f)分别为精河M_S6.6地震引起的新10井水位同震响应和相应的SMG台垂向分量位移波形。利用S变换法(Stockwell et al，1996)对2次地震的地震波和水震波进行了时频特征提取，图 3(c)、3(g)分别为利用S变换得到的九寨沟M_S7.0地震时新10井水震波与SMG台地震波垂向分量信号的时频特征图，图 3(d)、3(h)分别为利用S变换得到的精河M_S6.6地震时新10井水震波与水磨沟台地震波垂向分量信号的时频特征图。

由九寨沟M_S7.0、精河M_S6.6地震引起的新10井水震波形态(图 3(a)、3(b))与其相应的SMG台地震波垂向分量的形态(图 3(e)、3(f))可清晰地看到S波和面波。图 3(c)、3(g)的时频特征显示，九寨沟M_S7.0地震时，水震波和地震波均存在一个显著的频率段，即0.03~0.14Hz(7~33s)，而且在0.06~0.07Hz(14~16s)频段内水位的变幅度最大，幅度值为8mm(红色表示高值、蓝色表示低值)，说明在14~16s频段内新10井水位对地震波引起的地壳变形放大作用最强，而这也是面波(瑞利波)的主要周期成份。图 3(d)、3(h)的时频特征显示，精河M_S6.6地震时，水震波和地震波都存在2个频率段，即0.2~0.5Hz(2~5s)和0.04~0.12Hz(8~25s)，但相比地震波而言，水震波8~25s的长周期信号不明显，水位在0.25~0.35Hz(3~4s)频段内的变化幅度最大，幅度值为15mm，说明新10井水位在3~4s周期内地震波引起的地壳变形放大作用最强，而在8~25s(0.04~0.12Hz)周期内的放大作用很小。

综上分析可知，新10井水位对九寨沟地震和精河地震的水震波与地震波形态较为相似，整体来看，地震波幅度较大时，水震波引起的水位变化的幅度也越大。但也存在局部的差异性，新10井水位的水震波和地震波在九寨沟M_S7.0地震地震波作用时，二者主要表现为长周期成分，并且一致性较好(7~33s)；而在精河M_S6.6地震地震波作用时，二者的短周期成分较一致，水震波的长周期成分几乎没有(8~25s)，这可能是地震发生的地点距新10井较近，地震波作用于井水位的时间很短，水震波不发育所致。可见，水震波与地震波关系密切，由于井水位的同震响应是地表垂直运动和含水层孔隙压力波动共同作用的结果，因此，井水位的变化在一定程度上反映了地震波引起的地壳介质体应变的放大作用与井孔条件、含水层参数、地震波周期等密切相关(Cooper et al，1965；Sun et al，2015)。

3 地震波引起的含水层参数变化

从以上分析可知，地震波作用除了会引起地表介质的上下波动之外，还会引起含水层内孔隙压力的波动变化。水震波和地震波(垂向分量)在形态和频率上均高度相似，但不同周期的信号引起的水震波的响应幅度不同，这与含水层水文地质参数密切相关，特别是与储水系数和渗透系数有关(Cooper et al，1965)。而获取水文地质参数的方法可以分2大类，一类是水文地质实验法(抽水实验、注水实验、渗水实验等)，另一类是地下水动态观测资料反演。水文地质实验法对于井-含水层系统具有较大扰动，因此，不适合地震监测井这类需要连续无干扰观测条件的井孔。为获得地震前、震时及震后含水层渗透性的变化情况，本文将分别从以下2个方面获取含水层参数(储水系数和渗透系数)：一方面是利用水震波和地震波的形态特征反演地震波作用过程中的储水系数和渗透系数；另一方面是利用固体潮分析法，反演地震波作用之前、之后的储水系数和渗透系数，以期为进一步探究井水位同震变化机理提供理论依据。

3.1 地震波作用过程中的含水层参数

水位对地表运动的响应幅度与储水系数、含水层厚度和面波波数等因素密切相关，利用水位与垂向位移间的功率谱密度关系，可以求得井-含水层系统的储水系数(Shih，2009)

$ S = \frac{{1.836{\rm{\pi }}{b}}}{\lambda }{\left( {\frac{{{S_{hh}}}}{{{S_{ww}}}}} \right)^{ - 1/2}} $

(1)

$ S = \frac{{1.836{\rm{\pi }}{b}}}{\lambda }{\left( {\frac{{{S_{hh}}}}{{{S_{hw}}}}} \right)^{ - 1}} $

(2)

$ S = \frac{{1.836{\rm{\pi }}{b}}}{\lambda }{\left( {\frac{{{S_{wh}}}}{{{S_{ww}}}}} \right)^{ - 1}} $

(3)

式中，b为含水层有效厚度；λ为波长，由地震波波速和周期确定；S_hh、S_ww分别为水位h(t)和地面运动垂向位移w(t)的自功率谱密度，二者数值相等；S_wh、S_hw分别为h(t)和w(t)的互功率谱密度。式(1)、(2)、(3)分别为利用水位自功率谱与垂向位移自功率谱的比值、水位自功率谱与水位-垂向位移互功率谱的比值、水位-垂向位移互功率谱与垂向位移自功率谱的比值得到的含水层的储水系数，依据式(1)~(3)可分别得到S_hh、S_ww以及S_wh和S_hw 的值，从而估算出井-含水层系统平均储水系数。

图 4(a)、4(b)分别为九寨沟M_S7.0、精河M_S6.6地震引起的水位、垂向位移的自功率谱密度和互功率谱密度。由图 4可见，九寨沟地震在0.04~0.07Hz频段内、精河地震在0.5~0.6Hz频段内水位，垂向位移的功率谱密度最大，且一致性较高，因此，可用该显著频段的信号来估算储水系数(Shih et al，2013)。

新10井含水层在埋深16~24m处(图 1(b))，含水层有效厚度b为8m，计算波长λ=v/ω时所需的波速v值由震中距除以面波到时与发震时刻之差获得，频率ω取决于图 4中的功率谱密度的最高值。地表运动(地震波垂向分量)的质点位移幅度由速度幅度转化得到，即s_max=v_max·τ/2π(τ为地震波周期)。依据式(1)~(3)，利用九寨沟M_S7.0、精河M_S6.6地震引起的水位、垂向位移的功率谱密度关系，分别估算出新10井含水层平均储水系数为2.10×10^-5、3.08×10^-5(表 1)。由表 1可见，精河地震后新10井含水层的储水系数增大。

在承压性的井-含水层系统中，井水位对孔隙压力和地表垂向运动的响应可用放大因子(或幅度响应比)来表示(Cooper et al，1965)

$ \left. \begin{array}{l} A = {\left[ {{{\left( {1 - \frac{{{\rm{\pi }}r_{\rm{w}}^2}}{{T\tau }}K{e_{\rm{i}}}{\alpha _{\rm{w}}} - \frac{{4{{\rm{\pi }}^2}{H_{\rm{e}}}}}{{{\tau ^2}g}}} \right)}^2} + {{\left( {\frac{{{\rm{\pi }}r_{\rm{w}}^2}}{{T\tau }}K{e_{\rm{r}}}{\alpha _{\rm{w}}}} \right)}^2}} \right]^{ - 1/2}}\\ A' = A \cdot \frac{{4{{\rm{\pi }}^2}{H_{\rm{e}}}}}{{{\tau ^2}g}},{\alpha _{\rm{w}}} = {r_{\rm{w}}}{(\omega {S/T})^{1/2}} \end{array} \right\} $

(4)

式中，A为水位对含水层孔隙压力波动的放大因子；A′为水位对地表垂向位移的放大因子；r_w为井孔内径；H_e为观测井有效水柱高度；T为导水系数；S为储水系数；τ为地震波周期；ω为地震波频率；α_w为中间变量；g为重力加速度；Ke_i和Ke_r为零阶开尔文函数的虚部和实部。依据Cooper等(1965)、向阳等(2017)和式(4)，可反演得到含水层等效厚度和渗透系数(图 5)。依据式(4)得到渗透系数K值，其中，井径r_w=54mm，含水层厚度d=8m，有效水柱高度H_e=18m，综上所述，九寨沟地震时新10井储水系数S为2.10×10^-5，渗透系数K为61m/d；精河地震时新10井储水系数S为3.08×10^-5，渗透系数K为147m/d(表 1)。

基于以上分析，利用水震波与地震波特征反演的九寨沟M_S7.0、精河M_S6.6地震新10井含水层参数表明，精河地震时孔隙度和渗透系数都出现了不同程度的增大，尤其是渗透系数增大了1倍左右。说明九寨沟地震和精河地震地震波作用期间新10井的渗透性明显增强。

3.2 地震波作用前后含水层参数

以上利用水震波和地震波响应特征获得的是地震波作用过程中新10井的含水层参数，而且依据Cooper等(1965)的理论模型，反演的是水平向的含水层参数，即含水层中应力变化与井孔中水位变化之间的关系。而地震波作用前后新10井-含水层系统水文地质参数的变化则可利用潮汐响应(相位滞后与响应幅度)来计算，即利用水位对潮汐的响应幅度反演含水层储水系数，利用相位滞后来反演渗透系数(Hsieh et al，1987；Roeloffs，1996)。由于储水系数不存在方向性，而渗透系数有方向性，因此，在潮汐参数相位差为负值时，可利用Hsieh等(1987)的水平流模型得到水平向含水层参数；在潮汐参数相位差为正值时，可利用Roeloffs(1996)的垂向流模型得到垂直向含水层参数(Shi et al，2016)。本文则采用Venedikow调和分析法(Venedikov et al，2003)得到新10井水位的相位差为正值，因此，采用了Roeloffs(1996)的垂向流模型计算了2次地震前后含水层垂向参数的变化，如导水系数T、储水系数S和扩散系数D。

由于九寨沟M_S7.0、精河M_S6.6地震相继发生，为了分析新10井水位在地震前后井-含水层水文参数的变化，本文选取了九寨沟地震前1个月内(2017年7月8日~8月7日)和地震之后1个月内(2017年8月10日~9月7日)新10井水位数据，反演得到了2次地震前后的含水层参数变化。由于这2次地震之前的2017年7月18日科曼多尔群岛海域发生了M_S7.8地震，为了探讨该地震对含水层是否也有影响，选取了该地震前1个月的数据(2017年6月18日~7月17日)加以分析。图 6所示为Roeloffs(1996)的垂向流模型中扩散系数-相位滞后理论曲线与实际值的对比(扩散系数D与渗透系数K、储水系数S、含水层有效厚度b间的关系为：D=Kb/S)。由图 6可见，九寨沟M_S7.0、精河M_S6.6地震前的相位滞后要比地震之后的大，扩散系数要比地震之后的小，而2017年7月18日科曼多尔群岛海域M_S7.8地震对其影响不大。

利用固体潮分析法反演的新10井-含水层系统各参数见表 2。由表 2可见，在含水层储水系数不变的情况下，其扩散系数和渗透系数明显增大，而且数值减小了近1倍。

4 讨论

通过对比表 1(水震波和地震波特征)、表 2(固体潮响应参数)2种方法反演的新10井水位在地震波作用期间和地震波作用前后的含水层参数变化可以看出：①利用Cooper等(1965)的水平流模型反演得到的含水层参数，其渗透系数在精河6.6级地震时比九寨沟7.0级地震时大，表明在地震波作用下含水层的水平向渗透性明显增大；②利用Roeloffs(1996)的垂向流模型反演得到的2次地震前后的含水层参数，其渗透系数在2次地震前比地震后明显偏低，表明地震波作用导致了含水层垂向的渗透性有所减小。而用2种方法反演得到的含水层储水系数并未发生明显变化，表明地震对储水系数的影响较小。另外，由表 1、2还可见，新10井-含水层系统储水系数值均在10^-5量级内，而渗透系数在水平向和垂向存在量级上的差异，水平向为61~147m/d(表 1)，垂向为20~50×10^-5m/d(表 2)，表明新10井-含水层系统内的水流方向主要以水平流为主。

由前人研究结果可知，地震波作用可导致含水层渗透性增强(Shi et al，2016；Sun et al，2015)，本文分析结果与此相同。如图 7所示，新10井-含水层系统既有水平向的井孔-含水层间的水流交换，也有垂向的含水层-弱透水层间的水流交换，地震波作用之后，水平向的渗透系数增加，垂向的渗透系数减小，其井孔水位也出现上升变化(图 2(a))。而水平向的渗透性增加主要是由于地震波作用过程中，周期性的张压作用导致孔隙内的水流发生快速移动，使孔隙内的填充物重新排列或移动，也即当地震波穿过含水层系统时，周期性的水流作用使水流通道从堵塞变为疏通，最终导致含水层的渗透性增大。

由于新10井含水层内的水流作用主要以水平向的为主，地震波作用导致的渗透性增加也以水平向为主，在含水层储水系数变化不大的情况下，水平向的水流作用增强也会导致垂向水流交换的减弱。所以，如图 7所示，2次地震之后水平向的水流交换增强(水平向渗透系数增大)，而垂向水流交换减弱(垂向渗透系数减小)。从图 2水位曲线亦可见，精河6.6级地震后水位出现了更为明显的上升变化，表明精河6.6级地震引起了新10井含水层水平向渗透性明显增强。

储水系数和渗透系数是含水层系统中2个最重要的参数。储水系数是含水层水头变化1个单位时，从底面积为1个单位、高等于含水层厚度的柱体中所释放(或储存)的水量；渗透系数是当水力坡度为1时，水在介质中的渗流速度，它是描述介质渗透能力的参数。水震波和地震波方法(表 1)与固体潮方法(表 2)得到的储水系数在量级上一致，而渗透系数明显不同。已有研究表明，含水层渗透系数并不是恒定不变的参数，它会随外力作用于井-含水层系统程度的不同而表现出不同的变化(Candela et al，2014)。本文采用水震波与地震波的方法反演得到的渗透系数明显大于采用固体潮方法反演得到的渗透系数，这是因为地震波引起的水震波振荡周期较短，而且短时间内含水层与井孔之间的水量交换大，水流速度快，因此，反演得到的这一时段的平均渗透系数偏大；而固体潮作用于水位的周期很长，变幅小，水流速度基本比较稳定，含水层与井孔之间水量交换缓慢，所以，得到的渗透系数偏小。这说明周期性的外力作用不会改变含水层系统的储水能力，但会影响含水层系统的渗透能力，即渗流速度。另外，由于新10井成井时没有保存抽水实验数据，故无法计算其实际的渗透系数，因此，只能通过反演的方法(水震波、固体潮响应)得到其等效值。

5 结论

本文基于新10井水位对九寨沟7.0级、精河6.6级地震的同震响应认识，分别利用水震波与地震波以及固体潮调和分析法，反演了新10井水位在地震波作用期间和地震波作用前后的含水层参数，并在此基础上，分析和讨论了新10井水位同震变化特征及其可能的响应机理，得出以下初步结论。

(1) 新10井水位对九寨沟地震、精河地震的同震响应形态不同，前者表现为水位振荡—平稳，后者表现为水位振荡—上升。

(2) 利用水平流模型反演所得的新10井含水层渗透系数，在九寨沟地震时为61m/d，而在精河地震时为147m/d，表明地震波作用导致水平向的渗透性增强。

(3) 利用垂向流模型反演得到的新10井含水层渗透系数，在2次地震之前约为49×10^-5 m/d，之后约为18×10^-5 m/d，表明地震波作用导致垂直向的渗透性减弱。

致谢: 新疆维吾尔自治区地震局监测中心为本项研究提供了水磨沟台地震波形数据，这对本文研究及取得的认识至关重要，在此表示衷心的感谢。