Seismic Spectral Parameters Fitting Analysis of Reservoir Area in Western China

Abstract

In this study, the seismic events (>M_L 2.0 earthquake motion records recorded by China earthquake net stations) in Western China are calculated and the source parameters of several hundreds of seismic events are obtained, such as corner frequency f_c, zero frequency limit Ω₀, seismic moment M₀, stress drop ∆σ, and apparent stress σ_app. The 95% confidence interval of each parameter of each earthquake is obtained. The sample area is about 100 km × 400 km. Then, the calculation results are compared with statistical results, which shows that the correctness and feasibility of the fitting analysis are feasible and correct. Finally, several conclusions and prospects are presented and discussed.

1. Introduction

For the actual location of earthquake occurrence and energy release, the study of source is an important part of understanding the earthquake mechanism. With the deepening of the study on the source mechanics process, the parameters needed to describe the source model gradually increase, such as seismic moment, fault scale, fracture velocity, stress drop, etc. These parameters are introduced to describe the characteristics of the source population. At present, the study of the earthquake source mainly depends on the most direct record of the earthquake wave. On the other hand, there is a distance between the seismic station and the earthquake source; thus, the seismic wave recorded by the station is not a direct reflection of the seismic source, but includes the medium information in the propagation process and the response of seismometer instruments. To this end, it is necessary to remove the influence of seismic recording instruments and seismic propagation paths on seismic waves to obtain source information from seismic waves, which are obtained from earthquake net stations. Since the instrument response of seismic records is, in fact, known, how to accurately obtain the seismic source parameters also comes down to how to effectively remove the seismic propagation path effect from the seismic wave records.

Obviously, it is easier to deduct the influence of the propagation medium and seismograph instrument from seismic wave records in the frequency domain. The deduction of the influence of the recording instrument can directly remove the frequency response of the instrument from the seismic spectrum record, while the deduction of the influence of the propagation medium needs to consider different methods according to different stations and seismic characteristics. When there is a large number of near-station distribution, the inelastic attenuation of the medium can be ignored and the theoretical geometric attenuation model R-1 for seismic wave propagation in full space is used to deduct the path effect [1,2]. When the epicentral distance is far, the geometric attenuation model and inelastic attenuation model with different epicentral distances should be considered, such as the piecewise attenuation model [3,4]. In addition, if the seismic station is not erected on bedrock rock, the seismic station response should be considered [5]. the distribution of earthquakes is concentrated and the distance between earthquakes is far less than the distance between earthquakes and stations. It is considered that the propagation paths of earthquakes are the same and small-magnitude earthquakes can be used as the empirical Green’s function to deduct the propagation path effect.

1.1. Source Parameter Characteristics’ Study

The study of source parameter characteristics of medium and small earthquakes has always been an important content in seismology research, including earthquake self-similarity, source scaling characteristics, etc. [6,7]. Whether the earthquake self-similarity exists universally is still one of the most controversial scientific issues. Many seismologists have studied the relationship among earthquake source parameters (i.e., calibration relationship) to reveal whether the rupture mechanism of large earthquakes and small earthquakes contains the same physical process, and whether the energy radiation efficiency of earth earthquakes is the same as that of small earthquakes. Some scholars believe that earthquakes are self-similar ([8,9,10,11,12,13]) and the results of other scholars support that earthquakes are not self-similar [14,15,16,17,18,19]. On the other hand, the stress drop and apparent stress in the source parameters are related to the stress level of the seismogenic environment physically; thus, the study of their spatial and temporal distribution characteristics is one of the means to explore the precursory information of strong earthquakes, to calculate the apparent stress of the global shallow earthquakes with a magnitude greater than 5.8 from 1986 to 1991 by using the far-field P-wave, and to discuss the distribution of the global apparent stress and the magnitude of the apparent stress of different types of earthquakes in different tectonic environments. Prieto et al. [20] systematically calculated the source parameters of a large number of small and medium earthquakes in California, the temporal and spatial variation process of stress drop, and the correlation between this change and strong earthquakes. Allmann and Shearer [21] found that before the M_L6.0 earthquake in December 2004, the stress drop in the focal area was significantly higher than that in other areas on the fault, while after the M_L6.0 earthquake, the stress drop in the focal area showed a significant decline. Hardebeck and Aron [22] discussed the relationship between the distribution of high-stress drop and the locking section of fault and rock strength, and pointed out that the distribution area of high-stress drop on the fault represents that the medium here is stronger or bears higher external shear stress. That is to say, the area where the high-stress drop is concentrated may be the potential source nucleation area of moderate and strong earthquakes. Liu et al. [23] and Stankova-Pursley et al. [24] calculated the apparent stress of small- and medium-sized earthquakes in a subduction zone and analyzed the regions where the plate coupling is strengthened and weakened.

