Determination of Safety Monitoring Indices for Roller-Compacted Concrete Dams Considering Seepage–Stress Coupling Effects

Abstract

Analyzing the working conditions of a dam using safety monitoring indices (SMIs) is a relatively intuitive and effective method for dam safety evaluation. Therefore, a reasonable and accurate method for determining the SMIs of a dam is of vital importance for dam safety assessment. However, the current methods for determining the SMIs of dams, especially roller-compacted concrete (RCC) dams, have many shortcomings, such as ignoring the construction process of the dam, the coupling effect among multiple physical fields, etc. In this paper, a novel SMI determination method considering the seepage–stress coupling effects was proposed for RCC dams with the assistance of a constructed seepage and stress coupling model so as to address the deficiency of existing RCC dams in determining SMIs. The coupled mathematical model was developed in COMSOL Multiphysics to establish a finite element analysis model of an RCC gravity dam in Henan Province, China. Moreover, the seepage anisotropy of the RCC construction layers was also considered in the model. Finally, the seepage, stress, and deformation characteristics of the RCC dam were analyzed based on the model, and the seepage and deformation SMIs of the dam were determined and compared with traditional methods. The results show that seepage, stress, and displacement fields are distributed similarly for both coupled and uncoupled models. However, in contrast to the uncoupled model, the hydraulic head contour distribution is more dispersed in the coupled model. Additionally, the stress and displacement simulated by the coupled model increase at different rates, with a more pronounced stress concentration near the dam heel. Comparing the seepage and stress SMIs of RCC dam obtained from different methods, it was found that the indices of dam seepage discharge and crest displacement that are calculated by considering the seepage–stress coupling effect and anisotropic characteristics of RCC construction layers are 34.78% and 31.98% lower than results obtained by ignoring these two effects, respectively. Therefore, it is crucial to consider the seepage–stress coupling effect and the anisotropic characteristics of RCC when determining the SMIs for RCC dams.

1. Introduction

Due to the multiple effects of external environmental factors such as calcium leaching, sulfate attack, and freeze–thaw cycles, roller-compacted concrete (RCC) dams inevitably experience aging problems in the long-term operation process [1,2,3,4,5,6,7,8]. The deterioration of the dam from normal to pathological is a gradual (i.e., performance decreases over time) and abrupt (i.e., sudden change in the serviceability of the dam) process. Under the action of internal and external loads, some potential factors may cause the operational state of the dam to change abruptly, thus threatening the safety of the dam project [9]. Therefore, it is necessary to adopt corresponding methods to predict the operational state of the dam, such as establishing a finite element model that conforms to the engineering reality to simulate the performance evolution of the dam during its whole-life-cycle service process, and then conduct a safety evaluation of the dam [10,11,12,13]. In addition, the service state of the dam can also be evaluated by using the actual engineering monitoring data, but this requires the help of predetermined safety monitoring indices (SMIs). Using SMIs to analyze the working state of dams is a direct way, and the key lies in the determination of dam SMIs. During the long-term operation of dams, their bearing capacity will evolve due to changes in material properties, making the formulation of monitoring indices complex. Additionally, the seepage problem of RCC dams is more complex because of the involvement compaction in its construction process in a series of layers [14]. Therefore, the determination of seepage and deformation SMIs for RCC dams is one of the difficulties in current dam safety monitoring research.

There are many methods for determining the SMIs of dams. For concrete dams, the typical methods for determining the indices of the service state monitoring effects are as follows: the confidence interval method, the limit state method, and the small probability method [15]. The confidence interval method is a method that uses statistical theory or the finite element calculation method to establish a mathematical model between the monitoring effect and the load based on the existing dam monitoring data, and the model is used to calculate the SMIs of the monitoring quantity under various loads [16,17]. The limit state method, according to the different methods of calculating the total effect S and the resistance R of the critical load combination, is divided into the safety factor method, the first-order moment extreme value state method, and the second-order moment extreme value state method [18]. The SMIs calculated by this method are mainly the extreme values of the effect quantity. The typical small probability method combines the monitoring effects produced by the load combinations that are unfavorable for strength and stability and determines the SMIs of dams based on existing observation data. It is a more mature method for determining dam SMIs. Currently, the small probability method has been widely used in determining SMIs such as dam deformation, seepage, tensile strength, and temperature [17,19]. This method qualitatively links the effects produced by the load combinations that are unfavorable for strength and stability and estimates the monitoring indices based on previous measured data [17]. Therefore, it is more reasonable and accurate than empirical values. However, the small probability method only considers dam deformation occurring under the previous unfavorable loads to estimate the extreme value, which lacks the ability to identify abnormal dam behavior under conventional loads [17,20].

With the continuous development of methods for determining SMIs for dams, some emerging and efficient methods have been proposed in recent years [21,22,23,24,25,26,27,28,29,30] which have greatly promoted the development of dam safety monitoring theory. Currently, in the process of determining the indices for seepage and deformation safety monitoring of concrete dams, the methods can be generally divided into two categories: statistical analysis and model analysis. Statistical analysis methods are usually simple and convenient to operate and are therefore widely used in engineering. However, statistical analysis methods only start from the actual monitoring data, without considering the causes and mechanisms of dam state changes, nor do they relate to the grade and type of dams. As a result, the physical probability is not clear enough [17]. The model analysis method is based on a mathematical model established via the finite element method, which uses numerical calculation to simulate the seepage and deformation of dams and combines the measured data to invert the seepage and deformation parameters of dam body and foundation. The inverted parameters are then substituted into the numerical model to determine the SMIs of dams. The model analysis method has clear physical concepts and can realize the simulation of dam seepage and deformation under various adverse conditions, as well as the prediction of the dam state throughout the whole life cycle. By doing so, it effectively solves the problem of short monitoring time and short data series of dam measurement. Therefore, the model analysis method is an effective method to determine the comprehensive SMIs of dams.

