Your new experience awaits. Try the new design now and help us make it even better

ORIGINAL RESEARCH article

Front. Earth Sci., 07 January 2026

Sec. Georeservoirs

Volume 13 - 2025 | https://doi.org/10.3389/feart.2025.1704200

This article is part of the Research TopicNatural Fracture Characterization for Geoenergy ManagementView all 3 articles

Groundwater regime evolution in an excavation damaged zone of underground water-sealed petroleum storage caverns

Updated
Xianzhong YaXianzhong Ya1Xiaoli Liu,
Xiaoli Liu2,3*Jing ZouJing Zou1Fang WangFang Wang4Sijing Wang,Sijing Wang2,5Xinlei Zhang,Xinlei Zhang2,6
  • 1Sinopec Engineering Incorporation, Beijing, China
  • 2State Key of Laboratory of Hydroscience and Engineering, Tsinghua University, Beijing, China
  • 3College of Water Conservancy and Civil Engineering, Xizang Agriculture and Animal Husbandry University, Linzhi, China
  • 4Architectural Design & Research Institute of Tsinghua University Co., Ltd., Beijng, China
  • 5Institute of Geology and Geophysics, Chinese Academy of Sciences, Beijing, China
  • 6School of Civil Engineering, Liaoning Technical University, Fuxin, China

The excavation of underground engineering structures induces damage in the surrounding rock, resulting in the formation of an excavation damaged zone (EDZ) around the excavation profile. The evolving EDZ will alter the permeability of the surrounding rock mass and thus lead to a changing seepage field. In this article, underground water-sealed petroleum storage caverns are taken as study objects. Based on the measured fractures in the field, a random fracture network simulation technique is used to study the initial permeability characteristics. In addition, the occurrence of the EDZ and its evolution are studied. Combining the initial permeability with the EDZ-induced permeability alteration, and the rainfall process considered, the changing 3D seepage field is analyzed numerically. Water inflow into the cavern working face is predicted, and the influence of excavation on the groundwater table is examined. The results show that the simulation values of the groundwater table agreed well with the monitored ones, considering the effect of cavern excavation. As groundwater conditions are critical to the effectiveness and safety of petroleum storage, investigating the evolution of the groundwater regime during excavation is essential for optimizing construction sequencing and designing reliable water-sealing conditions during operation.

1 Introduction

Energy is experiencing rapid economic development. Petroleum storage is an effective approach to alleviating or avoiding the influence of the oil crisis on the economy, society, and security, through balancing supply and demand and stabilizing oil prices and strategic oil reserves for unexpected emergencies. At present, many countries have built a petroleum strategic reserves system with a large and increasing amount of oil storage. Petroleum storage includes conventional storage on the ground and emerging underground storage. Because of the large amount of storage with very little land usage, underground storage plays an increasingly important role in oil reserves. Because of oil’s lower density than water and its poor intermiscibility with water, underground water-sealed petroleum storage caverns have usually been built at a certain depth below the phreatic line, where the oil can be sealed by underground water. Those deep rock caverns are natural oil reservoirs supported by the rock bearing capacity. Compared with other oil reserve approaches, underground water-sealed petroleum storage caverns have many advantages, such as good safety performance, low surface land use, low investment, long usage period, and little environmental pollution (Zhang et al., 2019; Qian et al., 2025).

The first underground water-sealed petroleum storage caverns were built in Sweden in 1948. From then on, increasing numbers of underground water-sealed petroleum storage caverns have been developed and applied in other parts of the world. At present, more than 300 underground water-sealed petroleum storage facilities have been built in more than 30 countries, which are mainly located in the Scandinavian Peninsula, South Korea, Japan, Germany, France, and Saudi Arabia (Chung et al., 2009; Xue et al., 2015). In China, the underground water-sealed petroleum storage caverns project is still in its infancy and requires further comprehensive development.

The excavation of underground rock caverns is accompanied by the unloading of surrounding stresses and the release of energy. The surrounding rock is unloaded in the direction normal to the cavern surface, and the excavation damaged zone (EDZ) forms gradually, induced by the rebound deformation of the surrounding rocks on the excavation boundary. In the EDZ, the mechanical and hydraulic properties of rocks greatly change as a result of rock damage and disturbance. The EDZ controls the cavern stability. As a result of the external loading and unloading, the EDZ evolves, and the permeability coefficient in the EDZ increases sharply, which greatly endangers the stability of the rock cavern. Consequently, understanding the evolution of the mechanical and hydraulic properties of the EDZ is of great significance in engineering practice, with respect to supporting optimization, safe construction, and stability evaluation.

