Abstract

To explore the distribution and evolution of the surrounding rock of underground caverns of hydropower stations, the dynamic compression tests for granite sample were carried out and the dynamic constitutive relation suitable to computer programming was established in terms of strain-based static elastic-plastic damage constitutive relation. By inputting five earthquakes of various frequency spectra, the fully nonlinear dynamic analysis was performed for some large underground caverns in the gorge and mountain region in Sichuan Province in China with the explicit difference method. The simulation results show that, as the input seismic waves become strong gradually, the damage zone generally extends deeply from the middle of side walls. The damage zone reaches the maximum value in areas immediately after the seismic acceleration peak arrives and keeps constant thereafter. According to the initial analysis, if the dominant frequency of input wave approaches the natural frequency of underground caverns, the maximum area of damage zone will increase. Those findings may provide practical data for earthquake-resistant design and construction.

1. Introduction

In the past, it was believed that the damage caused by earthquakes was not serious for underground facilities compared to surface structures, and minimal attention was focused on the seismic study of underground facilities until the 1990s when earthquakes occurred in Kobe (Japan, 1995), Chi-Chi (Taiwan, 1999), and Kocaeli (Turkey, 1999), destroying selected underground structures [1]. The collapse of the Daikai subway station was the first collapse of an urban underground structure due to earthquake forces rather than ground instability. During the Ms 8.0 Wenchuan earthquake which occurred on May 12, 2008, eleven tunnels of the Dujiangyan to Wenchuan highway suffered different degrees of damage [2], and the safe operation of the two hydropower stations in Yingxiuwan and Gengda was threatened [3]. In addition, the surrounding rock of underground caverns of Yingxiuwan Hydropower Station in Southern China deformed seriously during the 2008 Wenchuan earthquake, the lining cracked, and the bottom concrete plate was displaced, which paralyzed the hoisting machine service [4]. Hence, it is essential to investigate the seismic dynamic response of the engineering rock mass of underground caverns.

Three methods, namely, theoretical study, experimental study, and numerical modelling, are commonly used for studying this subject. The theoretical study does not easily consider complicated geological cases [57]. The experimental study is sometimes limited by the existing test techniques and scale effects [8, 9]. The limitations of these two methods motivated the rapid development of the numerical modeling which provides an economical, convenient, and efficient approach.

In recent years, many researchers have investigated the damage mechanism and influence factors of tunnels or caverns under seismic loading [1014]. Zhang applied the modified DDA method to study the seismic response of the underground houses of the Dagangshan Hydropower Station in Western China, and valuable results were obtained [15]. Cui et al. [1619] analyzed the seismic performance of a rock joint with a modified continuously yielding model by DEM software and also studied the control effect of a large geological discontinuity on the seismic response and stability of underground rock caverns. Fu et al. [20, 21] used the discontinuous deformation analysis method to analyze the seismic wave propagation across a rock mass of a large rock cavern complex.

Furthermore, the reasonable seismic input is vital to the dynamic response of underground rock caverns in the mountainous region. Most earthquakes applied at the numerical or analytic models were obtained directly or indirectly from the seismic records measured on or near the ground. However, some seismic data measured at various depths show that the peak acceleration in deep position is generally less than that on the ground. The Japan code stipulates that the acceleration at 20 m depth is 1/2-1/3 of that on the ground. Hu and Xie studied the seismic propagation in Hosokura mining region in Japan and found that the average horizontal peak acceleration at about 400 m depth is approximately 0.65 times that on the ground [22]. In addition, the mountain, especially the mountain peak, can amplify the acceleration. Fortunately, most actual mountainous seismic records were measured at mountain foot, which is little influenced by the amplification effect. In mountainous region, the so called “ground” is assumed to be at the river surface based on the fact that most seismic records were measured at the mountain foot. Thus, for many underground caverns of hydropower stations in mountainous region in Southern China, the numerical model bottom is always hundreds of meters deep from the river surface, and the depth attenuation should be considered in order to simulate the dynamic response under seismic loading in a relatively reasonable way.