Despite the extensive work conducted by scholars on the determination of SMIs for dam seepage and deformation, there are still some unresolved issues. These are as follows: (1) The existing SMIs for dams are mostly aimed at dam deformation, and they mostly use mathematical statistics methods. However, these methods cannot consider the causes of dam state changes, and they do not relate to different dam types and dam grades. (2) During the operation of dams, the seepage and stress of the dam body actually affect each other, but the existing SMIs for seepage and deformation are mostly considered from a single perspective; they do not consider the coupling effect of the two fields. (3) Due to the speciality of the construction technology of RCC dams, the rolling layer surface produced by them has different seepage characteristics from conventional concrete dams. Therefore, it is necessary to consider the influence of the rolling layer surface when determining SMIs for RCC dams. For this purpose, this paper considers the influence of the anisotropy of the seepage characteristics of the dam body and the roller layer surface on the seepage and stress fields of the dam. Based on this, a three-dimensional finite element model of an RCC gravity dam is established, and a seepage–stress fully coupled analysis method is constructed and applied in COMSOL Multiphysics commercial finite element software (version: COMSOL Multiphysics 5.2a; COMSOL Co., Ltd.: Stockholm, Sweden) to analyze the seepage and stress fields of the dam. Furthermore, the seepage and deformation SMIs under coupled and uncoupled conditions were proposed and compared with traditional methods.

2. Methodology

2.1. Governing Equations for Seepage and Stress Fields

For the seepage field, the dam body concrete, and the dam foundation rock of the RCC dam can be regarded as porous medium materials. Then the seepage field control equation can be expressed as [31] follows:

$$\{\begin{array}{l}\nabla \cdot \left[-{\delta }_K{k}_s/\left({\rho }_{w}g\right)\left(\nabla p+{\rho }_{w}g\Delta z\right)\right]={\delta }_Q{Q}_s\\ H=p/\left({\rho }_{w}g\right)+z\end{array},$$

(1)

where δ_K and δ_Q are the flux and source proportionality coefficients, usually taken as 1.0; k_s is the permeability of the roller-compacted concrete or rock, m²; ρ_w is the density of water, kg/m³; g is the gravitational acceleration, m/s²; p is the pore water pressure, Pa; z is the elevation, m; Q_s is the source term of the fluid; H is the total hydraulic head, m; and ∇ is the Laplace operator.

For the stress field, the concrete and rock of the RCC dam body and foundation can be regarded as ideal elastic bodies, and the stress component σ_ij can be expressed as follows:

$${\sigma }_{ij}=2G{\epsilon }_{ij}+\lambda {\delta }_{ij}{\delta }_{kl}{\epsilon }_{kl}-\alpha {\delta }_{ij}p ,$$

(2)

where G is the shear modulus, Pa, G = E/[2(1 + ν)]; E is the elastic modulus; ν is the Poisson’s ratio; λ is the Lame parameter, λ = Eν/[(1 + ν)(1 − 2ν)]; ε_ij is the strain component; α is the Biot coefficient; and δ_ij is the Kronecker delta.

Based on the theory of elasticity, the stress and strain should have the following relationship:

$${\epsilon }_{ij}=\frac{1}{2}\left(u_{i,j}+u_{j,i}\right).$$

(3)

According to the static equilibrium condition, it can be known that

$${\sigma }_{ij,j}+F_i=0.$$

(4)

By substituting Equations (2) and (3) into Equation (4), the modified Navier equilibrium equation with displacement as the basic unknown and containing coupling terms is obtained:

$$G{u}_{i,jj}+\left(G+\lambda \right)u_{j,ji}-\alpha p_i+F_i=0 \left(i=x, y, z\right),$$

(5)

where u_i is the displacement in the i direction; and F_i is the volume force in the i direction.

By combining Equations (1) and (5), the mathematical model for the seepage–stress analysis of the RCC dam by the fully coupled method is obtained. Then, by applying the corresponding boundary conditions and initial conditions, the fully coupled problem of seepage and stress of the RCC dam can be solved. Since COMSOL Multiphysics commercial finite element software is used to model and solve the seepage–stress fully coupled problem of RCC dam in this paper, the coupled partial differential equation group can be expressed as follows:

$$\{\begin{array}{ll}\nabla \cdot \left(D_{\alpha }^F\cdot \partial U^F/\partial t\right)+D_{\alpha }^F\cdot \partial U^F/\partial t+\nabla \cdot {\Gamma }^F=F^F& \mathrm{on} \Omega \\ G^F+{\left(\partial R^F/\partial U^F\right)}^T{\mu }^F=-n\cdot {\Gamma }^F& \mathrm{on} \partial \Omega \\ R^F=0& \mathrm{on} \partial \Omega \end{array},$$

(6)

where

$$\{\begin{array}{l}{U}^F={\left\{u^{F_1}, u^{F_2}, \cdots , u^{F_N}\right\}}^{\prime } \\ {\Gamma }^F={\left\{{\mathsf{\Gamma }}^{F_1}, {\mathsf{\Gamma }}^{F_2}, \cdots , {\mathsf{\Gamma }}^{F_N}\right\}}^{\prime } \\ F^F={\left\{{\mathrm{F}}^{F_1}, {\mathrm{F}}^{F_2}, \cdots , {\mathrm{F}}^{F_M}\right\}}^{\prime } \\ G^F={\left\{{\mathrm{G}}^{F_1}, {\mathrm{G}}^{F_2}, \cdots , {\mathrm{G}}^{F_M}\right\}}^{\prime }\\ R^F={\left\{R^{F_1}, R^{F_2}, \cdots , R^{F_M}\right\}}^{\prime }\\ {\mu }^F={\left\{{\mu }^{F_1}, {\mu }^{F_2}, \cdots , {\mu }^{F_M}\right\}}^{\prime }\end{array}$$

(7)

$$\{\begin{array}{l}{\Gamma }^{F_{lj}}=c^{F_{lkji}}\frac{\partial u^{F_k}}{\partial x_i}-{\alpha }^{F_{lkj}}{u}^{F_k}+{\gamma }^{F_{lj}} \\ F^{F_{lj}}=f^{F_{lki}}\frac{\partial u^{F_k}}{\partial x_i}-{\alpha }^{F_{lk}}{u}^{F_k}\\ G^{F_l}=g^{F_l}-q^{F_{lk}}{u}^{F_k}\\ R^{F_m}=r^{F_m}-h{F}_m{}^{F_{ml}}{\mu }^{F_l}\end{array},$$

(8)

where Ω is the computational domain; ∂Ω is the outer boundary of Ω; n is the outward normal direction of Ω; t is time; N is the number of unknown variables, l = 1~N; M is the number of constraint conditions, m = 1~M; F^Fl, G^Fl and R^Fm are scalars; Γ^Fl is a vector; μ^F1, μ^F2, …, μ^FM are Lagrange multipliers; ${\mathit{D}}_{\alpha }^F$ is a 3rd-order tensor, $d_{\alpha }^{F_{131}}=d_{\alpha }^{F_{142}}={\rho }_w/{\rho }_0$ (the rest are 0); ${\mathit{D}}_a^F$ is a 2nd-order tensor, $d_{\alpha }^{F_{11}}=\varphi \alpha$ (the rest are 0); and ϕ is porosity.

The expressions of each parameter in Equation (8) are as follows:

$$\left[\begin{array}{cc}{c}^{F_{1111}}& c^{F_{1112}}\\ c^{F_{1121}}& c^{F_{1122}}\end{array}\right]=\frac{{\rho }_w}{{\rho }_0}\frac{1}{{\mu }_w}\left[\begin{array}{cc}{k}_{xx}& k_{xy}\\ k_{yx}& k_{yy}\end{array}\right]$$

(9)

$$\begin{array}{l}\left[\begin{array}{cc}{c}^{F_{1331}}& c^{F_{3312}}\\ c^{F_{3321}}& c^{F_{3322}}\end{array}\right]=\left[\begin{array}{cc}{\varsigma }_1& 0\\ 0& {\varsigma }_2\end{array}\right], \left[\begin{array}{cc}{c}^{F_{3411}}& c^{F_{3412}}\\ c^{F_{3421}}& c^{F_{3422}}\end{array}\right]=\left[\begin{array}{cc}0& {\varsigma }_1\\ {\varsigma }_2& 0\end{array}\right]\\ \left[\begin{array}{cc}{c}^{F_{4311}}& c^{F_{4312}}\\ c^{F_{4321}}& c^{F_{4322}}\end{array}\right]=\left[\begin{array}{cc}0& {\varsigma }_2\\ {\varsigma }_1& 0\end{array}\right], \left[\begin{array}{cc}{c}^{F_{4411}}& c^{F_{4412}}\\ c^{F_{4421}}& c^{F_{4422}}\end{array}\right]=\left[\begin{array}{cc}{\varsigma }_2& 0\\ 0& {\varsigma }_1\end{array}\right]\end{array}$$

(10)

$${\gamma }^{F_{11}}=0, {\gamma }^{F_{12}}=-{\rho }_{w}g, {\gamma }^{F_{2j}}={\gamma }^{F_{3j}}={\gamma }^{F_{4j}}=0 \left(i=1, 2\right)$$

(11)

$$f^{F_1}=Q_s, f^{F_2}=0, f^{F_3}=F_x={\gamma }^{F_{4j}}=0 \left(i=1, 2\right)$$

(12)

$$\left\{\begin{array}{c}{\beta }^{F_{331}}\\ {\beta }^{F_{332}}\end{array}\right\}={\rho }_w{c}_w\frac{k}{{\mu }_w}\left(\nabla p-{\rho }_{w}g\right)$$

(13)

$${\beta }^{F_{311}}={\beta }^{F_{412}}=1-\frac{K_A}{K_B},$$

(14)

where ς₁ = E(1 − ν)/3(1 − 2ν); ς₂ = E/2(1 + ν); K_A is the drainage bulk modulus, Pa; and K_B is the effective bulk modulus of the dam concrete, Pa. Each term in α^Fijk and α^Flk is 0; the rest of the terms in β^Fijk and c^Fijkl are 0.

2.2. Coupling Variable

Once the seepage–stress fully coupled model established, it is further necessary to determine the cross-coupling terms between the two sets of equations. Among them, the relationship between the permeability tensor and the stress state is the key to analyzing the coupling of seepage and stress fields in geotechnical materials. On the one hand, the permeability tensor k_ij of geotechnical materials is affected by the stress component σ_ij, (i.e., k_ij = g(σ_ij)); on the other hand, the stress state σ_ij of geotechnical materials is also affected by the permeability tensor due to the existence of seepage water pressure F_d (i.e., σ_ij = f(k_ij)). Louis [32] established an empirical relationship between the isotropic average permeability coefficient k_f of the rock mass and the effective normal stress σ based on the borehole pressure test results at different depths in a dam site rock mass:

$$k_f=k_0{e}^{-\alpha \prime \sigma },$$

(15)

where k₀ is the initial permeability coefficient of the material, m/s; σ is the effective normal stress, Pa, σ = γH − p, γ is the bulk density of the rock mass, N/m³, H is the thickness of the rock mass, m, p is the pore water pressure, Pa; and α’ is an empirical coefficient, which is related to the integrity of the rock mass material.

For the RCC dam body, its permeability coefficient anisotropy is mainly determined by the tangential and normal permeability coefficients (i.e., k_t and k_n) of the layer surface. Due to the action of the dam self-weight and the upstream water pressure, the hydraulic gap width of the layer surface changes; thus, the tangential permeability coefficient is affected by the stress state and changes. With the help of Formula (15), the material permeability coefficient considering the anisotropy of the layer seepage can be expressed by the following formula:

$$\left[k\right]=k_0\left[\begin{array}{ccc}\exp \left(\lambda {\sigma }_1^{\prime }\right)& & \\ & \exp \left(\lambda {\sigma }_2^{\prime }\right)& \\ & & \exp \left(\lambda {\sigma }_3^{\prime }\right)\end{array}\right],$$

(16)

where λ is the influence coefficient, which can be determined by experiments. In this paper, both concrete and rock layer are taken as 0.35; and σ_i^’ is the effective principal stresses, i = 1, 2, 3.

2.3. Monitoring Indices Determination Method

The most widely used method for determining the service state SMIs of concrete dams is the confidence interval estimation method. This method takes a significant level α, then P_α = α is a small probability event. In the determination of SMIs for RCC dams, by establishing the model between the monitoring effect quantity and the load of the dam service state, and through calculation and analysis, the difference between the measured value E and the effect quantity $\stackrel{⏜}{E}$ generated by each load when acting is obtained. It is assumed that a sample with a sub-sample number of n can be obtained by model calculation [17]:

$$E=\left\{E_{\mathrm{m}1},E_{\mathrm{m}2},\cdots ,E_{\mathrm{m}n}\right\}$$

(17)

$$\overline{E}=\frac{1}{n}{\displaystyle \sum _{i=1}^n{E}_{\mathrm{m}i}}$$

(18)

$${\sigma }_E=\sqrt{\frac{1}{n-1}\left({\displaystyle \sum _{i=1}^n{E}_{\mathrm{m}i}^2-n{\overline{E}}^2}\right)}.$$

(19)

Using the statistical verification method to verify its distribution, the distribution function of the probability density function can thus be determined. The normal distribution function is the most commonly used variable distribution form. Assuming E_m is the extreme value of the environmental effect quantity, when E > E_m, the dam will be in danger, and its probability expression is as follows:

$$P_{\alpha }=P_{\left(E>E_{\mathrm{m}}\right)}={\displaystyle \underset{E_{\mathrm{m}}}{\overset{\infty }{\int }}f\left(E\right)\mathrm{d}E}.$$

(20)

Once the distribution of E_m is obtained, it is necessary to estimate E_m by determining the failure probability α, which is typically determined based on the importance of the dam. Based on the dam’s significance, once the failure probability α is determined, E_m can be directly obtained from the distribution function of E_m:

$$E_{\mathrm{m}}=F^{-1}\left(\overline{E}, {\sigma }_E, \alpha \right).$$

(21)

For the calculation of a normal distribution, let $t=(\delta -\overline{\delta })/\sigma$ (t~N(0, 1)), and the calculation formula for a specific failure probability σ is then transformed into the calculation formula for the standard normal distribution:

$$f\left(t\le t_{\mathrm{m}}\right)=1-\sigma ={\displaystyle \underset{-\infty }{\overset{t_{\mathrm{m}}}{\int }}\frac{1}{\sqrt{2\mathsf{\pi }}}}{e}^{-\frac{1}{2}{t}^2}\mathrm{d}t.$$

(22)

For a specific failure probability, it can be obtained by referring to the standard normal distribution table (as seen in

Table 1. The standard normal distribution for a specific failure probability.

Failure Probability δm/%	Standard Normal Distribution tm
1	2.33
5	1.65

). Based on the retrieved information from the table, by substituting t_m into the calculation, δ_m can be determined.

3. Engineering Example

3.1. Project Overview

A reservoir is located at the outlet of the Beimang River in Henan Province, China. It controls a drainage area of 94 km², accounting for 56% of the total watershed area. The total storage capacity of the reservoir is 10.75 million m³, and it is a medium-sized Grade III water conservancy project that is comprehensively utilized. The development objectives of the reservoir include flood control, improvement of local groundwater environment, urban water supply, irrigation, aquaculture, and tourism. The design criteria for the reservoir are a 50-year flood design and a 500-year flood check. The dam type is an RCC gravity dam, with the outer part of the dam constructed using conventional concrete and the inner part using RCC. The dam consists of non-overflow dam sections, overflow dam sections, flood discharge bottom outlets, and water supply and irrigation dam sections (Figure 1a). The basic profile of the dam is triangular, with a maximum dam height of 77.6 m, a dam crest elevation of 317.60 m, and a dam crest length of 224.5 m. It is divided into 10 dam sections, including a left bank impervious dam section of 63.50 m (No. 1–3 dam sections), a diversion dam section of 15.0 m (No. 4 dam section), an overflow dam section of 68.00 m (No. 5–7 dam sections), and a right bank impervious dam section of 74.0 m (No. 8–10 dam sections). The overflow weir of the dam adopts an open-channel free overflow with no gate control, arranged in five openings, with a net width of 12.0 m per opening. The crest elevation of the overflow weir is 313.0 m, and the type of the overflow weir is a practical WES-type weir. The energy dissipation at the end of the overflow dam is achieved through a trajectory bucket-type energy dissipation.

Based on the geological conditions and layout characteristics of the project, the reservoir safety monitoring system focuses on the dam body, dam foundation, and anti-seepage works, while also considering other hydraulic structures. The dam is equipped with three monitoring sections, namely, Section 2 (Station 0 + 33.50), Section 6 (Station 0 + 104.00), and Section 8 (Station 0 + 150.50). The main monitored indices include dam deformation, seepage, stress–strain, and environmental parameters. The monitoring of the maximum dam height section in the overflow dam section is of particular importance. Therefore, the analysis is conducted using the actual monitoring data from the designated key monitoring section: Section 6 (No. 6 dam section) at Station 0 + 104.00 (refer to Figure 1b). Within the dam body, five piezometers (P1–P5) are installed, all located on the compaction surface. In addition, four piezometers (P6–P9) are positioned at the dam foundation (indicated by red dots in Figure 1b). Furthermore, nine displacement meters are distributed at different locations on the dam foundation, as shown by the blue dots in Figure 1b. For further details on the monitoring data, reference can be made to the literature [33].

3.2. Establishment of Finite Element Model

A finite element model for seepage–stress coupling analysis of an RCC gravity dam was established using the COMSOL Multiphysics commercial software to simulate the seepage, stress, and deformation of the dam body and foundation. Figure 2 shows the finite element model of the main overflow section (i.e., Section 6) of the RCC gravity dam. The calculation coordinate system in Figure 2a is defined as follows: the x-axis is positive in the downstream direction and negative in the upstream direction; the y-axis is positive in the vertical upward direction and negative downward; the z-axis is positive from the left bank to the right bank and negative vice versa. According to the actual situation of the reservoir dam and the distribution of foundation rock layers, the calculation range of 20 m (i.e., z direction) is selected on Section 6; 175 m and 175 m (x direction) are cut from the heel and toe of the dam, respectively, in the upstream and downstream directions; and the vertical direction is taken down to the elevation of 100 m (i.e., y direction). The model is discretized using tetrahedral mesh, and the mesh division result is shown in Figure 2b. After discretization, the dam body and foundation have a total of 27008 nodes and 18535 elements. Among them, a total of 1406 tetrahedral elements were generated in the anti-seepage curtain, as well as 2457 grid nodes.

The boundary conditions of the model should be pre-determined before performing numerical calculations. To better illustrate the setting of the boundary conditions, a schematic of the boundary setting for a 2D RCC dam is shown in Figure 3. For the seepage field, the boundaries below the upstream water level (i.e., AHI) and the downstream channel boundary (i.e., DE) were set as constant head boundaries. The surrounding intercept boundaries (i.e., AB, BC, and CD) as well as the boundary at elevations above the upstream water level (i.e., GH), were set as no-flow boundaries. The boundaries at elevations below the upstream water level and above the downstream channel (i.e., EF and FG) were defined as the permeable layer boundary. The permeable layer boundary is a hybrid boundary, a detailed description of which can be found in the literature [33]. For the stress and deformation calculation, the displacements in the three directions of the model bottom surface (i.e., BC) were constrained, and only the lateral displacements in the x and z directions of the other side surfaces were constrained (i.e., AB and CD). In the model calculation process, the loads mainly consider the seepage pressure and the self-weight.

3.3. Computational Parameters

The accuracy of the constructed seepage–stress coupling model for the RCC dam largely depends on the actual physical and mechanical parameters of the project. For such a large-scale model with monitoring data, the inversion analysis method is a simple and reliable way to determine the engineering parameters. In our previous study, we used the advanced HMPSO–RBFNN model to invert the main parameters (i.e., the permeability coefficient and elastic modulus) involved in the seepage and stress deformation calculations of the project, and we obtained relatively reliable results [33], as shown in

Table 2. Calculation parameters for seepage and stress coupling model.

Material	Permeability Coefficient (m/s)	Density (kg/m3)	Elastic Modulus (GPa)	Poisson’s Ratio
Strongly weathered rock	1.16 × 10−3	2650	20.16	0.20
Weakly weathered rock	1.74 × 10−5	2650	21.22	0.20
Fresh rock mass	8.10 × 10−7	2650	32.81	0.167
Vertical RCC (kz)	3.02 × 10−9	2300	20.74	0.167
Horizontal RCC (kx, ky)	8.57 × 10−8	2300	20.74	0.167
Conventional concrete	2.71 × 10−9	2300	24.98	0.167
Anti-seepage curtain	2.55 × 10−7	2400	20.53	0.20

.

Since the inversion analysis calculation model only considers the data provided by engineering geology, it does not account for potential unknown karst caves and faults present at the actual reservoir site. Therefore, it is assumed that the model and the permeability coefficient obtained through inversion incorporate the “averaged” permeability effects caused by these factors into the permeability coefficient of the “faults and fractures” section, considering their influence on reservoir seepage. In this model, the drainage pipes are treated as homogeneous entities with a permeability coefficient of 1 m/s and mechanical parameters similar to concrete. It should be noted that only the permeability coefficient and elastic modulus are inverted parameters because they have been proven to be the most important parameters affecting the seepage and stress deformation calculations. Density and Poisson’s ratio, as insensitive parameters for the output of seepage and stress deformation calculations, are taken based on engineering experience in this paper. Since the mechanical parameters in this paper were obtained from the literature [33], and the reasonableness of the seepage and deformation parameters were verified by comparing numerical results with the actual monitoring data. Therefore, the validity of the model in this paper is not repeated, and more details about the validation of the model can be found in the literature [33]. It should be noted that the finite element analysis in this paper is performed assuming large displacements.

3.4. Characteristics of Seepage and Stress Fields

Figure 3 shows the distribution of hydraulic head in the dam body and foundation under normal reservoir level conditions during operation, considering and not considering the seepage–stress coupling effect. As can be seen from Figure 4, the distribution pattern of hydraulic head contours for the dam considering and not considering the seepage–stress coupling effect is basically consistent. By comparison, it can be seen that the hydraulic head contours are dense at the concrete face slab and curtain of the RCC dam, indicating that the face slab and curtain have good anti-seepage ability. When considering the seepage–stress coupling effect, the stress is small in the shallow rock and large in the deep rock. Therefore, in the vicinity of the shallow heel of the dam, the hydraulic head contour distribution under coupled and uncoupled conditions is relatively close. In contrast, the seepage–stress coupling effect in the deep rock is obvious, and the difference in hydraulic head contour distribution is significant, and the hydraulic head contours considering the seepage–stress coupling effect are biased towards the upstream.

Figure 5 shows the uplift pressure at the dam foundation obtained by the model under two scenarios: with and without consideration of the seepage–stress coupling. It can be seen from Figure 5 that the uplift pressure calculated by considering the seepage–stress coupling is slightly higher than that without considering the seepage–stress coupling, and the maximum difference between the two calculations can reach 4 m, which occurs near the downstream side of the dam. This is mainly because when the seepage–stress coupling is considered, the fractures on the upstream side of the dam body open due to the effects of upstream water pressure, dam self-weight, and seepage volume force, resulting in an increase in permeability, while the fractures on the downstream rock mass close, resulting in a decrease in permeability, thus leading to an increase in dam uplift pressure.

Figure 6 shows the contour maps of seepage gradient distribution on the concrete face slab of the dam under the coupled and uncoupled scenarios. The results shown in Figure 6 indicate that the seepage gradient distribution pattern on the concrete face slab is basically consistent under the coupled and uncoupled scenarios, and their maximum seepage gradients are both concentrated at the connection between the dam bottom surface and the curtain. When considering the seepage–stress coupling, the maximum seepage gradient on the concrete face slab is 25.3; when ignoring the coupling effect, the maximum seepage gradient on the concrete face slab is 24.1.

The contour maps of seepage gradient in the roller-compacted concrete of the dam under the coupled and uncoupled scenarios are shown in Figure 7. The pore water pressure in the dam body is relatively low, and the pore water pressure gradient is also small, as shown in Figure 6, which is mainly due to the drainage effect of the dam body drainage pipes. The seepage gradient in the dam body is therefore at a low level. The maximum seepage gradient in the dam body is 0.97 when ignoring the coupling, and 0.98 when considering the coupling effect.

The contour maps of the seepage gradient in the grouting curtain under the coupled and uncoupled scenarios are given in Figure 8. It can be seen that the seepage gradient distribution pattern in the curtain is basically consistent under the coupled and uncoupled scenarios, mainly characterized by a gradual decrease in seepage gradient with the decrease in curtain elevation. When considering the seepage–stress coupling effect, the maximum seepage gradient in the curtain is 28.74, while it is 26.38 when ignoring the coupling effect. The seepage gradient in the anti-seepage curtain is larger at its top and lower at its bottom, which is mainly related to the influence of the stratum. The upper part of the foundation rock is strongly weathered, and the permeability of the material is larger, and most of the water head is reduced by the anti-seepage curtain, so the seepage gradient in the curtain in the strong weathering zone is larger. For the anti-seepage curtain in the fresh rock mass area, its permeability is close to that of the rock mass, and a large part of the water head is reduced by the fresh rock mass, so the seepage gradient at the bottom of the curtain is smaller. In addition, the model considering the seepage–stress coupling effect can consider the decrease in rock mass permeability due to rock mass compression, so the seepage gradient in the curtain under the uncoupled scenario is larger than that under the coupled scenario.

Figure 9 and Figure 10 show the stress distribution in the RCC dam under the coupled and uncoupled scenarios, respectively. The results shown in the figures indicate that the overall stress distribution pattern in the dam is basically the same under the coupled and uncoupled scenarios, and the stress distribution conforms to the general rule, i.e., the maximum compressive stress occurs at the heel of the dam, and the maximum tensile stress occurs at the toe of the dam. There are obvious differences in the stress magnitude results of the dam body between when considering and ignoring the coupling effect. When considering the coupling effect, the maximum tensile stress at the heel of the dam is 0.36 MPa; while it is 0.39 MPa when ignoring the coupling effect.

The displacement distribution in the RCC dam under coupled and uncoupled conditions is shown in Figure 11. The vertical displacement distribution pattern of the dam is basically the same under both conditions (Figure 11a,b), with large vertical displacements at the upper part and small vertical displacements at the lower part. The maximum vertical displacement is 15.5 mm when considering the coupling effect, and 22.3 mm when ignoring the coupling effect. The horizontal displacement distribution pattern of the dam is also similar under both conditions (Figure 11c,d). However, the horizontal displacement contour lines at the bottom of the dam body are biased towards the downstream when ignoring the coupling effect. The maximum horizontal displacement is 2.85 mm when considering the coupling effect, and 2.06 mm when ignoring the coupling effect, and it occurs at the top of the dam in both cases.

3.5. Determination of Safety Monitoring Indices

Based on the results of the seepage–stress coupling analysis of the RCC dam in Section 3.4, and considering the monitoring equipment installed in the overflow dam section, this section intends to further propose SMIs for the dam. Since the reservoir project is a medium-sized project, the failure probability α = 1% is taken when using the small probability method to formulate the SMIs, and the following formula is used for calculation:

$$f\left(\delta \le {\delta }_{\mathrm{m}}\right)=1-\mathsf{\sigma }={\displaystyle \underset{-\infty }{\overset{{\delta }_{\mathrm{m}}}{\int }}\frac{1}{\sqrt{2\mathsf{\pi }}}}{e}^{-\frac{1}{2}{t}^2}\mathrm{d}t.$$

(23)

3.5.1. Seepage Safety Monitoring Indices

Seepage is a major factor affecting the safety of RCC dams. In dam safety monitoring, uplift pressure, seepage discharge and around-dam seepage monitoring are essential monitoring items. Among them, for RCC gravity dams, the foundation uplift pressure is an important factor affecting the stability of the dam. The monitoring of the foundation uplift pressure usually uses the uplift pressure coefficient α as its monitoring index.

Table 3. Uplift pressure coefficient under coupled and uncoupled conditions.

Monitoring Point	Water Level under Uncoupled Condition (m)	Water Level under Coupled Condition (m)	Water Level Difference (m)	α under Uncoupled Condition	α under Coupled Condition	δm (α = 1%)
P6	312.86	312.90	0.04	0.17	0.19	0.46
P7	242.45	244.38	1.93
P8	234.15	236.03	1.88
P9	230.00	232.08	0.08

gives the uplift pressure monitoring indices for the overflow dam section under coupled and uncoupled conditions.

As can be seen from Table 3, the foundation uplift pressure considering the seepage–stress coupling effect is higher than that under the uncoupled condition. The uplift pressure coefficient of the foundation under the coupled condition is 0.19, while that under the uncoupled condition is 0.17, both of which are smaller than δ_m. According to the design situation, the uplift pressure of the dam meets the requirements. However, the foundation uplift pressure considering the coupling effect is greater than that which ignores the coupling effect. Therefore, in order to conduct more accurate safety monitoring and safety evaluation of the dam, it is necessary to consider the interaction of seepage and stress.

According to the provisions of the Technical Specification for Safety Monitoring of Concrete Dams in China (DL/T5178-2003) [34], the monitoring of the seepage pressure of the dam body for first- and second-level concrete gravity dams can be carried out as needed, and third-level dams are not monitored.

Table 4. Total head under coupling and uncoupling conditions.

Monitoring Point	Uncoupled Condition (m)	Coupled Condition (m)	Head Difference (m)	δm (α = 1%)
P1	292.15	292.75	0.60	335.42
P2	275.99	280.20	4.21	320.05
P3	274.42	276.00	1.58	305.58

presents the seepage pressure at some measuring points of the concrete face slab and the dam body under coupled and uncoupled conditions.

A comparison of the total water head at each measurement point under coupled and uncoupled conditions in Table 4 shows that the total hydraulic head at each point is higher under the seepage–stress coupled condition than under the uncoupled condition. Among them, the total hydraulic head at the P2 measurement point is 4.21 m higher under the coupled condition than under the uncoupled condition. Therefore, the stress field in the dam body has a significant effect on the pore water pressure of the face slab and the dam body, and the coupling effect between the two fields should be considered when determining SMIs.

Song [35] suggested that the monitoring index of the total seepage discharge of the dam foundation can be established based on the design-specified permeability q₁, the mean effective water head P₁ of the drainage holes, and the total length L₁ of the drainage holes. The safety factor of permeability q₁ is set to 10, and the computed Q₁ is used as the monitoring index of the seepage discharge of the dam foundation:

$$Q_1=0.1q_1{P}_1{L}_1.$$

(24)

For the seepage discharge of the dam body, it can be expressed as follows:

$$Q_2=q_2{P}_{2}H,$$

(25)

where q₂ is the permeability of the dam concrete; P₂ is the average water head on the upstream dam face; and H is the length of the dam.

Table 5. Seepage discharge monitoring indices.

Seepage Discharge	Uncoupled Condition (L/s)	Coupled Condition (L/s)	Empirical Formulae (L/s)	δm (α = 1%)
Dam foundation	15.4	10.7	33.5	35.14
Dam body	2.83	1.19	2.23	3.17

presents the seepage discharges of the dam body and foundation, considering both the coupled and uncoupled models, as well as the results obtained using Formulas (24) and (25). Comparing the results from different methods, it it found that the seepage discharge of the dam body, calculated using the coupled model, is lower than the seepage discharge obtained without considering the coupling, with a difference of approximately 6.34 L/s (i.e., 34.78%). When considering the coupling effect of seepage and stress, the seepage discharge of the dam body is also lower than the monitoring indices derived from empirical formulae and traditional small probability methods, with differences of approximately 23.84 L/s (i.e., 39.00%) and 26.42 L/s (i.e., 68.96%), respectively. Compared to other approaches for determining SMIs, the proposed method for monitoring the dam’s seepage discharge, considering seepage–stress coupling effects, exhibits a higher level of conservatism. Therefore, from a safety monitoring perspective, it is necessary to consider the influence of coupling effects between the two domains when determining the SMI for dam seepage discharge.

3.5.2. Deformation Safety Monitoring Indices

Based on the seepage–stress coupling analysis results in Section 3.4, and considering the distribution of observation points, the deformation monitoring indices for the RCC gravity dams are proposed.

Table 6. Displacement monitoring indices.

Monitoring Point	Coupled Condition (mm)	Uncoupled Condition (mm)	δm (α = 1%)
D1	13.40	19.70	21.73
D2	13.14	19.35	20.87
D3	12.80	18.88	18.35
D4	13.53	19.75	19.73
D5	13.20	19.34	20.82
D6	12.78	18.81	20.20
D7	13.28	19.31	21.77
D8	12.95	18.92	19.42
D9	12.54	18.42	20.81
Top of the dam	21.40	26.51	27.63

shows the dam deformation monitoring indices under the conditions of seepage–stress coupling and uncoupling. As seen from Table 6, under the coupled and uncoupled conditions, the displacements of the monitoring points decrease gradually with the decreasing elevation. When the seepage–stress coupling effects are considered, the displacements of the dam are smaller than those obtained without considering the coupling effect, with a maximum difference of approximately 6.30 mm (i.e., 31.98%). Furthermore, when the coupling effects are considered, the displacements of the dam are also smaller than the deformation monitoring indices derived from the traditional small probability method, with a maximum difference of approximately 8.49 mm (i.e., 39.00%). Compared to alternative approaches for establishing SMIs, the proposed determination of the dam’s deformation monitoring indices, considering the coupling effects, exhibits a higher level of conservatism. From a safety monitoring perspective, it is necessary to consider the deformation data that account for the coupling effects when determining the deformation monitoring indices for dam safety.

4. Conclusions

In this work, a seepage and stress fully coupling model was proposed for RCC dams to overcome the deficiency of existing RCC dams in determining SMIs without considering the coupling effects of the two fields. Based on consideration of the seepage anisotropy of the RCC construction layers, the proposed model is applied to analyze the seepage and stress field characteristics of an RCC gravity dam in Henan Province, China. Furthermore, the SMIs for seepage and deformation of the RCC dam, considering the seepage–stress coupling effects, are determined by the proposed model. The main conclusions are as follows:

(1)

The seepage field distribution pattern of the RCC dam is basically consistent with or without the seepage–stress coupling. However, by considering the seepage–stress coupling, the model can simulate the reduction of the rock permeability due to the pore closure caused by the reservoir water pressure. This phenomenon is reflected in the seepage field as the model simulates more dispersed equipotential lines of water head than those without the coupling effect; the equipotential lines on the upstream side of the anti-seepage wall tend to move upstream; and those on the downstream side tend to move downstream. Moreover, the seepage gradient of the concrete face slab also increases slightly, but the seepage gradient of the dam body and curtain decreases.

(2)

The stress and deformation distribution patterns of the RCC dam are also basically consistent under the conditions of considering and not considering the coupling, with the maximum compressive stress occurring at the heel and the maximum tensile stress occurring at the toe, but there are obvious differences in the stress magnitude. When considering the seepage–stress coupling effect, the maximum tensile stress at the heel is 0.36 MPa, while when not considering the coupling effect, the maximum tensile stress at the heel is 0.39 MPa. When considering the coupling effect, the maximum downstream displacement is 2.85 mm; while, when not considering the coupling effect, the maximum downstream displacement is 2.06 mm. The maximum downstream displacement occurs at the top of the dam in both cases.

(3)

By analyzing the seepage and stress fields of the RCC dam under coupled and uncoupled scenarios, this paper proposed the SMIs of uplift pressure, pore water pressure, seepage discharge, and deformation for RCC dams that account for the seepage–stress coupling effect. The proposed indices were compared with those obtained by ignoring the coupling effect and using empirical formulae and the traditional small probability method, and it was found that the seepage and deformation SMIs considering the coupling effect were more conservative. Therefore, when determining the seepage and deformation SMIs for RCC dams, the seepage–stress coupling effect should be considered.