1.2. Recent Research Status

In recent years, with the improvement of reservoir digital seismic networks, abundant digital seismic waveform data have been used to study the source parameter characteristics of reservoir earthquakes. Hua et al. [25,26] found that the relationship between stress drop and earthquake size is consistent with the results of Nuttli’s ISD model [27]. The apparent stress of earthquakes increases with the increase of magnitude, which means that large earthquakes are more efficient radiation bodies of seismic energy than small earthquakes. Compared with earthquakes, the stress drop of the two reservoirs is obviously lower, about 10 times smaller. It is considered that this may be due to the increase of pore pressure of underground media caused by reservoir impoundment, which leads to a reservoir-induced earthquake under a relatively low tectonic stress. Qiao and Zhang [28] found that the stress drop of the earthquake source parameters near the Pubugou reservoir area was significantly lower than that of Zipingpu reservoir, while the corner frequency of blasting events in the reservoir area was significantly lower than that of natural earthquake source parameters.

Since the study range’s observation data are very new and have not been open published, the seismic spectral parameters fitting analysis of the area base on the new data have not been further studied. This article focuses on this aspect, which will have to be of practical significance.

2. Analysis Range in Reservoir Area of Lifeline

The study range comprises the earthquake network stations in river areas in Western China, as shown in Figure 1; the selected area includes several reservoir inundation areas and their adjacent areas. In the study area, according to the reservoir water storage area, we selected more than 500 earthquakes (M_L ≥ 3.0) within 400 km from the reservoir water area from January 2007 to December 2015 to research. The seismic network of reservoirs in the lower reaches of the river is distributed along the river. With the progress of network construction in different periods, we use these waveform records to calculate and analyze the source parameters of small- and medium-sized earthquakes in the study area.

In the selection of calculation method, considering that the stations are densely distributed along several rivers and the earthquakes are evenly distributed on both sides of the reservoir, and the seismic records are all near earthquake waveforms, the inelastic depletion of the medium can be ignored when the propagation path effect is removed; and only the elastic attenuation can be deducted by using the theoretical geometric attenuation model R-1.

3. Analysis Method and Result of Seismic Spectrum Parameters

The seismic stations of the two reservoirs are both fixed stations. According to the requirements of site selection in the seismic design code of the people’s Republic of China (GB 50011-2018, part specification for the construction of seismic stations: seismic stations), the sites are selected on bedrock rock mass. Therefore, it is considered that the site response of the two reservoir network stations should be near type I. The specific calculation method process is as follows: (1) the seismic waveform is transformed into the frequency domain by fast Fourier transform; and (2) the integration is performed in the frequency domain; that is, the power spectral integration S_V and S_D of velocity and ground motion displacement are obtained.

$$S_{\mathrm{D}}=2{{\displaystyle \int \left[D\left(f\right)\right]}}^{2}df$$

(1)

$$S_{\mathrm{V}}=2{{\displaystyle \int \left[V\left(f\right)\right]}}^{2}df$$

(2)

Therefore, the limit values of corner frequency and zero frequency of the two parameters of the seismic spectrum are as follows:

$$f_0=\frac{1}{2\pi }\sqrt{\frac{S_{\mathrm{V}}}{S_{\mathrm{D}}}}$$

(3)

$${\Omega }_0^2=4S_{\mathrm{D}}\sqrt{\frac{S_{\mathrm{D}}}{S_{\mathrm{V}}}}$$

(4)

Accordingly, other source parameters can be obtained from the above two values. Among them, the seismic moment M₀ is expressed:

$$M_0=4\sqrt{\frac{5}{2}}\pi \rho {\beta }^3{\Omega }_0$$

(5)

The average of three components of a station is taken as the seismic moment M₀ of the station.

$$M_0=\sqrt{M_{0\mathrm{Z}}^2+M_{0\mathrm{N}\mathrm{S}}^2+M_{0\mathrm{E}\mathrm{W}}^2}$$

(6)

where M_0Z, M_0NS, M_0EW are the square of seismic moment M₀ in the vertical direction, north-south direction, and east-west direction, and M₀ is the radius γ_a of the root rupture.

$${\gamma }_a=\frac{2.34\beta }{2\pi f_0}$$

(7)

$$\Delta \sigma =\frac{7M_0}{16{\gamma }_a^3}$$

(8)

In Equation (8), the ∆σ is stress drop.

In Equations (5) and (7), β is S wave velocity taken as 3.5 km/s; ρ is the density of the medium taken as 2.7 × 10³ kg/m³.

$${\sigma }_{\mathrm{app}}=\eta \overline{\sigma }=\mu \frac{E_s}{M_0}$$