From the early 1980s, EDZs have received much attention, especially for the safety evaluation of nuclear waste treatment plants. Many theoretical and experimental works have been conducted to classify EDZs. In 1982, an experimental and monitoring study on the EDZ was carried out in the Underground Research Laboratory of Atomic Energy of Canada Limited (Read, 2004). Autio et al. (1998) classified the EDZ into the crushed zone, the fractured zone, and the fissured zone. Bäckblom and Martin (1999) found the EDZ comprised a failed zone, a damaged zone, and a disturbed zone. Sato et al. (2000) conducted field tests about EDZs and classified them as a fracturing zone, a desaturation zone, and a stress redistribution zone. Sheng et al. (2002) found a damage zone, an affected zone, and a slightly affected zone in the EDZ. Martino and Chandler (2004) believed that, after excavation, the surrounding rock consisted of an excavation damaged zone and an excavation disturbed zone. Cai and Kaiser (2005) proposed a quantitative evaluation of the degree and the range of EDZ with a micro-mechanical model, using parameters of crack density and crack distribution. Tsang et al. (2005) summarized some existing findings and proposed the concepts of the excavation influence zone and excavation damaged zone based on water flow and seepage properties. The seepage flow state and its properties changed dramatically in the excavation damaged zone, but not in the excavation influence zone. Based on the mesomechanics, Souley et al. (2001) studied the evolution of the damage variable and the seepage property in the EDZ. Kwon et al. (2006) investigated the engineering characteristics of the EDZ in an engineering project. Li et al. (2022) and Qiao et al. (2022) developed an identification method for dominant seepage channels in fractured rock masses based on water curtain system tests and discrete element modeling and established discriminant criteria. Cao et al. established the multi-parameter nonlinear dynamic seepage equation of fractured rock masses and studied the structural stability of the seepage system in fractured rock masses (Cao et al., 2024; Jia et al., 2024; Cao et al., 2025). In addition, the evolution of the seepage properties of damaged rock has been extensively studied (Salari et al., 2004; Blümling et al., 2007; Levasseur et al., 2013; Deshpande et al., 2025; Sun et al., 2025).

Existing studies have predominantly focused on the qualitative classification and static characterization of the EDZ, with limited attention paid to the dynamic coupling mechanisms between damage evolution and transient groundwater seepage during construction. Moreover, the quantitative relationship between damage variables and permeability evolution in the EDZ has not been systematically established, and coupled numerical models suitable for large-scale storage projects remain underdeveloped. The validation of existing models against real-world engineering data is insufficient, constraining their reliable application in practice. To address these gaps, in this article, three-dimensional unstable seepage evolution with excavation induced damage surrounding underground water-sealed petroleum storage caverns was studied numerically, considering the initial permeability characteristics in the surrounding rock and the rainfall process. The damage constitutive model of an EDZ was also established to investigate its damage evolution during cavern excavation. The method and findings were applied to the large-scale, ongoing Huangdao underground water-sealed petroleum cavern storage project, taking a project of underground water-sealed petroleum cavern storage as a case study. The seepage evolution law and its correlation with damage properties in the EDZ were determined. A comparison between numerical modeling results and field measurements was performed to verify the numerical method. In addition to academic merits in rock mechanics, hydro-mechanical coupling, and underground excavation, the findings are also of engineering significance for real-time adjustments during the construction period and sequence during cavern construction and to effectively design the water-sealed parameters after construction. The objectives of this article are to establish an efficient numerical method to investigate the seepage properties in the EDZ and its correlation with damage evolution during and after the excavation of water-sealed petroleum caverns and to verify its applicability and reliability through comparison with field measurements of an ongoing engineering project.

2 Project profile

As the national oil strategic reserve project, the project of underground water-sealed petroleum caverns storage studied in this article is the first extra-large-scale domestic project at present, which is 300 × 104 m3 of design capacity. The location of cavern storage is superior with convenient traffic. It is located in the northern part of the Tailong region in the south of the Jiaodong Peninsula with hilly landforms. The caverns are in Longque Mountain and trend from the east to the west, with an altitude of 280–350 m. There is an escarpment on the north side of the ridge and the abrupt slope on the south side, which is approximately 35–55°. The gully on the mountain foot extends in the north-south or northeast direction. The main part of the cavern storage is located in the southwest part of Longque Mountain, with an average altitude of 220 m. The highest point is Dadingzi, at an altitude of 350.9 m. The lowest point is at the position of the mouth of the well of the zk012 drill hole, which is at 97.50 m altitude. The surface of the cavern storage location is covered by a large area of forest, including pine trees and locust trees, by certain areas of exposed rock, and by a few areas of farmland near the mouth of the well at the mountain foot.

The cavern storage area is located on a join of two plates, the North China plate and the Yangtze plate, which is the Jiaonan–Weihai orogenic belt with a ductile shear zone and a brittle faulting structure as a main characteristic and a poorly developed fold structure. It is also located on the south rim of the Mouping–Jimo fault zone, trending mainly from the north to the east or nearly from the east to the west. On the west side of the location is the Laojunta faulting zone lying in the northeast direction, and on the east side is the Sunjiagou faulting zone in the same direction, while on the north side is the front Maliangou faulting zone trending from east to west direction. All the faulting zones mentioned above are out of the cavern storage location, which influences the siting of cavern storage. A structural map of the cavern storage is shown in Figure 1.