In addition, there are abundant microfractures, pores, or cracks within surrounding rock of underground caverns, so the growth of the damage zone under seismic loading should be also considered. Based on the static damage accumulation model provided by Zhu et al. [23], a dynamic damage accumulation model is proposed which is assumed to be proportional to the static one. This paper established a 3D numerical simplified model to analyze the damage zone evolution under seismic loadings of various Fourier spectra.

2. Depth Attenuation in Mountainous Region

The Dagangshan Hydropower Station in mountainous Western China was taken as the case study, and the mountain peak is about 600 m higher than valley bottom. One horizontal seismic record measured in the Panzhihua region in 2008 Wenchuan earthquake was selected as the input earthquake wave for the numerical model, as shown in Figure 1. This seismic record is about 300 km away from the epicentre. Its Fourier spectrum is shown in Figure 2, and the frequencies more than 5 Hz are filtered. The design peak acceleration for the Dagangshan dam site under the exceedance probability of 3% in 50 years is 4.403 m/s2.

To approximately obtain the acceleration attenuation with depth, a simple model with no consideration of caverns was established. The mountain top was averagely flattened to simply simulate the mountain condition. The seismic wave was vertically input at the bottom of the model, with free-field boundary. The peak acceleration of Panzhihua earthquake is assumed to be 1 m/s2 for the purpose of simplification, and the damp ratio is 0.05. The FLAC3d was used to simulate the acceleration propagation, and the results are shown in Figure 3. The results show that the horizontal peak acceleration at the bottom is about 72 percent of that at the valley, which means the dynamic response will be overestimated if the depth attenuation is neglected.

3. Damage Accumulation Model

The damage constitutive model is the foundation for evaluating the safety of rock masses that surround underground caverns during an earthquake. Zhu et al. provided a static elastic-plastic damage constitutive model to describe the deformation and failure of brittle rocks such as granite [23]. The model is expressed aswhere is the initial static elastic modulus, is the unit matrix, and is the static damage third-order tensor. If the damage variables in three principal directions are independent of each other, then the damage evolution equation iswhere is the principal strain, , is the initially cracked static strain, and is the static failure strain.

Assume that the rock strain consists of the elastic strain and crack strain, namely,where is the axial strain, is the lateral strain, is elastic axial strain, is the elastic lateral strain, is the axial strain due to crack propagation, and is the lateral strain due to crack propagation. and can be obtained by Hooker’s law, and the crack strain can be determined by subtracting the crack strain from the elastic strain.

The and curves under the static test of 10−5·s−1 strain ratio for granite drilled from the Dagangshan hydropower site are shown in Figure 4, which indicates that the crack stain mainly comes from the lateral direction. Thus, it is assumed that the only elastic strain occurs along the direction, neglecting the axial brittle damage. In addition, the test results show that the lateral failure strain is not sensitive to the surrounding pressures and loading ratios. Thus, a simple and practical failure criterion was established as follows:

Assume the dynamic elastic-plastic damage constitutive model is similar to the static one, and thus is written aswhere is the initial dynamic elastic modulus.

According to the assumption above, the dynamic constitutive curve is similar to the static one except that their peak stresses and initial elastic modulus are different, as shown in Figure 5. Thus, it is inferred that there is some relationship between their damage variables.

Figure 5 shows the geometric relation asThus, , and then multiplied by , in terms of and . Thus, .

For the purpose of simplification, it is assumed here that the lateral failure strain of bristle rock is not affected by the strain ratio, and the strain amplification coefficients are identical for , and hence

The relationship among rock peak stress, initial elastic modulus, and strain ratio is established by the dynamic test results of Dagangshan granite.

The threshold stain for dynamic damage is assumed to not be affected by the strain ratio, and then it can be determined by the following equation:

The formula above can be used to numerically simulate the damage evolution of Dagangshan underground caverns by using the embedded FISH language of the FLAC3D software.

4. Numerical Model and Boundary Conditions

The three caverns of Dagangshan Hydropower Station all have a gate arch shape, with axis direction NE55. A three-dimensional model is established for three adjacent underground caverns by means of FLAC3D software [24], as shown in Figure 6. The model dimensions are 600 m × 420 m × 200 mm, and the caverns are about 390 to 520 m deep. The top of the model is applied by the corresponding overlying pressure in accordance with the model top depth and rock density. No reinforcement support is considered in this model.