(9)

In Equation (9), σ_app is the apparent stress, the shear modulus is μ = 3 × 10¹⁰ N/m³, and the radiation energy of the seismic wave is E_s = 4πρβS_V.

In the specific processing, we use the observation data of seismic velocity above M_L2.0 near the reservoir area recorded by the reservoir network in the lower reaches of the river. In order to minimize the impact of the propagation path on the calculation results and ensure that many stations participate in the calculation, the data of stations with an epicenter distance less than 80 km are selected in the initial stage of the network construction, and the epicenter distance is selected within 50 km in the middle and latter stages of the station data development work. At first, the instrument response parameters of each component of each station are obtained by seed format data. At the same time, the seed format or EVT format waveform data are converted into ASCII code data; and the mean value, trend, and instrument response are deducted, which can be called by the source parameter calculation program. In the process of calculation, we intercept the S wave as the research object. The interception standard is from the time when the S wave arrives until it decays to twice the noise level and then, uses the 4th order Butterworth filter to conduct a 0.5~24 Hz band-pass filter. The velocity spectrum is obtained by FFT on the filtered velocity record. For the complex velocity spectrum, the corresponding displacement spectrum can be obtained by dividing 2πfi. The power spectrum can be obtained by multiplying the velocity and displacement spectra by their conjugate complex numbers. Finally, S_V and S_D in Equations (1) and (2) can be obtained by integrating them. Each source parameter value can be obtained by these two values.

The final calculation result is the average of the results of multiple stations. The average result is not necessarily the actual value of the parameter. We need to test the confidence interval of the result. Since the calculation results of each station are independent of each other, we refer to the method proposed by Prieto et al. [13] and use the concept of jackknife variance [29] to estimate the 95% confidence interval of each parameter.

Assuming that the values of a source parameter calculated by K stations are X₁, X₂, ..., X_K, the probability characteristics can be described by the parameter θ.

$$\theta =\theta \left[X_1, X_2, \dots , X_K\right]$$

(10)

According to jackknife, throw away one of the values and see what happens to the parameter θ. After throwing away a value X_i,

$${\theta }_i=\theta \left[X_1, \dots , X_{i-1},X_{i+1} ,\dots , X_K\right]$$

(11)

For the K^th θ_i estimate value, the variance can be expressed as follows:

$$\mathrm{var}\left\{\theta \right\}={\sigma }^2=\frac{K-1}{K}{{\displaystyle \sum _{i=1}^K\left[{\theta }_i-{\overline{\theta }}_i\right]}}^2$$

(12)

where, $\overline{\theta }=\frac{1}{K}{\displaystyle \sum _{i=1}^K{\theta }_i}$.

Tukey [30] pointed out that $\left(\ln {\beta }_i-\ln {\overline{\beta }}_i\right)/\sigma$ is close to the student’s t distribution with a degree of freedom of K − 1 in a small sample size. Therefore, the confidence interval for both sides $1-\alpha$ is as follows:

$$\Psi e^{-t_{K-1}\left(1-\frac{\alpha }{2}\right)\sigma }<\Psi \le \Psi e^{t_{K-1}\left(1-\frac{\alpha }{2}\right)\sigma }$$

(13)

Among them, α is the bilateral confidence, σ is the square root of variance, t is the T value corresponding to the student t distribution in K − 1 degree of freedom, Ψ is any source parameter, which can be corner frequency f_c, zero frequency limit Ω₀, seismic moment M₀, stress drop Δσ, apparent stress σ_app, etc. The larger the confidence interval range, the greater the difference of source parameter values given by each station due to the influence of the propagation path, site response, source rupture directionality, and other factors. The smaller the confidence interval range, the smaller the difference between stations.

The results of the source parameters of more than 500 seismic events are obtained, including the corner frequency f_c, zero frequency limit Ω₀, seismic moment M₀, stress drop Δσ, apparent stress σ_app, and the 95% confidence interval of each parameter of each earthquake is obtained. According to the statistics of six source parameters (Figure 2), the statistical results of the number of earthquakes with a corner frequency of source spectrum show that the number of earthquakes decreases from 2.6 Hz to both ends and the span ranges from 1.0 Hz to 4.5 Hz. The average value of corner frequency of all the earthquakes is 2.686 Hz and the median value is 2.664 Hz. The calculation results of the zero frequency limit Ω₀ of the source spectrum are statistically analyzed. The results show that the number of earthquakes decreases with the increase of Ω₀. As a spectral parameter directly related to earthquake size, it is similar to the statistical results according to magnitude. The distribution range of Ω₀ ranges from log10^−2.5 to log10⁰, with an average value of 0.083 and a median value of 0.014. The statistical results of the scalar seismic moment show that when the seismic moment is log10^13.3, the number of earthquakes decreases rapidly along the decreasing direction of seismic moment and decreases gradually along the increasing direction of seismic moment. The average and median values are 1.92 × 10¹⁴ and 3.19 × 10¹³, respectively. The statistical results of the earthquake rupture radius show that the rupture scale of earthquakes is basically within 1 km, mainly concentrated in 0.4~0.6 km; the average fracture size is 0.52 km; and the median is 0.49 km. The calculation results of stress drop show that the stress drop of more than 500 earthquakes is statistically close to a normal distribution when log10 is taken; that is, it is mainly distributed near 1, and the number of earthquakes decreases rapidly in the direction of greater than 1 and less than 1, with an average value of 2.75 bar and a median value of 1.17 bar. The statistical results of the apparent stress show that the apparent stress of most earthquakes is between 0.1 bar and 1.0 bar; the statistical average value is 0.52 bar; and the statistical median value is 0.2 bar.