Figure 1
Map showing geological fault lines in the Xin’an Town area, highlighting the Laojun, Front Miliangou, and Shaoji gou faulting zones. A work area is marked. Green lines depict terrain, and red lines indicate faults. A compass rose is in the top right.

Figure 1. Schematic diagram of the structure of the cave reservoir area.

According to the different geologic ages, rock characteristics, and engineering properties, the stratum lithology is classified into four types: Quaternary remnants diluvial layer, early Cretaceous monzonitic granite, late Proterozoic granite gneiss, and early Cretaceous bright spot dike and diorite dikes. Based on the project overview and the survey data, the weathered surface area of the mountain is classified as Class IV and Class V rock masses according to the degree of weathering, and the excavation area of the cavern is classified as a Class III2 rock mass. The specific parameters are detailed in Table 1.

Table 1
www.frontiersin.org

Table 1. Rock mechanics parameters.

The hydrogeological investigation found that the groundwater at the storage location mainly consists of the pore water in the loose rock mass and the fissure water in the bedrock, including one fissure in the shallow rock strata with mesh-type distribution and one fissure in the deep rock strata with nervation-type distribution.

The underground petroleum storage project includes two construction roadways, nine main caverns, six vertical shafts, and four roadways for the water curtain, which is illustrated in Figure 2. The entrances of two construction roadways are located at the southwest of the storage with the designed altitude of 70 m (in the standard altitude system of Huanghai in 1956, the same below). The roadways continue to the north part of the storage and cross at the north end on both the east and west sides of the storage. There are three branches along the main caverns going to the south of the storage, with a total length of approximately 5819 m and the terminal designed altitude of −30 m. The average slope is approximately 13.3%. The cavern span is 8.3 m, and the height of the cavern is 7.5 m. The nine main caverns are parallel in the direction of 45° west by north. Each of the three caverns is connected to form a whole tank by the construction roadways between them, with three tanks in total. The designed bottom altitude of the main cavern is −50 m, with a designed length of 560–650 m, a designed cavern span of 20 m, and a designed height of 30 m. The total length of nine caverns is 5,537 m. The net distance between the adjacent construction roadways of the main caverns is 25 m, with the net distance between the adjacent main caverns of 30 m. There are two vertical shafts in each of three group-cavern tanks, with diameters of 3 m and 5 m respectively, mouth altitudes of 100 m, and depths of 110 m. At 30 m above the top of the main cavern (10 m altitude), there are four horizontal roadways for a water curtain, oriented vertically to the main cavern direction, with a cavern span of 5 m, a height of 4.5 m, and a total length of 2835 m.

Figure 2
Three-dimensional models showing a terrain and construction process simulation. The top model illustrates a landscape with varying elevations in a range of colors from green at the base to blue at the peaks, indicating height. The lower model includes a wireframe overlay depicting elevation changes, with horizontal brown lines and multicolored lines representing construction elements, set against an axis grid. The terrain model is labeled with height in meters.

Figure 2. Three-dimensional view of the locations and arrangement diagram of the petrol storage caverns.

3 The evolution of EDZ damage in the surrounding rock around the cavern during construction

3.1 The damage constitutive model of the EDZ of the surrounding rock

The excavation work in underground engineering results in stress redistribution, deformation induced by unloading, damage, and failure, which leads to the unrecoverable EDZ. Based on the damage mechanics, Shao et al. (1999) present an expression of a damage tensor as shown in Equation 1:

D=k=1Nmkak3a03a03nknk,(1)

where a0 is the average radius of initial random cracks (penny crack), ak is the extended crack radius, N is the total group number of cracks, mk is the density of the f kth group cracks, and nk is the unit vector of the kth group of cracks.

The deformation energy of EDZ is expressed as a linear equation of the damage variable as follows:

wε,D=gtrεD+λ2trε2+μtrεε+ αtrεtrεD+2βtrεεD,(2)

where λ, μ are the lame parameters of the undamaged material, g is a constant which represents the unrecoverable deformation, and α, β are the correction factors induced by damage in the free energy formula. ε, D are the deformation and damage tensors, respectively. The five parameters related to the constitutive model (λ, μ, g, α, β) can all be determined through standard laboratory tests. The testing methods are detailed in reference (Homand-Etienne et al., 1998; Shao et al., 1999) regarding the determination of these parameters.

The constitutive equation of rock in the damage zone can be derived from Equation 2 as shown in Equation 3:

σ=wε,Dε=gD+λtrεI+2με+ αtrεDI+trεD+2βεD+Dε,(3)

where σ represents the stress tensor.

3.2 The simulation of the damage evolution process in the EDZ of the surrounding rock

The numerical simulation of the EDZ evolution and the associated hydro-mechanical coupling process was conducted using FLAC3D. The model incorporated the damage constitutive model described in Section 3.1. For the boundary conditions, the bottom of the model was fixed, the lateral boundaries were constrained in the normal direction, and the ground surface was free. Hydraulically, the lateral and bottom boundaries were set as no-flow boundaries, while the ground surface was set as a flux boundary to account for rainfall infiltration. Based on this model, the EDZ formation and evolution resulting from the construction of the water-sealed petroleum storage caverns were simulated. The model of caverns is shown in Figure 3. The model boundary is a square with a length of 1,400 m.