The zone around the underground caverns is mainly composed of relatively fresh and intact granite, with several faults. However, if the faults and fractures are considered in the model, it is very hard to deal with the dynamic simulation due to the much time cost. Thus, the whole rock mass is considered intact rock mass, and the rock property parameters used for numerical simulation are listed in Table 1.

The reflection of seismic waves on the boundaries will heavily affect the dynamic simulation of underground caverns, so it is important to apply proper artificial boundaries around the model. Thus, the free-field boundary is used, and the seismic waves are applied at the model bottom [25]. Considering the depth attenuation, the peak input acceleration is 4.043 × 0.72 = 2.911 m/s2 under the exceedance probability of 3% in 50 years.

5. Input Seismic Wave

The underground cavern is often affected by low-frequency compositions of the earthquake [10], so the five earthquake records which contain rich low-frequency components are adopted as the input waves. In order to compare the dynamic response of underground caverns under seismic waves of various frequency spectra, Panzhihua earthquake record, Kobe earthquake, and three artificial earthquakes are selected for numerical simulation. The Kobe earthquake happened in Japan in 1995. The other three waves are artificially made by triangle series superposition method. The Panzhihua earthquake, Kobe earthquake, and three artificial earthquakes differ greatly in acceleration-time history and Fourier spectrum.

A seismic wave is usually a series of discrete acceleration points, and it should be filtered and rectified before used for numerical simulations [26]. Higher frequency will need smaller mesh dimension, larger grid number, and longer calculation time under the same computation accuracy condition. In this study, the frequencies larger than 8 Hz are removed from the seismic wave, which greatly reduces the number of grids and improves the computation efficiency.

The acceleration-time history of Panzhihua earthquake is shown in Figure 1, and the Fourier spectrum of Panzhihua earthquake is shown in Figure 2. The acceleration-time histories and Fourier spectra of other four earthquakes are shown in Figures 714.

Figures 714 show that the Panzhihua earthquake has the predominant frequencies of 2.3 Hz and 3.2 Hz, Kobe earthquake 1.7 Hz, first artificial earthquake 1.4 Hz and 2.9 Hz, second artificial earthquake 2.6 Hz and 4.3 Hz, and third artificial earthquake 3.4 Hz. The five input earthquakes have obviously different Fourier spectra.

6. Results and Discussion

The above five earthquakes were input into the numerical model for simulation separately. For each earthquake, the damage zone evolution of the surrounding rock of underground caverns was carried out under three typical moments which were selected based on the great change in acceleration.

6.1. Panzhihua Earthquake

The distributions of the damage zone of the surrounding rock under Panzhihua earthquake at the 1.522 s moment, 3.862 s moment, and 10.31 s moment are shown in Figures 1517, respectively. Figure 16 shows the extension of damage zones of powerhouse, transformer chamber, and tailrace surge chamber with time. At the 1.522 s moment, the damage zone of powerhouse is 778.1 m2 in area, that of transformer chamber 235.6 m2, and that of tailrace surge chamber 770.9 m2, totaling about 1785 m2. At the 3.862 s moment, the damage zone of powerhouse is 959.9 m2 in area, that of transformer chamber 355.9 m2, and that of tailrace surge chamber 913.9 m2, totaling about 2229 m2. At the 10.31 s moment, the damage zone of powerhouse is 1010 m2 in area, that of transformer chamber 380.1 m2, and that of tailrace surge chamber 1032.9 m2, totaling about 2422 m2. After the 10.31 s moment, the damage zones of the surrounding rock around three caverns remain constant. The evolution of damage zones around three caverns with time is shown in Figure 18.

Figures 1518 indicate that the damage zone near the middle of side walls of powerhouse and tailrace surge chamber slightly extends inward as the input acceleration increases. When the acceleration approaches the peak, the damage zone extends greatly. The damage zone is generally located near the tops and sides of the tailrace surge chamber and powerhouse as well as the bottom and sides of the transformer chamber. The damage zone area reaches the maximum when the peak acceleration arrives, and then remains constant.

6.2. Kobe Earthquake

