Hydrogeological Responses to Distant Earthquakes in Aseismic Region

Abstract

For the first time precise measurements of the groundwater level variations in the territory of the Mikhnevo geophysical observatory in an aseismic region (Moscow region, Russia) have been carried out since February 2008 at a sampling rate of 1 Hz. The groundwater level variations under quasi-stationary filtration are considered indicators of the dynamic deformation of a fluid-saturated reservoir represented by carbonate-terrigenous sediments. Both permanent (long-term) factors—atmospheric pressure, lunar-solar tides, and periodic (short-term) ones—seismic impacts from distant earthquakes, are used as probing signals for analyzing the filtration parameters of aquifers of different ages. Hydrogeological responses to the passage of seismic waves from earthquakes with magnitudes of 6.1–9.1 with epicentral distances of 1456–16,553 km was recorded in 2010–2023. Dependences of dynamic variations of the pore pressure in the upper weakly confined and lower confined aquifers on the ground velocity are approximated by different regression functions. Spectral analysis of hydrogeological responses made it possible to identify coseismic and postseismic effects from distant earthquakes. The postseismic effects in the form of an episodic increase in the pore pressure may be caused by a skin effect—clogging of microcracks nearby the wellbore by colloidal particles under intensive seismic impact.

1. Introduction

Precise hydrogeological monitoring is mainly carried out in seismically active regions. According to the data records it was established that the passage of seismic waves from earthquakes leads to coseismic and postseismic underground water level and pore pressure variations [1,2]. Poroelastic deformation of a fluid-saturated reservoir is considered the main mechanism of coseismic oscillations in the far field [3]. Postseismic hydrogeological effects manifest themselves as abrupt or gradual rises or falls in the water level combined with high-frequency oscillations [4,5]. Postseismic effects in the near field at distances comparable to the fault length (L) [6] with a seismic energy density of more than 10³ J/m³ are driven by the static deformation of the fluid-saturated reservoir. This includes the cracking of rock massif, and liquefaction of unconsolidated sediments and can manifest themselves as a long-term (from 1 month or more) step-like and gradual decrease or rise in the water level [7,8,9]. Clogging and unclogging of the main fractures [10,11,12], geochemical processes, and degassing [13,14,15] can be considered as the mechanisms of postseismic variations in the intermediate field of the earthquakes (at the distances (1–10 · L) with the seismic energy density of 10³ to 10⁻³ J/m³. Postseismic hydrogeological effects are also noted in the far field (at distances over 10 · L) for earthquakes with a seismic energy density of less than 10⁻³ J/m³ and manifest themselves in different ways in the host rock and in the fault damaged zones. For example, a decrease in pore pressure by 2.8 hPa was recorded during the earthquake in Pakistan on 24 September 2013 M7.7 in the Gome 1 well at 396 m depth close to the Dead Sea Fault Transform. An increase in pore pressure by 1.2 hP was noted in the Meizar 1 well at the depth of 1250 m located within the Arabian plate [16].

Estimation of the filtration properties of a fluid-saturated reservoir based on the analysis of the phase shift between the M₂ tidal wave identified in the ground displacement and the water level indicates a change in permeability during seismic impact [17,18,19,20]. The calculated permeability changes during the passage of seismic waves differ for wells drilled in the host rock and in the fault damaged zones. For example, the permeability changes for wells within the Tancheng Lujiang fault zone (a large continental-scale active strike-slip fault zone with a length greater than 3000 km) were greater compared to the values for wells in the host rock for six earthquakes in 2007–2016 [21].

Seasonal variations in the groundwater level can lead to a disruption of the quasi-stationary filtration due to changes in the hydrostatic pressure and transmissivity of the fluid-saturated reservoir [22,23]. It should also be taken into account that the change in permeability calculated from short time series (up to 10 days) cannot be considered reliable [24]. A significant contribution to the groundwater level variations and distortion of phase shift is also made by the anthropogenic impact associated with the exploitation of the aquifer [25].

This work is devoted to determining the range of background changes in the filtration properties of the fluid-saturated reservoir and studying its deformation in weakly confined and confined conditions under the seismic impact of distant earthquakes. In this study, we generalize the results of registration of hydrogeological responses in platform conditions (within the Mikhnevo observatory) during the passage of seismic waves from distant earthquakes in 2010–2023. At first, the background characteristics of the reservoir and their variations (the reaction to a permanent impact—barometry and tides) are estimated. These data are necessary to properly assess possible changes in the filtration properties of the reservoir under the quasiperiodic impact.

2. Hydrogeological Conditions of the Area under Investigation

Two observation wells #1 (115 m deep) and #2 (60 m deep) were drilled 30 m apart at the territory of the Mikhnevo (MHV) geophysical observatory of Sadovsky Institute of Geospheres Dynamics of Russian Academy of Sciences (IDG RAS). Glacioaqueous and moraine Quaternary sediments 10.2 m thick are distributed from the surface within the study area. Unconsolidated sediments are composed of sandy loam, red loam with gruss and gravel, and lie on the eroded surface of solid rocks. Solid rocks are represented by carbonate-terrigenous strata of the Middle and Lower Carboniferous age. Deposits of the Kashirsky and Vereya horizons of Middle Carboniferous, Protvinsky, Steshevsky, and Tarussko–Oksky complexes of Lower Carboniferous are widespread in the geological section from top to bottom (Figure 1).

The Kashirsky horizon of Middle Carboniferous consists of the Lopasny and Narskay strata represented by unevenly fractured gray limestones separated by the Khatunsky horizon composed of variegated clay. Contact zones (clay-limestone) and main cracks are accompanied by increased rock fracturing.

The Vereya horizon of Middle Carboniferous 19.5 m thick is composed of dense red clay with subordinate interlayers of clayey limestone and marl. The rocks of the Vereya horizon with erosion lie on the deposits of Lower Carboniferous, which include the Protvinsky and Aleksinsky horizons separated by the Steshevsky horizon.

The Protvinsky horizon of Lower Carboniferous is represented by brown limestone with clay interbeds. The Steshevsky horizon is composed of greenish-brown clay with interlayers of marl and clayey limestone. The Aleksinsky horizon is composed of gray limestone irregularly fractured.

Two aquifers are predominantly present within the study area. The Kashirsky aquifer (upper horizon) is weakly confined (WLU), underlain by a regional aquiclude—Vereya clays. The Aleksinsko-Protvinsky aquifer (lower horizon) is confined (WLC), and has a single level surface for the rocks of Lower Carboniferous due to the unsustainable thickness of the relative aquiclude—Steshevsky clay.

The main hydrogeological parameters were determined based on the results of test pumping and flow metering in the wells [26]. The transmissivity of the Kashirsky aquifer tapped in the interval of 44.0–56.2 m in well #2 is 15 m²/day, the hydraulic conductivity is 1.2 m/day, the hydraulic diffusivity and specific yield are 4.7·10² m²/day and 0.02, respectively. Increased rock fracturing was identified in the intervals of 43.5–43.9, 44.7–45.2 and 47.7–48.0 m according to the resistivity and caliper logs.

In well #1 that taps the Aleksinsko–Protvinsky aquifer in the interval of 92–115 m, the head reaches 24 m, the transmissivity does not exceed 3 m²/day, the hydraulic conductivity is 0.13 m/day, and the hydraulic diffusivity and storativity are 1.3·10⁴ m²/day and 2.3·10⁻⁴. Water inflow intervals were identified at the depths of 92–94 m and 99–100 m according to flow measurement, transmissivity for these intervals increases up to 5 m²/day.

The hydrographic network belongs to the basin of the Oka river located at a distance of 14 km south of the territory of the Mikhnevo observatory. The hydrological regime of the Oka river is of flat type and is characterized by a high spring flood according to the gauging station located near the town of Serpukhov (data provided by the Hydrology Department of the Central Administration of the Hydrometeorological Service of the Russian Federation) (Figure 2a). The climate is mildly continental, with cold winters and warm summers. The amplitude of annual variations of the atmospheric pressure varies from 43 to 90 hPa. The maximum atmospheric pressures up to 1040 hPa are observed in winter, and the minimum ones—up to 970 hPa—in summer.

The regime of the Kashirsky aquifer is disturbed due to aquifer exploitation. Local mini-depression cones are identified in the regional decrease trend (Figure 2b). Cones are associated with periodic water pumping (with a duration of ~3 h every ~6–7 days) from a technical well located 300 m downstream in the direction of the underground flow [27]. Nevertheless, an annual cyclicity is traced in the level variations in the form of intensive rises during the floods (in spring) with amplitudes up to 1.9–2.3 m.

The Aleksinsko–Protvinsky aquifer has longer periods of seasonal recharge. The maximum water level of the confined aquifer was established in June 2021. The minimum value was noted during 2011 and associated with the previous dry year 2010 and the lowest hypsometric position of the level in the Oka river (Figure 2a). The amplitude of annual variations of the confined aquifer varies from 1.1 to 3.1 m. These changes are caused by the hydraulic connection with the river. The river’s valley erodes the top of Lower Carboniferous deposits [22].

A common hydrodynamic regime of the fluid-saturated reservoirs demonstrating annual cyclicity is clearly identified in both water levels of the confined and weakly confined aquifers (Figure 2c).

The periods of spring rise and low-water decrease in the confined aquifer are traced 1–2 weeks ahead of the weakly confined aquifer. Relatively low transmissivity of the water-bearing rocks leads to smooth and gradual dynamic deformation of the fluid-saturated carbonate reservoir.

3. Methods

The instrumentation and measuring complex in the Mikhnevo observatory includes water level sensors LMP 308i (BD Sensors, Germany); atmospheric pressure sensor PAA-33X (Keller, Switzerland); broadband seismometers STS-2 (Streckeisen, Switzerland, sensitivity 1500 V/(m/s)) and CM-3-E (Russian Federation, sensitivity 2000 V/(m/s)). The registration of signals at a sample rate of 200 Hz is carried out by ADC E14-440 (L-Card, Russian Federation). Observation wells #1 and #2 are equipped with submersible level probes LMP308i. The level sensor in well #1 is operating since February 2008. Well #2 was equipped with the level probe in July 2013. In 2013–2023 the resonant frequencies of well #1 were 0.08–0.09 Hz (natural oscillation period 11–13 s) with level variations of 64.1–71.9 m, the resonant frequencies of well #2 were 0.15–0.19 Hz (natural oscillation period 5–7 s) with level variations 42.6–46.5 m. The atmospheric pressure was recorded from February 2008 to July 2019 by a digital weather station Vantage-Pro2 (Davis, USA). In August 2020 it was replaced by a PAA-33X sensor.

Hydrogeological and barometric data form a database covering the years 2008–2023. Processing the obtained data is aimed at monitoring the regime of aquifers of different ages. The barometric efficiency is calculated using the widely used linear regression method between the underground water level variations and the atmospheric pressure according to a geostatic model [28,29]. The ETERNA 3.0 software [30] was used to estimate the theoretical volumetric deformation and vertical ground displacement taking into account the coordinates of the observation point (51°16′58.7″ N 37° 35′09.7″ E). To isolate tidal waves in hydrogeological and barometric data, bandpass filtration was performed in the ranges of 12.32–12.52 h and 25.7–25.9 h after bringing the time step to 300 s.

To estimate the phase shift between the tidal component M₂ identified in the ground displacement and water level, a method based on the construction of phase trajectories in the coordinates “ground displacement—water level changes” in the form of ellipses was applied [22,31]. The obtained values of the phase shift were used to calculate the transmissivity according to the poroelastic model [3].

According to the model [3] equation for fluid flow in a confined aquifer (homogenous and isotropic) is:

$$\frac{{\partial }^{2}s}{\partial r^2}+\frac{1}{r}\frac{\partial s}{\partial r}-\frac{S}{T}\frac{\partial s}{\partial t}=0,$$

(1)

where s—drawdown, m; r—radial distance from the center of well, m; S—storage coefficient; T—transmissivity of an aquifer, m²/day; t—time, day.

Boundary conditions for periodic discharge from the aquifer to well with a volumetric rate $Q_0\exp \left(i\omega t\right)$ are:

$$\begin{array}{l}2\pi r_{w}T{\left(\frac{\partial s}{\partial r}\right)}_{r=r_w}=-Q_0\exp \left(i\omega t\right),\\ s\to 0,\hspace{0.17em}r\to \infty .\end{array}$$

(2)

The analytical solution for this equation allows us to express transmissivity (T) from phase shift (η) [3]:

$$\begin{array}{l}\eta =-t{g}^{-1}\left(\frac{F}{E}\right),\\ E\approx 1-\frac{\omega r_c^2}{2T}Kei\left({\alpha }_w\right),\\ F\approx \frac{\omega r_c^2}{2T}Ker\left({\alpha }_w\right),\end{array}$$

(3)

where ω = 2π/τ—fluctuation frequency, s⁻¹; τ—fluctuation period, s; r_c—casing radius, m; r_w—well radius, m; T—transmissivity, m²/day; Ker(α_w) & Kei(α_w)—Kelvin functions of 0th order; S—storage coefficient; α_w = (ωS/T)^1/2 ∙ r_w.

Figure 3 shows the relation between phase shift η and transmissivity T calculated for wells #1 and well #2.

We use casing radiuses r_c1 = 0.0635 m, r_c2 = 0.0585 m; well radiuses in the water intake interval r_w1 = 0.059 m, r_w2 = 0.056 m; storage coefficients S₁ = 2.3 × 10⁻⁴, S₂ = 2.9 × 10⁻⁴; watered intervals d₁ = 23 m, d₂ = 13.9 m.

From the curves in Figure 3, the ranges of transmissivity values T were determined, which correspond to the range of phase shift η ± σ calculated by the method of phase portraits [31], where σ is the standard deviation, which determines the calculation error.

The permeability of a fluid-saturated reservoir (k) can be estimated by the formula [32]:

$$k=\frac{\mu }{\rho gd}T,$$

(4)

where T is the transmissivity of rocks, m²/day; μ is the dynamic viscosity, μ = 1.78 × 10⁻³ Pa∙s; ρ is the water density, kg/m³; g is the acceleration of gravity, m/s²; d is the thickness of the reservoir, m.

At the next stage, data segments corresponding to the excess of the threshold water level change of 5 cm/day are excluded [22]. To estimate the background annual variations of the parameters, filtering was applied using the sliding window method (60 days, overlap 50%).

Additionally, the porosity of the fluid-saturated reservoir (n) was estimated in accordance with the relation [33]:

$$n=\frac{B_e{S}_s}{\beta \rho },$$

(5)

where B_e—barometric efficiency, S_s—specific storage, m⁻¹, S_s = 3.53 × 10⁻⁶ m⁻¹ for a confined aquifer, S_s = 4.97 × 10⁻⁶ m⁻¹ for a weakly confined aquifer [34], β—water compressibility, β = 4.6 × 10⁻¹⁰ m²/N, ρ—rock density, ρ = 2700 kg/m³.

A sample of hydrogeological and seismic records corresponding to the passage of seismic waves from distant earthquakes with a sampling rate of 1 Hz has been formed. The sample included earthquakes for which the peak-to-peak amplitude of the ground velocity is greater than 0.03 mm/s. For each earthquake, a segment lasting 6 h is created (3 h before and 3 h after the arrival of the seismic wave from the earthquake to the MHV station).

In the groups of surface waves identified in the ground velocity and water levels of the confined and weakly confined aquifers, the maximum amplitudes were determined. The amplitudes were measured between the successive maximum and minimum in a sliding window of 72 s with an overlapping of 50%. To remove the influence of local features of the measuring point, normalized spectra were constructed—the ratio of the spectra of the ground velocity and water levels constructed for 3-h intervals after the arrival of a seismic wave to the corresponding noise spectra for 3-h intervals before the arrival of a seismic wave to the MHV station.

The main parameters of seismic events—the arrival time of P-wave to the Obninsk station (OBN), the coordinates and depth of the earthquake source were taken from the catalog of the Federal Research Center “Geophysical Survey” of Russian Academy of Sciences (55.1146° N, 36.5674° E). Seismic station OBN, located at a distance of ~80 km from the Mikhnevo observatory, is equipped with an STS-2 seismometer, and it is a part of the international seismic network IRIS (http://ds.iris.edu/mda/II/OBN, accessed on 1 February 2023). The data recorded by the MHV and OBN seismic stations from remote earthquakes can be considered the results of measurements in a single observation point. For distant earthquakes that occurred from 2010 to 2019 the data from the MHV seismic station were used. For events in 2020–2023 (and four earthquakes that occurred on 11 April 2012, 14 May 2019, 25 June 2019, and 15 December 2019) the data of the OBN station were considered.

The probabilistic approach [35] was used to analyze the background characteristics of seismic, barometric, and hydrogeological data. The probabilistic approach is usually used to characterize the background seismic noise [35,36]. In this study, we apply it also to analyzing background variations in barometric and hydrogeological data. The initial data used for analysis were divided into intervals of 2-h duration each with a 1-h overlap. The linear trend was removed preliminarily (the least squares regression line was estimated and subtracted from the initial data). The power spectral density (PSD) was calculated for each 2-h interval with the Welch method in a running window of the length of 2/13 h with a 75% overlap. The obtained data of PSD were used to calculate the probability density for each frequency. For this purpose, at first the PSD was averaged in a running window of the length of 1 octave with an overlap of 1/16 octave. This allowed us to reduce frequencies and obtain an even sampling over frequency in the logarithmic scale. Then the probability densities of PSD distributions with the width of interval 1 dB were plotted for each frequency. Then resting on the calculated values of probability density for all the frequencies, the spectrogram of probability density function (PDF) of the occurrence of a certain value of PSD at each frequency was plotted. A statistical mode, the 10th and 90th percentile is used for further consideration. A comparative analysis of the PSD of seismic signals and water levels (pore pressure) for the most significant earthquakes with the background variations presented on the PSD diagrams was performed.

To compare the hydrogeological effects, the density of seismic energy (e) can be estimated as follows [37]:

$$M=2.7+0.69\lg e+2.1\lg r,$$

(6)

where r is the epicentral distance, km; M is the magnitude.

4. Results

4.1. Background Characteristics

Quasi-stationary factors—the atmospheric pressure and Earth’s tides, which manifest themselves in the range of periods from hours to several days, control the dynamics of poroelastic deformation of the fluid-saturated carbonate reservoir. The barometric efficiency was determined for long-term observations of the groundwater levels and atmospheric pressure. The barometric efficiency (Be) is −0.76 for the upper weakly confined aquifer associated with a fractured cavernous carbonate reservoir for 2013–2022. For the lower confined aquifer tapped in unevenly fractured limestones with subordinate clay interlayers it is less than −0.64 for 2010–2022.

The coherence functions between the atmospheric pressure and WLC, WLU aquifers clearly show a synchronous increase in the range of periods of 3.5, 6, 8, and 12 h (Figure 4a,b). The coherence at the period of ~12 h is also traced between the ground displacement and the water levels for semidiurnal tidal wave S₂ of the solar-type (Figure 4c,d). The coherence function between the ground displacement and the WLC aquifer is close to 1 for the period of 12.42 h, corresponding to the semidiurnal tidal wave M₂ of lunar type. For the ground displacement and the WLU aquifer, the coherence function at the periods of the M₂ wave is less pronounced and less than 0.2.

The carbonate fluid-saturated reservoir is affected by the Earth’s tides with different intensities in confined and weakly confined conditions [38]. The tidal factor (ratio of the amplitudes of the main types of tidal waves identified in the hydrogeological data to volumetric deformation) of the WLC aquifer is equal to 0.36 mm/nano, for the WLU—0.20 mm/nano. A technique based on the analysis of the phase shift between the M₂ tidal wave identified in the ground displacement and water levels was used to assess the filtration properties of the fluid-saturated reservoir (Figure 5). Additionally, the external factors (natural and anthropogenic) which can affect the regime of the confined and weakly confined aquifers, are taken into account.

At first, intervals with a non-quasistationary filtration with daily level variations exceeding ±5 cm/day (Figure 5) were excluded from the calculated phase shift series. A selection of intervals for data processing for the WLU aquifer was made corresponding to the level position relative to the clay interlayers, which are characterized by an increase in the amplitudes of the M₂ tidal wave up to the maximum values. Taking into account the factors noted above, the filtration parameters of the WLU aquifer (phase shift, transmissivity, permeability) can be determined only for four intervals (

Table 1. Background characteristics of the fluid-saturated carbonate-terrigenous reservoir.

Year	Month	φ°	T, m2/day	k × 10−13, m2	Month	φ°	T, m2/day	k × 10−13, m2
	Confined Aquifer				Weakly Confined Aquifer
2011	02–06	−3.3–−1.3	2.6–7.6	5.1–14
2012	05–08	−2.4–1.5	3.8–6.1	7.4–12
2013	06–08	−5.5–−2.9	1.5–3.1	2.9–6.0
2014	05–07	−7.5–−2.7	1.1–3.3	2.0–6.3
2015	05–08	−6.2–−4.5	1.3–1.9	2.5–3.7	07–09	−11–−2.5	0.8–4.0	1.6–7.9
2016	05–09	−7.6–−4.9	1.0–1.7	2.0–3.4	07–09	−9.5–−2.0	0.9–5.0	1.8–9.9
2017	04–08	−6.0–−7.1	1.1–1.4	2.2–2.6
2018	05–09	−6.8–−5.8	1.4–1.2	2.3–2.7
2019	03–08	−10.2–−7.6	1.1–1.2	1.5–2.0	07–08	−10.7–−2.0	0.9–5.0	1.8–9.9
2020	05–09	−10.4–−3.7	1.2–2.3	1.4–4.5	06–08	−4.0–−2.6	3.0–4.0	5.9–7.9
2021	05–09	−6.8–−5.0	1.2–1.7	2.3–3.2

). Red numbers in Figure 5b correspond to the earthquakes: 1–11 March 2011 Tohoku, Japan; 2–11 April 2012 off the western coast of North Sumatra; 3–16 September 2015 Chile; 4–8 September 2017 Mexico; 5–28 September 2018 Minahassa Peninsula, Sulawesi; 6–25 October 2018 Ionian sea; 7–22 February 2019 Peru-Ecuador; 8–22 July 2020 Alaska Peninsula. Numbers 3–8 in Figure 5b correspond to the numbers in the first column in

Table 2. Earthquakes with synchronous hydrogeological responses in the confined and weakly confined aquifers.

(Number in Figure 5)	Date and Arrival Time of P-Wave to OBN Station	Depth, km	Magnitude Mw	Epicentral Distance r, km	Seismic Energy Density e∙10−4, J/m3	Peak-to-Peak Vertical Ground Velocity Vp-p, mm/s	Peak-to-Peak Pore Pressure of WLU, Pa	Peak-to-Peak Pore Pressure of WLC, Pa	Region
	Mediterranean seismic belt
6	25.10.2018 22:59:29.6	16	7.1	2333	1.34	1.37	445	24	Ionian sea
-	06.02.2023 01:21:42.1	15	7.8	1934	24.49	2.18	554	59	Turkey
-	06.02.2023 10:28:45.9	12	7.7	1873	19.34	2.56	1447	78	Turkey
	West Pacific seismic belt
5	28.09.2018 10:15:24.8	12	7.3	9519	0.04	0.12	60	53	Minahassa Peninsula, Sulawesi
	East Pacific seismic belt
3	16.09.2015 23:11:29.5	17	8.3	14020	0.31	0.58	226	29	Chile
4	08.09.2017 05:02:50.4	15	8.2	11000	0.47	2.46	1084	78	Mexico
7	22.02.2019 10:31:16.46	121	7.5	11764	0.04	0.12	31	29	Peru-Ecuador
8	22.07.2020 06:23:48.7	37	7.8	7702	0.04	0.59	226	20	Alaska Peninsula

.

Figure 5. Diagrams of the daily level alteration (a), phase shift (b), transmissivity (c), permeability (d). Red numbers in (b) correspond to the earthquakes: 1–11 March 2011 Tohoku, Japan; 2–11 April 2012 off the western coast of North Sumatra; 3–16 September 2015 Chile; 4–8 September 2017 Mexico; 5–28 September 2018 Minahassa Peninsula, Sulawesi; 6–25 October 2018 Ionian sea; 7–22 February 2019 Peru-Ecuador; 8–22 July 2020 Alaska Peninsula. Red numbers 3–8 in Figure 5b correspond to the numbers in the first column in Table 2. Dotted lines in (b–d)—error bars.

Figure 5. Diagrams of the daily level alteration (a), phase shift (b), transmissivity (c), permeability (d). Red numbers in (b) correspond to the earthquakes: 1–11 March 2011 Tohoku, Japan; 2–11 April 2012 off the western coast of North Sumatra; 3–16 September 2015 Chile; 4–8 September 2017 Mexico; 5–28 September 2018 Minahassa Peninsula, Sulawesi; 6–25 October 2018 Ionian sea; 7–22 February 2019 Peru-Ecuador; 8–22 July 2020 Alaska Peninsula. Red numbers 3–8 in Figure 5b correspond to the numbers in the first column in Table 2. Dotted lines in (b–d)—error bars.

The transmissivity of the fluid-saturated carbonate reservoir in the upper part of the section varies from 0.8 to 5.0 m²/day, in the lower part of the section—from 1.0 to 1.9 m²/day for the same time intervals in Table 1. The reservoir permeability represented by Lower Carboniferous deposits has a lower value compared to the overlying stratum.

The porosity of two aquifers of different ages was estimated using the values of the barometric efficiency and specific storage, which is inversely proportional to the tidal factor, according to (5) [33]. The porosity of the upper part of the section of Middle Carboniferous is 31% and is consistent with the data of laboratory determinations of the porosity of limestones in well #2 (Figure 1). The porosity of the lower part of the section of Lower Carboniferous does not exceed 14%. The filtration parameters of the WLC aquifer have lower values in comparison with the WLU aquifer.

The reactions of aquifers of different ages tapped at different depth intervals in close wells (two wells 30 m apart) to quasi-stationary long-period factors (atmospheric pressure and Earth’s tides) differ. They depend on the filtration parameters of the aquifer.

4.2. Hydrogeological Responses to Distant Earthquakes

Hydrogeological responses to distant earthquakes recorded at the Mikhnevo observatory can be grouped according to the localization of the earthquake sources in seismic belts (Figure 6). Responses to 19 earthquakes that occurred at the epicentral distances of 1456–2575 km within the Mediterranean belt were recorded in 2011–2023. All earthquakes are shallow (with a depth of 12–26 km) with magnitudes of M_w 6.1–7.3.

The amplitudes of the pore pressure in the confined aquifer vary from 24 to 445 Pa, in the weakly confined one—from 24 to 84 Pa, with the peak ground velocity being 0.06–1.51 mm/s (Figure 7). Hydrogeological responses of the fluid-saturated reservoir to the two earthquakes in Turkey on 6 February 2023 are an exception and will be discussed in more detail below.

During the passage of seismic waves from 13 earthquakes with magnitudes of M_w 6.4–7.9, which occurred at the epicentral distances of 3206–5571 km within the Trans-Asian seismic belt, variations in the pore pressure of the WLU aquifer were 43–373 Pa at a ground velocity of 0.08–0.87 km/s. In one case a response of the WLC aquifer with an amplitude of 62 Pa was recorded during the earthquake in Pakistan on 16 April 2013 (M 7.7) with the peak ground velocity being equal to 1.01 mm/s (Figure 7).

The largest number of the hydrogeological responses (39 in WLU and 6 in WLC) was traced during earthquakes with magnitudes of M_w 6.2–9.1 at epicentral distances of 6492–16553 km within the West Pacific seismic belt. The amplitudes of pore pressure variations in the upper horizon did not exceed 20–225 Pa; the ground velocity was 0.04–0.58 mm/s. The maximum amplitude of the hydrogeological response of the lower aquifer up to 407 Pa was established at the ground velocity of 3.78 mm/s for the Tohoku earthquake, Japan, M 9.1 [39].

During earthquakes with magnitudes of M_w 7.0–8.3 located within the East Pacific seismic belt at the epicentral distances of 7067–15,079 km, variations in the pore pressure in the WLU aquifer vary from 22 to 226 Pa at the ground velocity of 0.06 to 0.63 mm/s; in the WLC aquifer—from 20 to 146 Pa at the ground velocity of 0.12–1.79 mm/s. The maximum amplitude of the water level variations of the WLU aquifer of 1.1 hPa was recorded at the ground velocity of 2.46 mm/s during the passage of seismic waves from the earthquake off the coast of Chiapas, Mexico with the magnitude of 8.2 on 8 September 2017 [40].

The analysis of hydrogeological responses indicates that the largest amplitudes of the pore pressure variations during the passage of seismic waves from distant earthquakes were recorded in the upper aquifer. Of particular interest is the study of the frequency range of dynamic deformation of the fluid-saturated carbonate-terrigenous reservoir, which is located in two different hydrogeological conditions (confined and weakly confined ones) under the same type of seismic impact. Simultaneous responses of the water level recorded at the Mikhnevo observatory were traced both in the WLC and WLU aquifers from eight earthquakes that occurred in the Mediterranean, West, and East Pacific seismic belts (Table 2).

The responses of the fluid-saturated reservoir to the earthquakes that occurred on 28 September 2018 near the Minahassa Peninsula, Sulawesi and on 22 February 2019 near the Peru-Ecuador border area are weakly expressed. Hydrogeological responses to a medium-focus earthquake compared to a small-focus earthquake are two times less at the same values of ground velocity of 0.12 mm/s, seismic energy density of 0.04∙10⁻⁴ J/m³, and oscillation period of 36 s. An increase in the intensity of the pore pressure variations in the normalized spectra was noted only during an earthquake near the Minahassa Peninsula, Sulawesi. The extremum of water level of the weakly confined aquifer corresponds to the frequency range of 0.04–0.06 Hz, of the confined aquifer—0.02 Hz.

The amplitude of the pore pressure in the upper aquifer does not exceed 226 Pa at the ground velocity of 0.58–0.59 mm/s recorded from earthquakes on 16 September 2015 near the central coast of Chile and 22 July 2020 near the Alaska Peninsula. The amplitude of the hydrogeological response of the lower aquifer to the earthquake near the Alaska Peninsula is 1.5 times less; the pressure head increased by 3.5 m over 5 years of observations (Figure 2). The maximum of the ground velocity and pore pressure in the WLU aquifer in the normalized spectra for the two earthquakes are similar. The peaks stand out in the range of 0.04–0.05 Hz corresponding to the arrival of body and surface waves during the earthquake that occurred on 22 July 2020. During the earthquake on 16 September 2015 two frequency ranges of 0.04–0.05 Hz and 0.007–0.009 Hz are traced. Identified frequency ranges correspond to the propagation of reflected and refracted PKIKP waves [41]. The extremes of the pore pressure variations in the WLC aquifer are shifted to the low-frequency area in the range of 0.004–0.006 Hz during the earthquake near Chile and to 0.004–0.005 Hz for the earthquake near Alaska.

For the earthquake in the Ionian Sea on 25 October 2018 with a ground velocity of 1.37 mm/s the peaks of the pore pressure variations in the WLU aquifer in the frequency range of 0.03–0.08 Hz coincide with an increase in the intensity of the ground velocity in the range of 0.03–0.08 Hz in the normalized spectra. An increase in the intensity of the hydrogeological response in the WLC aquifer was noted at the lower frequency of 0.02 Hz. The normalized spectra of the pore pressure variations in the upper and lower aquifers and ground velocities are similar during the catastrophic earthquake on 8 September 2017 near the coast of Chiapas, Mexico with the magnitude of M_w 8.2. Extreme values of the ground velocity appear at frequencies of 0.03 and 0.04 Hz, and extreme values of the water level—in the frequency range of 0.02–0.03 and 0.04 Hz.

Hydrogeological responses to two earthquakes that occurred in Turkey on 6 February 2023 with magnitudes of Mw 7.8 and 7.7 were traced in both confined and weakly confined aquifers. The amplitude of the pore pressure in the upper aquifer was 554 Pa with a period of 41 s, in the lower one—59 Pa with a period of 22 s during the passage of seismic waves from the first earthquake through the carbonate-terrigenous reservoir (Table 2). During the second earthquake, which occurred 9 h later, the amplitude of the pore pressure in the WLU aquifer increased to 1.4 hPa, and in the WLC to 78 Pa with the same oscillation period of 24 s (Figure 8). The water level response to the third earthquake that occurred in Turkey on 20 February 2023 with a magnitude of 6.4 and ground velocity of 0.09 mm/s was traced only in the WLU aquifer with an amplitude of 40 Pa.

To determine the postseismic effects, the response of a fluid-saturated reservoir was considered on a daily data series filtered from the influence of atmospheric pressure (Figure 8d). A smooth postseismic increase in pressure was observed in the lower confined aquifer after the first earthquake by 256 Pa and after the second one by 18 Pa against the clearly pronounced effect of Earth’s tides on pore pressure. On the contrary, Earth’s tides are not expressed in weakly confined aquifers. Only after the second earthquake in Turkey was a postseismic increase in pore pressure in the WLU aquifer of 65 Pa traced for 3 h.

The peaks of the ground velocity are identified in the ranges of 0.03–0.05 and 0.03–0.06 Hz on the normalized spectra constructed from the registration data of the first and second earthquakes, respectively (Figure 9). The peak of the pore pressure variations in the upper aquifer appears at frequencies of 0.08 and 0.05 Hz. In the lower aquifer, they are shifted to the low-frequency area and are recorded at frequencies of 0.05 and 0.03 Hz.

Differences in filtration properties of the fluid-saturated carbonate-terrigenous reservoir demonstrated in paragraph 4.1 can be considered as a “band pass” filter in the reaction of the two aquifers to seismic action. High-frequency responses (with a period of less than 60 s) are mainly registered in the upper weakly confined aquifer tapped by well #2. Low-frequency responses (with periods exceeding 60 s) are registered in the lower confined aquifer tapped by well #1. Previously the authors of this paper used this approach to typify the reaction of an aquifer to seismic actions [26]. The established difference in filtration properties of the reservoir also manifests at a higher scale level, in the reaction to the quasi-stationary factors. Hydrogeological responses to the Earth’s tides can be clearly detected in the lower confined aquifer, unlike in the upper one.

5. Discussion

An analysis of the results of the water level monitoring has shown a difference between the response of the WLU and WLC aquifers to the passage of seismic waves from distant earthquakes. The reaction of the fluid-saturated reservoir to the seismic impact recorded at the Mikhnevo observatory was detected in a wide range of seismic energy densities from 10⁻⁷ to 30.7 × 10⁻⁴ J/m³ (Figure 10). Seismic energy density is a general parameter that shows a relation between magnitude and epicentral distance.

In platform conditions, the energy density of more than 0.5 × 10⁻⁵ J/m³ and a ground velocity of 0.12 mm/s recorded at the Mikhnevo observatory are considered as the threshold of seismic impact. Hydrogeological responses to such action can be recorded in two aquifers and indicate poroelastic deformation of the fluid-saturated carbonate-terrigenous reservoir. However, the frequency extrema of hydrogeological responses of the upper and lower aquifers (WLC and WLU) differ for earthquakes with ground velocities from 0.58–0.59 to 2.5 mm/s. The extrema of the water level spectra of the WLC aquifer are shifted to the low-frequency area if compared to the WLU one. The frequency extrema of the hydrogeological responses for two aquifers are similar for strong seismic impacts with ground velocities exceeding 2.6 mm/s.

Two-hour data series of earthquakes that occurred in Turkey on 6 February 2023 were used to calculate the PSD of the ground velocity. The estimated value was compared to the probability density function (PDF) of PSD for OBN data from 1 December 2022 to 18 February 2023 (Figure 11a). Two extrema were identified in the mode of PDF of PSD in the range of ~0.06–0.08 Hz (12.5–20.0 s) and ~0.1–0.5 Hz (2–10 s). The peak of PSD at the frequencies of 0.1–0.5 Hz corresponds to secondary microseisms. Such signals are recorded everywhere. They are the main component of noise on Earth. The extremum in the range ~0.06–0.08 Hz (periods 12.5–16.7 s) is probably related to the propagation of seismic waves from earthquakes in Turkey on 6 February 2023. The maximum PSD of the ground velocity for earthquakes is detected in the frequency range of 0.07–0.08 Hz.

The PDF of the WLU aquifer has a wider range of PSD of 17 dB between 10 and 90 percentile relative to the WLC aquifer (Figure 11b,c). The mode of the WLU aquifer in the frequency range from 0.001 to 0.6 Hz varies from 15 to −12 dB, and the mode of the WLC aquifer—from 10 to −5 dB. The peak of PSD of WLC for the first earthquake in Turkey on 6 February 2023 is shifted to the low-frequency area of 0.02–0.03 Hz relative to the maximum PSD of WLU in the range of 0.06–0.07 Hz at the ground velocity of 2.2 mm/s. The hydrogeological responses in the WLC and WLU aquifers to the second earthquake in Turkey at the ground velocity of 2.5 mm/s are identified in the same frequency range as PSD.

The response of the fluid-saturated reservoir under confined and weakly confined conditions depends on hydrogeological parameters. The WLC aquifer is characterized by relatively low values of hydraulic conductivity and storativity compared to the WLU aquifer. The results of resistivity measurements confirm a 1.5 times difference in water inflow intensity between the identified intervals in the Kashirsky aquifer compared to the Aleksin–Protvinsky one according to geophysical logging [26].

Coseismic responses to the passage of seismic waves from distant earthquakes through an aquifer with low transmissivity are less pronounced. This effect corresponds to results obtained by Roellofs (1996) [42], who noted the lack of coseismic response of the fluid-saturated reservoir with low permeability to seismic action. Results of explorations performed worldwide have shown that amplitudes of hydrogeological responses depend on the filtration parameters of aquifers of different ages, as well as on the geological and structural conditions in the region [15,43], etc. Coseismic changes recover faster in an unconfined aquifer if compared to a confined one. Variations in pore pressure depend on the compressibility of rocks. The results of the hydrogeological analysis at seismic events that occurred on the same day can affect the permeability of the fluid-saturated reservoir. This is supported by studies evaluating the coseismic responses recorded in well X10 (south of Urumqi in the Xinjiang Autonomous Region) during the paired earthquakes on 8 December 2016 and 8 August 2017 [44].

The upper weakly confined aquifer in this study is characterized by higher fracturing compared to the lower confined one. The postseismic increase in the pore pressure by 45 Pa lasting 3 h in WLU was established only during the second earthquake that occurred in Turkey on 6 February 2023 at the ground velocity of 2.2 mm/s and seismic energy density of 2.4 × 10^–3 J/m³ (Figure 7). The postseismic effect in WLC was established at the ground velocity of 1.8–3.8 mm/s and seismic energy density of (0.1–3.1) × 10^–3 J/m³ from distant earthquakes on 2 February 2010 near the coast of Central Chile [41], 11 March 2011 Tohoku, Japan [39], 11 April 2012 off the western coast of North Sumatra [40], and 6 February 2023, Turkey.

The response of the fluid-saturated reservoir depends on the amplitude-frequency composition of seismic waves of an earthquake and affects the intensity of postseismic effects in the aquifer of carbonate rocks [12]. For example, after the Wenchuan earthquake M 7.9 on 12 May 2008 an increase in permeability was traced by the tidal data, which then returned to the initial value over a year. After the Tohoku earthquake M 9.0 on 11 March 2011, an analogous increase in permeability was observed in the same well, but no recovery occurred. Analyses of seismograms recorded nearby the well have shown that the main frequency of Rayleigh waves of Tohoku was five times lower than that of the Wenchuan. The authors suggest that Rayleigh waves of lower frequency promoted an efficient unclogging of large fractures in the carbonate reservoir.

At the same time in the far field, a difference in hydrogeological responses recorded in the host rock and in the fault zones was detected. An increase in pore pressure (a decrease in permeability) has been observed in the wells located mainly in host rock. On the contrary, a decrease in pore pressure (an increase in permeability) has been detected in fault zones [16,21], etc.

Laboratory observations also show that permeability can both increase and decrease under the action of dynamic stress [45,46], etc. Results of laboratory tests support the possibility for dynamic pulses to destroy low-permeable colloid barriers that emerge in rock fractures during the sedimentation of micro-particles [10].

Although Figure 5 shows no clear changes in the filtration properties of the aquifer, data processing demonstrates that only strong impacts lead to postseismic effects in the form of gradual rises of the pore pressure of the aquifer (Figure 8). For example, a postseismic increase in pore pressure with a duration of less than 1 h was observed in the lower confined aquifer after the two earthquakes in Turkey and in the weakly confined aquifer after the second earthquake. Such short-term changes cannot be detected according to the poroelastic model on the base of tidal analysis. The noted postseismic effects in the form of an episodic increase in the pore pressure may be caused by a skin effect—clogging of microcracks in the near-wellbore by colloidal particles under intensive seismic impact.

6. Conclusions

Studying the seismic impact on the filtration properties of the fluid-saturated reservoir, it is necessary to take into account the background variations of “input” parameters.

Here, we have considered two different aquifers belonging to the same fluid-saturated carbonate-terrigenous reservoir.

The lower confined aquifer is weakly fractured, characterized by lower permeability of (1.4–4.5) × 10⁻¹³ m² compared to the upper unevenly fractured, cavernous weakly confined aquifer with a permeability of (1.6–9.9) × 10⁻¹³ m² within the same time intervals

The responses of two aquifers of different ages tapped in the intervals of 42.7–56.6 m and 92–115 m to quasi-stationary and periodic impacts differ. The upper weakly confined aquifer subjected to the anthropogenic load is characterized by higher sensitivity to atmospheric pressure. The effect of Earth’s tides is weak. The lower confined aquifer is characterized by low values of barometric efficiency and increased sensitivity to Earth’s tides.

For the first time, the analysis of the seismic impact from distant earthquakes on the fluid-saturated carbonate-terrigenous reservoir was performed taking into account the alteration of filtration parameters over the section (with depth). The decrease in permeability of the fluid-saturated reservoir with depth determines the differences in the intensity of the hydrogeological responses of the WLC and WLU aquifers under the same seismic impact.

The boundary conditions of the dynamic deformation mode of the terrigenous-carbonate reservoir during the passage of seismic waves are determined. Firstly, dynamic deformation of the fluid-saturated carbonate-terrigenous reservoir was established at the ground velocity of more than 0.12 mm/s in the aseismic region. Secondly, the value of the maximum ground velocity (more than 2.6 mm/s) was determined. Actions exceeding this threshold result in a common poroelastic reaction of the carbonate-terrigenous reservoir under confined and weakly confined conditions. Thirdly, the values of the ground velocity and seismic energy density are identified, at which the postseismic effect is traced. An increase in pore pressure during seismic action can be caused by a complex mechanism: a combination of the dynamic deformation of the fluid-saturated terrigenous-carbonate reservoir with the skin effect—local clogging of microcracks by colloidal particles.