Figure 3
Two diagrams show 3D models of grid structures. The first diagram features a red grid embedded in blue terrain contours, depicting spatial placement within a landscape. The second diagram displays a standalone blue grid. Both illustrate structural concepts, likely for architectural or engineering purposes.

Figure 3. Model of the petrol storage caverns.

The modeled EDZ distribution after the cavern excavation is shown in Figure 4. The EDZ distribution is related to the mechanical properties of the surrounding rock and the cavern distribution (the distance between caverns and cavern shape).

Figure 4
3D rendering of a grid-like structure with a blue color gradient, indicating a form of stress or temperature map. Axes labeled X, Y, and Z are shown in the top right corner.

Figure 4. EDZ distribution after the petrol storage caverns were constructed.

As illustrated in Figure 5, the evolution of the EDZ during excavation exhibits distinct spatiotemporal characteristics. When the excavation face advanced to a distance of 10 m from the profile (Figure 5a), initial damage zones emerged at the cavern arch and sidewalls, with a maximum depth of approximately 1.5 m. This is primarily attributed to excavation unloading, which leads to a rapid reduction in radial stress and a concentration of tangential stress around the opening, thereby generating tensile and shear–tensile microcracks in the shallow surrounding rock. After the face advanced 50 m past the profile (Figure 5b), stress redistribution was largely complete, and the EDZ expanded significantly and stabilized. The damage zone depths at the arch and sidewalls increased to 3.2 m and 2.8 m, respectively. This evolutionary process quantitatively demonstrates that EDZ development primarily occurs within the influential zone ahead of the excavation face. Once this zone is passed, the damage range stabilizes, indicating the surrounding rock has reached a new equilibrium under a secondary stress state. This mechanism clarifies that timely support installation within this “disturbed zone” behind the face is most effective in curbing further EDZ development.

Figure 5
Cross-sectional diagrams show EDZ distribution at distances of ten meters and fifty meters from an excavation face. The diagram on the left illustrates the ten-meter distribution, while the one on the right shows the fifty-meter distribution. Both diagrams include contour patterns in colored zones, indicating varying levels of disturbance, with axes labeled for orientation.

Figure 5. EDZ distribution and its evolution during excavation of the petrol storage caverns. (a) EDZ distribution 10 m away from the excavation working face. (b) EDZ distribution 50 m away from the excavation working face.

Figure 6 quantitatively elucidates the crucial role of the support system in controlling the EDZ and maintaining water tightness. Comparing Figure 4 (unsupported) with Figure 6a (systematic rock bolting), it is evident that the rock bolts arched ceiling (φ25@L = 4∼6 m, 1.5 m × 1.5 m, arranged in an interlaced pattern. Side wall: φ25@L = 4∼5 m, 1.5 m × 1.5 m, arranged in an interlaced pattern), through its suspending and combined arch effects, significantly restrained the relaxation of the shallow rock mass. This reduced the average EDZ thickness around the cavern from approximately 3.5 m to 2.5 m, a reduction of 28.6%. As shown in Figure 6b, the application of anchor cables between adjacent caverns effectively controlled the deformation and damage interaction within the rock pillar, reducing the EDZ extent at the pillar center by an additional 40% compared to the unsupported case. This quantitative reduction in damage has significant hydro-mechanical implications: a smaller, less developed EDZ results in a more restricted zone of increased permeability, thereby significantly mitigating the expansion of the groundwater depression cone and establishing a sound rock mass foundation for the effectiveness of the water-sealing system during the operational period.

Figure 6
Diagram showing EDZ distribution with two scenarios. (a) EDZ distribution around a tunnel after rock mass anchoring, with support elements depicted. (b) EDZ distribution between adjacent caverns after anchoring, illustrating stress effects between them.

Figure 6. EDZ distribution after the surrounding rock mass was anchored. (a) EDZ distribution after the surrounding rock mass was anchored. (b) EDZ distribution after adjacent caverns were anchored.

4 The evolution of the seepage parameter in the construction process

4.1 The determination of the initial seepage parameter

Because of the anisotropy property of cracked rock, the seepage process induced by the underground cavern excavation is very complex. In this article, the saturated–unsaturated seepage process of underground water is studied under the influence of excavation and disturbing the cracked rock mass and considering rainfall infiltration.

The seepage governing equation of underground water is as shown in Equation 4.

nρtSρ·kρgμPρgz=Q(4)

Where S is the saturability, k is the permeability, P is the underground water pressure, g is the gravitational acceleration, and Q is the source sink term. In the equation, it is assumed that the water is incompressible, meaning that the density ρ does not change with time t.