4. Results and Discussions

According to the definition of local earthquake magnitude M_L, the magnitude is determined based on the records of the Wood–Anderson seismometer. The cut-off frequency of the Wood–Anderson seismometer is 1.2 Hz. When the corner frequency of an earthquake is higher than this frequency, its amplitude is directly proportional to the seismic moment. The statistical results in Figure 2 show that the corner frequencies of earthquakes are basically greater than 1.2 Hz. Therefore, the seismic moment in Figure 3 increases with the increase of magnitude. The scale coefficient of fitting results is 0.996, which is close to the theoretical value 1, and the fitting correlation coefficient is 0.899.

Figure 4 shows the fitting relationship between the seismic moment and the corner frequency of hundreds of earthquakes; however, the fitting relationship of the blue dotted reference line is not satisfied. Similarly, the relationship between the seismic apparent stress and the seismic moment in Figure 5 is not constant. The characteristics of the two spectral parameters show that the earthquakes near the study area do not satisfy the assumption of earthquake self-similarity.

From the change of seismicity rate, there is a good correlation between seismicity and water level near the reservoir area. Moreover, the seismicity is obviously affected by the change of reservoir water storage. Therefore, we compare and analyze the focal parameter characteristics of 148 earthquakes before impoundment (Oct., 2007.10~Nov., 2012) and 406 earthquakes after impoundment (Nov., 2012~Sep., 2015) near the study reservoir area.

The blue solid circle in Figure 6 is the specific focal parameter value of each earthquake, which includes the corner frequency F_c, zero frequency limit Ω₀, fracture radius r, stress drop Δσ, apparent stress σ_app, and the gray vertical line segment, being the 95% confidence interval of the source parameter value. The red solid line is the fitting line of the calibration relationship of each source parameter.

The corner frequency F_c decreases with the increase of the earthquake before and after impoundment. Although the F_c value before impoundment is relatively discrete, the linear fitting coefficient before and after impoundment is very close overall. The zero-frequency limit increases with the increase of the earthquake before and after the reservoir impoundment, and after taking the logarithm, it is in the linear fitting coefficient with the magnitude. Its slopes are 0.965, approximately equal to 1; that means an equal proportion change.

As several reservoirs in the area have a high earthquake intensity, there have been seismic activities before water storage, and these seismic activities have no relationship with the water storage of the reservoirs, which can be identified as tectonic earthquakes. Reservoir-induced earthquakes and tectonic earthquakes in the same study area are selected for comparison and analysis. It shows that there is no obvious difference in focal parameters between natural tectonic earthquakes and reservoir-induced earthquakes.

5. Conclusions and Prospects

(1)

Based on long-term experience, and data observation and analysis, the corner frequency f_c, zero frequency limit Ω₀, seismic moment M₀, stress drop ∆_σ, and apparent stress σ_app are certificated to be the effective source parameters to analyse reservoir earthquakes. Then, those parameters, which are chosen to simulate, are reasonable and meaningful to hazard alleviation and prediction.

(2)

According to existing observation data and conditions, based on the 95% confidence interval of those parameters, through a simulation of the sample area, the calculation results fit well with the statistical results, which shows that the approach is correct and feasible.

(3)

In this study, the results obtained are based on the seismic phase observation report obtained by a river reservoir seismic network from September 2007 to July 2015. Now it is already 2022 and the research can be continued in combination with new data. Nevertheless, the observation data are the most intensive typically before and after impoundment during the period from 2007 to 2015.

(4)

The reservoir-induced earthquakes and natural tectonic earthquakes are difficult to identify. According to the long-term observation data, since a large number of new data appear before and after impoundment, we regard them within 40 km as reservoir-induced earthquakes. The identification of the reservoir-induced earthquakes and natural tectonic earthquakes needs to be studied further.