The distributions of damage zone of the surrounding rock under Kobe earthquake at the 4.389 s moment, 4.782 s moment, and 7.115 s moment are shown in Figures 1921, respectively. The evolution of damage zones around three caverns with time is shown in Figure 22. It is shown in Figures 1921 that the damage zone near the middle of side walls of powerhouse and tailrace surge chamber deeply extends inward as the input acceleration increases. The damage zone is generally located near the tops and sides of the tailrace surge chamber and powerhouse as well as the bottom and sides of the transformer chamber. At the 4.782 s moment, a continuous damage zone occurs in the region between the transformer chamber and tailrace surge chamber, which does not happen under the action of other earthquakes.

6.3. Three Artificial Earthquakes

The evolution of damage zones around three caverns with time under the first artificial earthquake is shown in Figure 23. The damage zone near the middle of side walls of powerhouse and the tailrace surge chamber strongly extends inward as the input acceleration increases. However, the extension area under the first artificial earthquake is less than that under the Kobe earthquake. The damage zone is generally located near the tops and sides of the tailrace surge chamber and powerhouse as well as the bottom and sides of the transformer chamber.

The evolution of damage zones around three caverns with time under the second artificial earthquake is shown in Figure 24. The total damage zone areas of surrounding rock around three caverns at 1.979 s, 6.35 s, and 7.2 s moments are 1750 m2, 1799 m2, and 1805 m2, respectively. The damage zone extends little under the second artificial earthquake.

The evolution of damage zones around three caverns with time under the third artificial earthquake is shown in Figure 25. The damage zone gradually extends inward as the input acceleration increases. The damage zone area reaches the maximum immediately after the peak acceleration arrives and then keeps constant. Generally speaking, the damage zone extends a little under the third artificial earthquake. The damage zone is generally located near the tops and sides of the tailrace surge chamber and powerhouse as well as the bottom and sides of the transformer chamber.

6.4. Discussion

Figure 26 shows the maximum damage zone areas under various earthquakes and no earthquake, which indicates that the maximum damage zone areas vary greatly with Fourier spectrum in the case of the same peak acceleration. Among the five earthquakes, Kobe earthquake caused the greatest damage zone with area of 4378 m2, while the second artificial earthquake resulted in the least damage zone with area of 1805 m2.

The natural frequency of underground caverns is solved to be approximately 1.9 Hz. Figures 714 show that the Panzhihua earthquake has the predominant frequency of 2.3 Hz and 3.2 Hz, Kobe earthquake 1.7 Hz, first artificial earthquake 1.4 Hz and 2.9 Hz, second artificial earthquake 2.6 Hz and 4.3 Hz, and third artificial earthquake 3.4 Hz. Therefore, the maximum damage zone area will increase when the predominant frequency of underground caverns approaches the natural frequency of underground caverns and vice versa.

7. Conclusions

In this paper, taking the Dagangshan Hydropower Station in China as a study case, the dynamic response of underground caverns was numerically simulated under various Fourier spectra without consideration of reinforced supports such as anchor. The following results are drawn:(1)Even if the acceleration attenuation effect with depth was considered, a continuous damage zone occurs in the region between the transformer chamber and tailrace surge chamber under the Kobe earthquake, which may result in severe failure of caverns.(2)The damage zone near the middle of side walls of the powerhouse and tailrace surge chamber gradually extends inward as the input acceleration increases. The damage zone area reaches the maximum when the peak acceleration arrives and then keeps constant.(3)The maximum damage zone area varies greatly with Fourier spectrum in the case of the same peak acceleration. The maximum damage zone area will increase when the predominant frequency of underground caverns approaches the natural frequency of underground caverns and vice versa.

Data Availability

The data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest

The authors declare that they have no conflicts of interest.

Acknowledgments

The research was supported by the National Basic Research Program of China (no. 2015CB057905), the National Key R&D Program of China (no. 2016YFC0401803), the National Natural Science Foundation of China (nos. 51409263, 11472292, 51709113, and U1704243), the Science and Technology Project of Henan Province (no. 172102310475), and the Open Research Fund of State Key Laboratory of Geomechanics and Geotechnical Engineering, Institute of Rock and Soil Mechanics, Chinese Academy of Sciences (no. Z015011).