The permeability property can be evaluated using a field test or by the geological conditions of the construction location. It can be evaluated quantitatively with the fracture network diagram shown, as shown in Figure 7, using the simulation method of fissure seepage.

Figure 7
Diagram showing groundwater flow paths and fracture networks beneath a rocky surface. Arrows indicate the flow path through fractures, emphasizing two central voids where flow converges. A legend explains the symbols for groundwater flow paths and fracture networks.

Figure 7. Fracture network of rock mass.

Using the methods mentioned above and the data of the underground water level in the drill holes before construction, the initial condition and the boundary condition of seepage flow are inverted, and the seepage field in the construction location is simulated under the rainfall condition. The hydraulic head distribution of the initial seepage field is shown in Figure 8.

Figure 8
Two 3D surface plots with color gradients represent different data sets. The left plot uses colors from blue to red indicating height in meters, ranging from 45 to 135. The right plot displays a similar structure with blue to red colors indicating pressure in kilopascals, ranging from 0 to 500. Both plots include axes labels and legend scales.

Figure 8. Groundwater head distribution and water table contours.

4.2 The rules of seepage parameter evolution of the EDZ in the surrounding rock

Many researchers have experimentally explored the change rules and process of seepage parameters under stress or deformation. All the results show that the relationship between stress and seepage parameter is nonlinear. Based on the deformation analysis of the microcracks, Souley et al. (2001) proposed the formula as follows:

logkk0=Ca3a03lrat3,whenaa0>lrat,k=k0,whenaa0>lrat,(5)

where C is the material constant, k0 is the initial seepage parameter, and lrat is the ratio of the crack length when the seepage can start to occur to the initial crack length. The values of C, k0, and lrat have been determined and calculated through literature methods (Homand-Etienne et al., 1998; Shao et al., 1999).

The change of seepage parameter induced by the crack growth in EDZ can be evaluated quantitatively using Equation 5.

In addition, numerical simulation can be used in researching the change of seepage under the rock deformation, namely, the relationship between the seepage parameter and stress. Figure 9 shows the schematic diagram of the stress–seepage coupling process within the rock mass, and the relationship between seepage and stress is shown as shown in Equation 6:

Kf=Kinf1σ+1KiniKinf,(6)

where Kinf is the seepage parameter of rock mass under infinite stress, and Kini is the initial seepage parameter of the rock mass. The relationship between the seepage parameter and stress is shown in Figure 10.

Figure 9
Illustration shows two diagrams of a material testing process. The top diagram depicts a stressed rectangular block with random blue lines representing internal stress distribution. Arrows labeled with the Greek letter sigma indicate applied stress. The bottom diagram features a similar block placed inside a chamber, compressing the material while arrows show stress direction. A tube connects the chamber to a container.

Figure 9. Schematic diagram of stress–seepage coupling process.

Figure 10
Line graph showing the relationship between stress (in megapascals) and permeability coefficient (in 10^−7 meters per second). The permeability decreases sharply from 10 to around 2 as stress increases from 0 to 5 MPa, then gradually levels off as stress reaches 20 MPa.

Figure 10. Relationship between stress and seepage.

4.3 The seepage field analysis of the EDZ in the surrounding rock

The simulated 3D seepage field in Figure 11 clearly reveals the systematic disturbance of the regional groundwater flow induced by the cavern group excavation. Figures 11a,b show a large-scale groundwater depression cone centered on the storage caverns. The underlying mechanism for this phenomenon is that the excavated underground space acts as a sink for groundwater, while the surrounding EDZ acts as a highly permeable “seepage channel belt,” which greatly facilitates groundwater convergence toward the caverns. The cone is asymmetric; its spatial distribution closely aligns with the orientation of the main caverns, and the groundwater drawdown is relatively smaller above the elevation of the water curtain tunnels. This preliminarily verifies the positive role of the water curtain system in maintaining water pressure above the caverns. This 3D seepage field provides a visual and quantitative basis for assessing water-sealing safety, demonstrating that excavation has fundamentally altered the site’s original hydrogeological conditions.

Figure 11
Panel (a) shows a 3D seepage field distribution with a color gradient from blue to red, indicating varying seepage levels. Contour lines highlight elevation changes. Panel (b) illustrates the 3D seepage field distribution on the profile, with intersecting vertical planes emphasizing field variations. Both include a color bar scale from blue (low) to red (high), and axes for orientation.

Figure 11. Numerical modeling of the 3D seepage field. (a) Three-dimensional seepage field distribution. (b) Three-dimensional seepage field distribution on the profile.

Figure 12 provides a comparative analysis that quantitatively demonstrates the influence of considering rock mass damage on the seepage field assessment. In the ideal scenario without considering damage (Figure 12a), the hydraulic gradient around the caverns is gentle, and the water pressure distribution is relatively high. However, when the permeability enhancement effect of the EDZ is introduced (Figure 12b), the hydraulic gradient around the caverns steepens significantly, and the groundwater pressure generally decreases by 20%–30%, resulting in a more pronounced depression cone. The sectional view in Figure 12c offers a particularly clear comparison: the phreatic level in the model considering damage is significantly lower around the cavern arch and sidewalls than in the model without damage.

Figure 12
Three diagrams illustrate stress distribution and fatigue effects in a material. Diagram (a) and (b) show color-coded stress levels with arrows indicating stress direction. Diagram (c) compares stress lines in scenarios with and without fatigue consideration, highlighting their differences with labels and arrows. Units are in meters.

Figure 12. Seepage field with and without considering damage. (a) Underground water pressure distribution of storage caverns on the profile without considering rock damage. (b) Underground water pressure distribution of storage caverns on the profile, considering rock damage. (c) Underground water pressure distribution of storage caverns on the profile with and without considering rock damage.

5 Dynamic evolution rule of underground water

During the engineering investigation, the underground water levels in 15 boreholes in the cavern’s storage location were measured and found to be 0.18–143.00 m deep with an altitude range of 93.07–268.48 m controlled by the topography here. The groundwater levels in and around the storage location are in the altitude range of 39.00–124.35 m, with the lowest level of 39.00 m, which is the same as that of the Yinjiahe reservoir on the south side of the location.

During the construction process, the water levels in all the boreholes, including the 15 existing boreholes and the seven ones added later, are continually monitored, providing much data for the model evaluation. Because the data from different drill holes follow similar patterns, in this case, only the XZ2# borehole is taken as an example. The evolution rule of underground water under the disturbing influence of cavern excavation is studied over the evolution period of 65 days and is shown in Figure 13. The close agreement between the monitored and simulated groundwater levels not only validates the model’s reliability but also reveals the dynamic mechanism of the groundwater response to excavation disturbance. The initial decline in water level on the 23rd day corresponds to the excavation face advancing into the zone of influence of borehole XZ2#. The subsequent gradual decline (from day 23–55) aligns perfectly with the time-dependent process of EDZ propagation and the associated continuous increase in permeability as the excavation progresses. This demonstrates that the groundwater drawdown is not instantaneous but is controlled by the damage evolution process. The stabilization of the water level after the 55th day indicates that the damage and permeability changes in the surrounding rock in this area had essentially stabilized, and the seepage field reached a new dynamic equilibrium.

Figure 13
A graph showing rainfall and water levels over a 60-day period. The y-axis on the left measures rainfall in millimeters, while the y-axis on the right measures water level elevation. Rainfall is marked with black circles, water level monitoring values with black squares, and calculated water level values with blue diamonds. Rainfall peaks sharply around days 20, 30, 50, and 60. Water levels decrease gradually after an initial high, aligning closely between monitored and calculated values.

Figure 13. A comparison of the calculated values and the monitored values of groundwater level in XZ2# borehole.

6 Conclusion

The evolution law of underground water during cavern excavation is analyzed based on underground water-sealed petroleum cavern storage engineering. The research results show that the underground water is affected by rainfall and cavern excavation. The groundwater level in the borehole begins to decline on the 23rd day and decreases gradually until the 55th day, eventually stabilizing. This trend reflects the significant impact of cavern excavation during this period, with groundwater fluctuation closely correlating with the construction progress. The consistency between simulated results and field monitoring data confirms the validity of the proposed model.

The progression of excavation continuously extends the rock damage zone, modifying seepage properties and affecting the local seepage field. The dynamic evolution of the seepage law in the storage location is closely related to the EDZ distribution. Studying the EDZ distribution regularity under different construction scheme conditions can provide important guidance for the water-sealed cavern storage in adjusting engineering construction and designing water-sealed conditions in the operating period.

Data availability statement

The original contributions presented in the study are included in the article/Supplementary Material; further inquiries can be directed to the corresponding author.

Author contributions

XY: Conceptualization, Writing – original draft, Writing – review and editing. XL: Funding acquisition, Software, Writing – original draft, Writing – review and editing. JZ: Data curation, Investigation, Writing – original draft, Writing – review and editing. FW: Data curation, Writing – original draft, Writing – review and editing. SW: Methodology, Supervision, Writing – original draft, Writing – review and editing. XZ: Methodology, Writing – original draft, Writing – review and editing.

Funding

The author(s) declared that financial support was received for this work and/or its publication. This work was supported by the National Key R&D Program of China (No. 2023YFB4005500).

Conflict of interest

Authors XY and JZ were employed by Sinopec Engineering Incorporation. Author FW was employed by the Architectural Design & Research Institute of Tsinghua University Co., Ltd.

The remaining author(s) declared that this work was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

The author XL declared that they were an editorial board member of Frontiers at the time of submission. This had no impact on the peer review process and the final decision.

Correction note

A correction has been made to this article. Details can be found at: 10.3389/feart.2026.1788688.

Generative AI statement

The author(s) declared that generative AI was not used in the creation of this manuscript.

Any alternative text (alt text) provided alongside figures in this article has been generated by Frontiers with the support of artificial intelligence and reasonable efforts have been made to ensure accuracy, including review by the authors wherever possible. If you identify any issues, please contact us.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

References

Autio, J., Siitari-Kauppi, M., Timonen, J., Hartikainen, K., and Hartikainen, J. (1998). Determination of the porosity, permeability and diffusivity of rock in the excavation-disturbed zone around full-scale deposition holes using the 14C-PMMA and He-gas methods. J. Contam. Hydrology 35 (1-3), 19–29. doi:10.1016/s0169-7722(98)00112-0

CrossRef Full Text | Google Scholar

Bäckblom, G., and Martin, C. D. (1999). Recent experiments in hard rocks to study the excavation response implications for the performance of a nuclear waste geological repository. Tunn. Undergr. Space Technol. 14 (3), 377–394. doi:10.1016/s0886-7798(99)00053-x

CrossRef Full Text | Google Scholar

Blümling, P., Bernier, F., Lebon, P., and Martin, C. D. (2007). The excavation damaged zone in clay formations time-dependent behaviour and influence on performance assessment. Phys. Chem. Earth 32 (8-14), 588–599. doi:10.1016/j.pce.2006.04.034

CrossRef Full Text | Google Scholar

Cai, M., and Kaiser, P. K. (2005). Assessment of excavation damaged zone using a micromechanics model. Tunn. Undergr. Space Technol. 20 (4), 301–310. doi:10.1016/j.tust.2004.12.002

CrossRef Full Text | Google Scholar

Cao, Z. Z., Zhang, S. Y., Li, Z. H., Du, F., Huang, C. H., and Wang, W. Q. (2024). Numerical research on disastrous mechanism of seepage instability of karst collapse column considering variable mass effect. Sci. Rep. 14 (1), 13900. doi:10.1038/s41598-024-63344-w

PubMed Abstract | CrossRef Full Text | Google Scholar

Cao, Z. Z., Xiong, Y., Xue, Y., Du, F., Li, Z. H., Huang, C. H., et al. (2025). Diffusion evolution rules of grouting slurry in mining-induced cracks in overlying strata. Rock Mech. Rock Eng. 58 (6), 6493–6512. doi:10.1007/s00603-025-04445-4

CrossRef Full Text | Google Scholar

Chung, I. M., Kim, T., Lee, K. K., and Han, I. Y. (2009). Efficient method for estimating groundwater head in the vicinity of the underground gas storage caverns in fractured media. J. Hydrologic Eng. 14 (3), 261–270. doi:10.1061/(asce)1084-0699(2009)14:3(261)

CrossRef Full Text | Google Scholar

Deshpande, A. S., Rahatgaonkar, A. M., and Singh, B. N. (2025). Mode-I macro-damage evolution analysis in functionally graded plates using localizing gradient continuum damage model. Mech. Adv. Mater. Struct., 1–21. doi:10.1080/15376494.2025.2486745

CrossRef Full Text | Google Scholar

Homand-Etienne, F., Hoxha, D., and Shao, J. F. (1998). A continuum damage constitutive law for brittle rocks. Comput. Geotechnics 22 (2), 135–151. doi:10.1016/s0266-352x(98)00003-2

CrossRef Full Text | Google Scholar

Jia, Y. L., Cao, Z. Z., Li, Z. H., Du, F., Huang, C. H., Lin, H. X., et al. (2024). Nonlinear evolution characteristics and seepage mechanical model of fluids in broken rock mass based on the bifurcation theory. Sci. Rep. 14 (1), 10982. doi:10.1038/s41598-024-61968-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Kwon, S., Cho, W. J., and Han, P. S. (2006). Concept development of an underground research tunnel for validating the Korean reference HLW disposal system. Tunn. Undergr. Space Technol. 21 (2), 203–217. doi:10.1016/j.tust.2005.06.008

CrossRef Full Text | Google Scholar

Levasseur, S., Collin, F., Charlier, R., and Kondo, D. (2013). On micromechanical damage modeling in geomechanics: influence of numerical integration scheme. J. Comput. Appl. Math. 246, 215–224. doi:10.1016/j.cam.2012.05.022

CrossRef Full Text | Google Scholar

Li, Y. T., Zhang, B., Wang, L., Wu, Y., Wang, H. X., and Peng, Z. H. (2022). Identification of dominant seepage channels in fractured rock masses of underground water-sealed oil storage: a case study. Bull. Eng. Geol. Environ. 81 (9), 357. doi:10.1007/s10064-022-02841-6

CrossRef Full Text | Google Scholar

Martino, J. B., and Chandler, N. A. (2004). Excavation-induced damage studies at the underground research laboratory. Int. J. Rock Mech. Min. Sci. 41 (8), 1413–1426. doi:10.1016/j.ijrmms.2004.09.010

CrossRef Full Text | Google Scholar

Qian, R. P., Liu, X. L., Ma, Q., Feng, G. R., Bai, J. W., Guo, J., et al. (2025). Effect of water intrusion on mechanical behaviors and failure characteristics of backfill body and coal pillar composite specimens under uniaxial compression. J. Clean. Prod. 502, 145388. doi:10.1016/j.jclepro.2025.145388

CrossRef Full Text | Google Scholar

Qiao, L. P., Li, W., Wang, Z. C., Zhong, S. R., and Liu, J. (2022). Assessment of inter-cavern containment property of underground water-sealed oil storage caverns combining discrete fracture network analysis with graph theory. Environ. Earth Sci. 81 (18), 442. doi:10.1007/s12665-022-10576-6

CrossRef Full Text | Google Scholar

Read, R. S. (2004). 20 years of excavation response studies at AECL's underground research Laboratory. Int. J. Rock Mech. Min. Sci. 41 (8), 1251–1275. doi:10.1016/j.ijrmms.2004.09.012

CrossRef Full Text | Google Scholar

Salari, M. R., Saeb, S., Willam, K. J., Patchet, S. J., and Carrasco, R. C. (2004). A coupled elastoplastic damage model for geomaterials. Comput. Methods Appl. Mech. Eng. 193 (27-29), 2625–2643. doi:10.1016/s0045-7825(04)00077-5

CrossRef Full Text | Google Scholar

Sato, T., Kikuchi, T., and Sugihara, K. (2000). In-situ experiments on an excavation disturbed zone induced by mechanical excavation in Neogene sedimentary rock at Tone mine, central Japan. Eng. Geol. 56 (1-2), 97–108. doi:10.1016/s0013-7952(99)00136-2

CrossRef Full Text | Google Scholar

Shao, J. F., Hoxha, D., Bart, M., Homand, F., Duveau, G., Souley, M., et al. (1999). Modelling of induced anisotropic damage in granites. Int. J. Rock Mech. Min. Sci. 36 (8), 1001–1012. doi:10.1016/s1365-1609(99)00070-2

CrossRef Full Text | Google Scholar

Sheng, Q., Yue, Z. Q., Lee, C. F., Tham, L. G., and Zhou, H. (2002). Estimating the excavation disturbed zone in the permanent shiplock slopes of the three Gorges project, China. Int. J. Rock Mech. Min. Sci. 39 (2), 165–184. doi:10.1016/s1365-1609(02)00015-1

CrossRef Full Text | Google Scholar

Souley, M., Homand, F., Pepa, S., and Hoxha, D. (2001). Damage-induced permeability changes in granite: a case example at the URL in Canada. Int. J. Rock Mech. Min. Sci. 38 (2), 297–310. doi:10.1016/s1365-1609(01)00002-8

CrossRef Full Text | Google Scholar

Sun, H., Liu, X. L., Ye, Z. N., and Liu, C. (2025). Mean square flow of fluid in the fractures during rock failures. Geomechanics Geophys. Geo-energy Geo-resources 11 (1), 29. doi:10.1007/s40948-025-00936-4

CrossRef Full Text | Google Scholar

Tsang, C. F., Bernier, F., and Davies, C. (2005). Geohydromechanical processes in the excavation Damaged Zone in crystalline rock, rock salt, and indurated and plastic clays - in the context of radioactive waste disposal. Int. J. Rock Mech. Min. Sci. 42 (1), 109–125. doi:10.1016/j.ijrmms.2004.08.003

CrossRef Full Text | Google Scholar

Xue, Y. G., Li, S. C., Qiu, D. H., Wang, Z. C., Li, Z. Q., Tian, H., et al. (2015). A new evaluation method for site selection of large underground water-sealed petroleum storage depots. Sci. China-Technological Sci. 58 (6), 967–978. doi:10.1007/s11431-015-5825-0

CrossRef Full Text | Google Scholar

Zhang, B., Shi, L., Yu, X., and Qi, S. W. (2019). Assessing the water-sealed safety of an operating underground crude oil storage adjacent to a new similar cavern - a case study in China. Eng. Geol. 249, 257–272. doi:10.1016/j.enggeo.2019.01.008

CrossRef Full Text | Google Scholar

Keywords: evacuation damaged zone, fracture network, groundwater regime evolution, water inflow, water-sealed petroleum storage caverns

Citation: Ya X, Liu X, Zou J, Wang F, Wang S and Zhang X (2026) Groundwater regime evolution in an excavation damaged zone of underground water-sealed petroleum storage caverns. Front. Earth Sci. 13:1704200. doi: 10.3389/feart.2025.1704200

Received: 12 September 2025; Accepted: 08 December 2025;
Published: 07 January 2026; Corrected: 26 January 2026.

Edited by:

John Hooker, University of the Incarnate Word, United States

Reviewed by:

Zhengzheng Cao, Henan Polytechnic University, China
Wenqiang Mu, Northeastern University, China

Copyright © 2026 Ya, Liu, Zou, Wang, Wang and Zhang. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

*Correspondence: Xiaoli Liu, eGlhb2xpLmxpdUB0c2luZ2h1YS5lZHUuY24=

Disclaimer: All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article or claim that may be made by its manufacturer is not guaranteed or endorsed by the publisher.