# Experimental Study on Seepage Anisotropy of a Hexagonal Columnar Jointed Rock Mass

Authors: Yanxin He; Zhende Zhu; Wenbin Lu; Yuan Tian; Xinghua Xie; Sijing Wang
Journal: Shock and Vibration
 (2021)

---

## Abstract

Many columnar jointed rock masses (basalt) are present at the Baihetan hydropower dam site, and their seepage characteristics have a significant impact on the project’s safety and stability. In this study, model samples consisting of material similar to the columnar jointed rock mass with different inclination angles (0°–90°) were prepared and laboratory triaxial seepage tests were performed to study the seepage characteristics of the columnar jointed rock mass under maximum axial principal stress. The experimental results showed that the similar material model samples of columnar jointed rock mass showed obvious seepage anisotropy. The nonlinear seepage characteristics were well described by the Forchheimer and Izbash equations, and the fitting coefficients of the two equations were in good correspondence. The curves describing the relationship between the inherent permeability and the stress of the samples with different dip angles were U-shaped and L-shaped, and a one-variable cubic equation well described the relationship. The 45° angle specimen had the highest sensitivity to the maximum principal stress, and its final permeability increased by 144.25% compared with the initial permeability. The research results can provide theoretical support for the stability evaluation of the Baihetan hydropower station.


---

## Body

## 1. Introduction

A columnar jointed rock mass is formed by rapid cooling and shrinkage of high-temperature lava. It is a discontinuous, multiphase, and anisotropic geological body. These types of rock masses are widely distributed globally. Numerous densely developed slender columnar jointed rock masses are exposed in the dam foundation and underground powerhouse area of the Baihetan hydropower station in China (Figure1) [1, 2].Figure 1
Columnar jointed rock mass in the dam site area of the Baihetan hydropower station.Current research on columnar jointed rock masses has primarily focused on the mechanical properties and anisotropic deformation and failure [3–6]. However, numerous practical engineering projects have shown that the damage and deformation of rock masses and changes in the seepage field caused by engineering disturbance are the dominant reasons for large-scale engineering instability and engineering geological disasters [7, 8]. Statistics have shown that more than 90% of rock slope failures are related to groundwater permeability, 60% of coal mine shaft damage is related to groundwater, and 30%–40% of dam failures of water conservancy and hydropower projects are caused by seepage [9, 10]. Therefore, it is crucial to study the seepage characteristics of columnar jointed rock masses under stress in water conservancy and hydropower projects for practical engineering applications.Research on rock mass seepage has mainly focused on the seepage characteristics of a single joint surface. Louis et al. verified the applicability of the cubic law of laminar flow in a seepage test of smooth parallel plate fractures [11, 12]. However, the structure of a natural joint surface hardly meets the assumption of parallel plate fracture. Therefore, Neuzil et al. put forward a modified cubic law by incorporating different definitions of joint surface roughness [13–16]. However, the variation of the joint surface roughness under stress is highly complex, resulting in a wide range of the joint gap width, complicating the theoretical calculation of the seepage-stress coupling of joints [17–19]. Therefore, it is feasible and effective to conduct indoor seepage-stress coupling tests on rock samples to explore the seepage characteristics under stress. Liu et al. investigated the nonlinear flow characteristics of a fluid at the fracture intersection using an indoor permeability test; they established a fracture network model to simulate the fluid flow state and analyze the influence of the equivalent hydraulic gap width and hydraulic gradient on the nonlinear seepage characteristics [20]; Liu and Li analyzed the relationship between water permeability and sand particle size using a self-designed seepage circuit [21]; Xia et al. studied the nonlinear seepage characteristics of a single joint surface with different roughness values using an indoor shear seepage test [22, 23]; Liu et al. performed conventional indoor triaxial seepage tests on fractured sandstone and analyzed its permeability characteristics [24]; Chao carried out a conventional triaxial cyclic loading and unloading seepage test on a similar material model of a columnar jointed rock mass with inert gas as the measuring medium to determine the relationship between permeability and confining pressure [25]. Few reports were published on the seepage characteristics of multijointed rock masses, especially columnar jointed rock masses.A natural columnar jointed rock mass has large columns. However, indoor tests are limited by the sample size, and the grinding of the natural rock core does not reflect the actual column structure. Therefore, the most effective and reliable method to analyze the permeability characteristics of a columnar jointed rock mass in indoor tests is the use of similar material for the samples. Based on this, according to the natural columnar joint rock mass size and the similarity principle of the model test, this study uses ordinary Portland cement and river sand to fabricate hexagonal prism cylinders to simulate a natural columnar joint rock mass and uses white cement which has the same grade as cementation material to make the similar material model sample of hexagonal prism columnar jointed rock mass. A true triaxial electrohydraulic servo seepage test machine is used to conduct a seepage-stress coupling test using distilled water as the measuring medium. The volume velocity of samples with different dip angles is determined under the maximum principal stressσ1 and different water pressures σs (0.2–0.8 MPa), and the change in permeability and the stress sensitivity of the columnar jointed rock mass is analyzed at different dip angles.


## 2. Sample Preparation and Test Equipment

### 2.1. Sample Preparation

Ordinary Portland cement, river sand, and water with a mass ratio of 1 : 0.5 : 0.4 are used as the column materials of the similar material model to ensure that the physical and mechanical properties of the test sample are similar to that of the natural rock [26]. An early-strength water-reducing agent of 0.4% cement quality was added to ensure the workability and cohesiveness of the mortar and improve the strength of the column in the initial setting. The bonding material between the pillars consisted of a mix of ordinary white Portland cement with the same quality and water with a ratio of 1 : 0.4.A regular hexagonal prism with a side length of 15 mm and a height of 160 mm was created to ensure that the column samples were similar to those of the naturally occurring columnar jointed rock mass at the dam site of the Baihetan hydropower station. A self-designed and manufactured cylinder mold (Figure2(a)) was used to create the sample columns. The mold consisted of several resin plates. The two sides of the mold were attached to long bolts, and the middle was fastened with a clamp to ensure easy assembly and demolding after the column had solidified.Figure 2
The mold for creating the column and the standard curing box.
(a)(b)First, lubricating oil was evenly applied to the surface of the mold. The cement mortar was then poured in layers, and the mold was continuously vibrated to minimize bubble generation and mortar settlement in the column. The mold was then placed in a constant-temperature box, and the mold was removed after 24 h. The demolded column (Figure2(a)) was placed in a standard constant-temperature and constant-humidity curing box (temperature 20 ± 1°C; humidity 95 ± 1.5%, Figure 2(b)) for 28 d. Subsequently, the columns were bonded together with the bonding material and were cured for 28 d in the curing box. The bonded and cured material was cut and ground into 100 × 100 × 100 mm samples of the cubic columnar jointed rock mass with column inclination angles of 0°, 15°, 30°, 45°, 60°, 75°, and 90°, as shown in Figure 3.Figure 3
Similar material samples of the columnar jointed rock masses with different dip angles.


### 2.2. Test Equipment

The test equipment was a true triaxial electrohydraulic servo test machine to determine the rock fracture seepage (Figure4). The test machine had a closed-loop control system with an electrohydraulic servo valve to conduct the permeability test under triaxial loading. Three directions were loaded independently to meet the requirements of the rock permeability characteristics (including seepage deformation and seepage failure).Figure 4
The true triaxial test system for determining the rock fracture seepage.The seepage water weighing system is a crucial part of the equipment and simulates the seepage velocity and water quantity of the rock fracture under different water pressures. The system consisted of a permeable water cylinder, servo valve, pressure gauge and joint, high-precision electronic balance, water collection container, water pipe, and other components. The servo valve controls the pressure of the permeable water cylinder and ensures that the sample has the required osmotic water pressure. The sample was encased by O-rings placed on the upper and lower permeable plates and was wrapped in a latex membrane (Figure5) to ensure that the seepage water could only pass through the cracks of the rock mass and did not flow in the transverse direction. The water seeped into the lower head through the rock fracture and into the water container on the electronic balance through a water pipe. The electronic balance was connected to a computer through a data interface to collect data automatically. The relationship between the rock sample’s triaxial stress and the seepage pressure, seepage velocity, and water quantity was determined.Figure 5
Sample installation.


### 2.3. Test Method and Test Scheme

Typical methods of the single-phase liquid seepage test in the laboratory include the constant water head, variable water head, flow pump, and unstable pulse methods [27]. In this test, we used the constant water head method. A water pressure difference was created at both ends of the rock sample, and the flow velocity and the fluid volume exiting the rock sample were determined.Before the test, the model sample was placed into a vacuum saturation device that used distilled water (Figure6). The samples were wrapped with latex films and loaded into the test machine. During the seepage test, the fluid was assumed to be incompressible; the test was conducted at a constant temperature of 23°C. During the test, the maximum principal stress σ1 was applied in the Z-axis direction, and the intermediate principal stress σ2 and the minimum principal stress σ3 were applied in the transverse direction (Figure 7). The loading scheme was as follows:(1)
The loading mode of the three-dimensional principal stress: force-controlled loading was used. First, the specimen was loaded toFx = Fy = Fz = 1 kN. Then, the stress of σ1 = σ2 = σ3 = 1 MPa was applied in the three directions at a speed of 30 kN/min. σ2 and σ3 remained unchanged at 1 MPa, and σ1 was loaded to 6 MPa with a step size of 1 MPa.(2)
The loading method of the seepage water pressure: after each stage of stress loading was stable, the seepage water pressureσs was applied at 0.2 MPa and was then increased to 0.8 MPa with a step size of 0.1 MPa.Figure 6
The vacuum saturation device.Figure 7
The schematic diagram of stress and hydraulic loading.


## 2.1. Sample Preparation

Ordinary Portland cement, river sand, and water with a mass ratio of 1 : 0.5 : 0.4 are used as the column materials of the similar material model to ensure that the physical and mechanical properties of the test sample are similar to that of the natural rock [26]. An early-strength water-reducing agent of 0.4% cement quality was added to ensure the workability and cohesiveness of the mortar and improve the strength of the column in the initial setting. The bonding material between the pillars consisted of a mix of ordinary white Portland cement with the same quality and water with a ratio of 1 : 0.4.A regular hexagonal prism with a side length of 15 mm and a height of 160 mm was created to ensure that the column samples were similar to those of the naturally occurring columnar jointed rock mass at the dam site of the Baihetan hydropower station. A self-designed and manufactured cylinder mold (Figure2(a)) was used to create the sample columns. The mold consisted of several resin plates. The two sides of the mold were attached to long bolts, and the middle was fastened with a clamp to ensure easy assembly and demolding after the column had solidified.Figure 2
The mold for creating the column and the standard curing box.
(a)(b)First, lubricating oil was evenly applied to the surface of the mold. The cement mortar was then poured in layers, and the mold was continuously vibrated to minimize bubble generation and mortar settlement in the column. The mold was then placed in a constant-temperature box, and the mold was removed after 24 h. The demolded column (Figure2(a)) was placed in a standard constant-temperature and constant-humidity curing box (temperature 20 ± 1°C; humidity 95 ± 1.5%, Figure 2(b)) for 28 d. Subsequently, the columns were bonded together with the bonding material and were cured for 28 d in the curing box. The bonded and cured material was cut and ground into 100 × 100 × 100 mm samples of the cubic columnar jointed rock mass with column inclination angles of 0°, 15°, 30°, 45°, 60°, 75°, and 90°, as shown in Figure 3.Figure 3
Similar material samples of the columnar jointed rock masses with different dip angles.


## 2.2. Test Equipment

The test equipment was a true triaxial electrohydraulic servo test machine to determine the rock fracture seepage (Figure4). The test machine had a closed-loop control system with an electrohydraulic servo valve to conduct the permeability test under triaxial loading. Three directions were loaded independently to meet the requirements of the rock permeability characteristics (including seepage deformation and seepage failure).Figure 4
The true triaxial test system for determining the rock fracture seepage.The seepage water weighing system is a crucial part of the equipment and simulates the seepage velocity and water quantity of the rock fracture under different water pressures. The system consisted of a permeable water cylinder, servo valve, pressure gauge and joint, high-precision electronic balance, water collection container, water pipe, and other components. The servo valve controls the pressure of the permeable water cylinder and ensures that the sample has the required osmotic water pressure. The sample was encased by O-rings placed on the upper and lower permeable plates and was wrapped in a latex membrane (Figure5) to ensure that the seepage water could only pass through the cracks of the rock mass and did not flow in the transverse direction. The water seeped into the lower head through the rock fracture and into the water container on the electronic balance through a water pipe. The electronic balance was connected to a computer through a data interface to collect data automatically. The relationship between the rock sample’s triaxial stress and the seepage pressure, seepage velocity, and water quantity was determined.Figure 5
Sample installation.


## 2.3. Test Method and Test Scheme

Typical methods of the single-phase liquid seepage test in the laboratory include the constant water head, variable water head, flow pump, and unstable pulse methods [27]. In this test, we used the constant water head method. A water pressure difference was created at both ends of the rock sample, and the flow velocity and the fluid volume exiting the rock sample were determined.Before the test, the model sample was placed into a vacuum saturation device that used distilled water (Figure6). The samples were wrapped with latex films and loaded into the test machine. During the seepage test, the fluid was assumed to be incompressible; the test was conducted at a constant temperature of 23°C. During the test, the maximum principal stress σ1 was applied in the Z-axis direction, and the intermediate principal stress σ2 and the minimum principal stress σ3 were applied in the transverse direction (Figure 7). The loading scheme was as follows:(1)
The loading mode of the three-dimensional principal stress: force-controlled loading was used. First, the specimen was loaded toFx = Fy = Fz = 1 kN. Then, the stress of σ1 = σ2 = σ3 = 1 MPa was applied in the three directions at a speed of 30 kN/min. σ2 and σ3 remained unchanged at 1 MPa, and σ1 was loaded to 6 MPa with a step size of 1 MPa.(2)
The loading method of the seepage water pressure: after each stage of stress loading was stable, the seepage water pressureσs was applied at 0.2 MPa and was then increased to 0.8 MPa with a step size of 0.1 MPa.Figure 6
The vacuum saturation device.Figure 7
The schematic diagram of stress and hydraulic loading.


## 3. Test Results and Analysis

### 3.1. The Analysis of the Seepage Characteristics with the Forchheimer Equation

#### 3.1.1. The Relationship between Volume Velocity and Pressure Gradient

Figure8 shows the relationship between the volume velocity Q and the pressure gradient dP/dL of the similar material samples with different inclination angles under maximum principal stress loading.Figure 8
The relationship between the volume velocityQ and the pressure gradient dP/dL for samples with different dip angles. (a) SS0°, (b) SS15°, (c) SS30°, (d) SS45°, (e) SS60°, (f) SS75°, and (g) SS90°.
(a)(b)(c)(d)(e)(f)(g)Significant differences are observed in the relationship between the volume velocity and the pressure gradient for the samples with different dip angles. Atα = 0°–60°, the fitting curves show nonlinear characteristics, and the degree of nonlinearity is different. At the same stress, the flow rate increases with an increase in the pressure gradient, indicating an increase in the permeability of the sample. In contrast, at α = 75°–90°, linear Darcy flow characteristics are observed.For the nonlinear flow relationship in fractures and porous media, Forchheimer proposed a zero-intercept quadratic equation to describe the nonlinear flow process from a macroscopic perspective [28, 29]: (1)−∇p=aQ+bQ2,(2)a=μkA,b=βρA2,where ▽p = dP/dL is the ratio of the pressure difference P between the inlet and outlet and the length L of the sample; Q is the overall volume velocity at the outlet; a and b are the model fitting coefficients, representing the specific gravity of the water pressure drop caused by the linear term and nonlinear term in the seepage test, respectively. When b = 0, the Forchheimer function degenerates to Darcy’s law; μ is the dynamic viscosity coefficient of the fluid; k is the intrinsic permeability of rock; A is the flow area of the rock; β is the nonlinear coefficient; ρ is the density of the fluid.Table1 lists the fitting coefficients a and b, the coefficients of determination R2, the values of the intrinsic permeability k, and the nonlinear coefficient β of the fractured rock obtained from equation (2). The R2 values of all test conditions are greater than 0.98, indicating a good agreement between the experimental and simulated values. Thus, the Forchheimer equation is suitable to describe the relationship between the volume velocity and pressure gradient of similar material samples of the columnar jointed rock mass.Table 1
The fitting results of the nonlinear Forchheimer equation and the calculated values ofk and β
Sample numberMaximum principal stressσ1 (MPa)Fitting coefficienta (Pa·s·m−4 × 1014)Fitting coefficientb (Pa·s2·m−7 × 1022)Coefficient of determinationR2Intrinsic permeabilityk (m2)Nonlinear coefficientβ (m−1)SS0°18.078E − 01−1.361E − 020.99931.579E − 15−6.765E + 1228.821E − 01−1.270E − 020.99851.446E − 15−6.312E + 1239.157E − 01−1.049E − 020.99911.393E − 15−5.214E + 1249.852E − 01−1.543E − 020.99911.294E − 15−7.669E + 1259.792E − 01−2.216E − 020.99891.302E − 15−1.101E + 1368.043E − 01−1.372E − 020.99761.585E − 15−6.819E + 12SS15°17.929E − 01−1.533E − 020.98831.608E − 15−7.619E + 1229.268E − 01−2.202E − 020.98771.376E − 15−1.094E + 1338.912E − 01−2.152E − 020.99331.431E − 15−1.070E + 1348.455E − 01−1.856E − 020.99401.508E − 15−9.225E + 1257.246E − 01−1.476E − 020.99461.760E − 15−7.336E + 1264.928E − 01−6.380E − 030.99352.588E − 15−3.171E + 12SS30°18.718E − 01−1.875E − 020.99821.463E − 15−9.319E + 1221.024E + 00−2.772E − 020.99901.246E − 15−1.378E + 1331.113E + 00−3.456E − 020.99601.146E − 15−1.718E + 1341.027E + 00−3.233E − 020.99631.242E − 15−1.607E + 1358.599E − 01−2.280E − 020.99841.483E − 15−1.133E + 1366.127E − 01−1.117E − 020.99872.081E − 15−5.552E + 12SS45°19.247E − 01−2.320E − 020.99331.379E − 15−1.153E + 1329.551E − 01−2.243E − 020.99521.335E − 15−1.115E + 1338.284E − 01−1.956E − 020.99731.539E − 15−9.722E + 1245.011E − 01−5.140E − 030.99902.545E − 15−2.555E + 1254.406E − 01−3.710E − 030.99962.894E − 15−1.844E + 1263.786E − 01−2.860E − 030.99833.368E − 15−1.421E + 12SS60°18.4222E − 01−1.8600E − 020.99651.514E − 15−9.245E + 1221.0840E + 00−2.3160E − 020.99791.176E − 15−1.151E + 1331.1346E + 00−1.9720E − 020.99611.124E − 15−9.801E + 1241.0895E + 00−2.1190E − 020.99631.170E − 15−1.053E + 1359.9456E − 01−2.4650E − 020.99911.282E − 15−1.225E + 1367.7904E − 01−1.3570E − 020.99921.637E − 15−6.745E + 12SS75°15.3832E − 011.9100E − 030.98672.369E − 159.493E + 1126.5698E − 012.1100E − 030.99461.941E − 151.049E + 1237.7073E − 01−1.0400E − 030.99501.655E − 15−5.169E + 1148.9612E − 01−6.2300E − 030.99161.423E − 15−3.096E + 1251.0872E + 00−9.0500E − 030.99631.173E − 15−4.498E + 1261.4753E + 00−2.6690E − 020.99358.643E − 16−1.327E + 13SS90°15.4456E − 011.6000E − 030.99072.342E − 157.952E + 1127.4407E − 01−6.7000E − 030.99621.714E − 15−3.330E + 1238.9512E − 01−8.2300E − 030.99751.425E − 15−4.091E + 1249.2963E − 01−1.6400E − 030.99751.372E − 15−8.151E + 1151.1876E + 00−2.8810E − 020.97791.074E − 15−1.432E + 1361.3268E + 00−3.7120E − 020.98029.611E − 16−1.845E + 13Table1 indicates that the variation trend of the linear term coefficient a in the fitting equation of the samples with different dip angles with an increase in the maximum principal stress σ1 is different. The trends are shown in Figure 9. It can be seen that at α = 0°–60°, the linear term coefficient a increases first and then decreases, at α = 75°–90°, and the linear term coefficient a increases monotonically.Figure 9
The relationship between the linear coefficienta and the stress for specimens with different dip angles. (a) SS0°, (b) SS15°, (c) SS30°, (d) SS45°, (e) SS60°, (f) SS75°, and (g) SS90°.
(a)(b)(c)(d)(e)(f)(g)Figure10 shows the relationship between the nonlinear coefficient β and the maximum principal stress σ1 for samples with different dip angles. It can be seen that only the nonlinear coefficient β of the samples with dip angles of 75° and 90° is positive when the stress has been initially loaded and the nonlinear coefficient β of other samples is negative. This phenomenon has also appeared in many scholars’ seepage tests; Liao et al. studied and discussed the possibility and causes of negative values of β [27, 30–32]. In this test, the author believes that the negative values of the nonlinear coefficients β reflect the test attributes of the similar material samples of the columnar jointed rock mass. There are two potential reasons: (1) A boundary layer exists on the surface of porous media. As the pressure gradient increases, the thickness of the boundary layer decreases, and the cross-sectional area of the fluid and the permeability increase, resulting in a negative value of β. (2) White cement was used as the filling medium between the columns of the samples. As the pressure gradient increases, the maximum stress the filling medium can bear is exceeded, resulting in the formation of new fractures or the connection of the original fractures in the medium. The pore structure of the medium is changed, increasing the porosity and permeability and resulting in a negative value of β. However, an in-depth investigation is required to reveal the underlying mechanism.Figure 10
The relationship between the nonlinear coefficientβ and the maximum principal stress σ1.


#### 3.1.2. The Relationship between the Intrinsic Permeability and the Maximum Principal Stress

Figure11 shows the relationship between the intrinsic permeability and the maximum principal stress of the samples with different dip angles. The samples of the columnar jointed rock mass show obvious seepage anisotropy, and there are significant differences in the relationship between the natural permeability k and the maximum principal stress σ1 for samples with different dip angles. The relationship can be summarized into two types: In the Type I curves (α = 0°–60°), the intrinsic permeability decreases initially and then increases with an increase in the maximum principal stress, showing a U-shape. Define the stress when the permeability decreases and start to increase as the inflection point stress, and the inflection point stress is 5, 3, 4, 3, and 4 MPa for the samples with different dip angles. In the Type II curves (α = 75°–90°), the intrinsic permeability exhibits a monotonous decrease with an increase in the maximum principal stress, and the curve is L-shaped.Figure 11
The relationship between the intrinsic permeabilityk and the maximum principal stress σ1 for samples with different dip angles. (a) SS0°, (b) SS15°, (c) SS30°, (d) SS45°, (e) SS60°, (f) SS75°, and (g) SS90°.
(a)(b)(c)(d)(e)(f)(g)According to the variation characteristics of the permeability-stress curve, try to fit the test data with one-variable cubic polynomial function:(3)k=a+bσ1+cσ12+dσ13,where k is the intrinsic permeability of the rock; σ1 is the maximum principal stress; a, b, c, and d are the fitting parameters. The R2 values of all fitting curves are greater than 0.96, indicating an excellent goodness-of-fit of the cubic polynomial function describing the relationship between the intrinsic permeability and the maximum principal stress of similar material samples of a columnar jointed rock mass.In addition, the values of the intrinsic permeability of the samples with different dip angles also show obvious anisotropy. The initial permeabilityks, the minimum permeability kmin, the maximum permeability kmax, and the final permeability kf of the samples with different dip angles are plotted in Figure 12.Figure 12
Theks, kmin, kmax, and kf of samples with different dip angles.The maximum decrease in the minimum permeability relative to the samples’ initial permeability is 63.51% for the specimen with a 75° dip angle, and the minimum decrease is 3.18% for the specimen with a 45° dip angle, representing a 20-fold difference. The maximum permeability of the 45° angle specimen increases by 144.25% compared with the initial permeability. The permeability of the 75° and 90° dip angle samples decreases monotonically, and the maximum permeability is the same as the initial permeability. The maximum permeability of the 0° and 60° dip angle samples is basically equal to the initial permeability. The final permeability increases by 0.43%, 60.90%, 42.28%, 144.25%, and 8.11% for the 0°–60° dip angle samples, whereas decreases of 63.51% and 58.96% are observed for the 75° and 90° dip angle samples, respectively (compared with the initial permeability). The sensitivity of the permeability of the samples with different dip angles to the maximum principal stress isα = 45° > α = 75° > α = 15° > α = 90° > α = 30° > α = 60° > α = 0°.


### 3.2. Analysis of the Seepage Characteristics Using the Izbash Equation

The Izbash equation, which is a power function, has also been widely used to characterize the nonlinear seepage characteristics of fluids [33]:(4)−∇p=λQm,where λ and m are empirical coefficients. When m = 1, the Izbash equation degenerates to the linear Darcy’s law. When m = 1-2, equation (4) represents the nonlinear seepage characteristics caused by a significant inertia effect. Values of m = 0-1 are often used to characterize nonlinear seepage in low-permeability rock media caused by a significant solid-liquid interface effect.The Izbash equation is used to fit the test data. The fitting coefficientsλ and m and the R2 values for different maximum principal stress are listed in Table 2. The results show that the Izbash equation is suitable for characterizing the nonlinear seepage characteristics of similar material samples of the columnar jointed rock mass well. The R2 values are all greater than 0.98. The order of magnitude of the fitting coefficient λ is 1013–1014, and the range of the fitting coefficient m is between 0.68 and 1.02.Table 2
The fitting parameters of the Izabsh equation.
Sample numberMaximum principal stressσ1 (MPa)Fitting coefficientλFitting coefficientmCoefficient of determinationR2SS0°19.326E + 130.850.997529.744E + 130.890.999239.798E + 130.920.999541.072E + 140.890.999551.140E + 140.820.998469.529E + 130.840.9993SS15°19.895E + 130.800.989221.142E + 140.790.989831.139E + 140.770.991441.067E + 140.780.994551.006E + 140.750.996667.177E + 130.770.9988SS30°11.043E + 140.810.993021.229E + 140.780.992831.353E + 140.760.988841.395E + 140.690.983351.266E + 140.680.984769.560E + 130.710.9921SS45°11.140E + 140.780.986421.117E + 140.820.990031.100E + 140.750.990846.255E + 130.850.999255.451E + 130.870.999364.901E + 130.860.9993SS60°11.069E + 140.780.995421.202E + 140.860.998331.226E + 140.890.997441.207E + 140.870.998351.188E + 140.800.998969.357E + 130.830.9989SS7515.369E + 131.020.986326.519E + 131.020.994437.757E + 130.990.995049.339E + 130.950.991951.100E + 140.960.995961.496E + 140.930.9928SS90°15.618E + 131.000.990328.247E + 130.910.997839.540E + 130.930.998249.515E + 130.980.997751.364E + 140.830.985761.497E + 140.830.9877The trends of the coefficientsλ with the maximum principal stress σ1 obtained from the Izbash equation for samples with different dip angles are plotted in Figure 13. It can be seen that the trend of the coefficients λ is basically consistent with that of the linear coefficient a obtained from the Forchheimer equation; there is a good corresponding relationship between the two. In addition, there is also a corresponding relationship between the two equations in describing the nonlinear seepage characteristics of similar material samples of columnar jointed rock mass; that is, the fitting coefficient m < 1 of the Izbash equation corresponds to the nonlinear coefficient β < 0 in the Forchheimer equation.Figure 13
The relationship between the fitting coefficientλ of the Izbash equation and the stress for samples with different dip angles. (a) SS0°, (b) SS15°, (c) SS30°, (d) SS45°, (e) SS60°, (f) SS75°, and (g) SS90°.
(a)(b)(c)(d)(e)(f)(g)Figure14 shows the R2 values of the Forchheimer equation and Izbash equations for fitting the test data. Both the zero-intercept quadratic function (Forchheimer equation) and the power function (Izbash equation) well describe the nonlinear flow of the similar material samples of the columnar jointed rock mass, and the difference between them is not significant. However, the Forchheimer equation has a clear theoretical basis and can be simplified by the Navier-Stokes equation if the flow channel does not change significantly. In contrast, the Izbash equation is purely empirical without clear physical meaning.Figure 14
TheR2 values of the Forchheimer and Izbash equations under all working conditions.


## 3.1. The Analysis of the Seepage Characteristics with the Forchheimer Equation

### 3.1.1. The Relationship between Volume Velocity and Pressure Gradient

Figure8 shows the relationship between the volume velocity Q and the pressure gradient dP/dL of the similar material samples with different inclination angles under maximum principal stress loading.Figure 8
The relationship between the volume velocityQ and the pressure gradient dP/dL for samples with different dip angles. (a) SS0°, (b) SS15°, (c) SS30°, (d) SS45°, (e) SS60°, (f) SS75°, and (g) SS90°.
(a)(b)(c)(d)(e)(f)(g)Significant differences are observed in the relationship between the volume velocity and the pressure gradient for the samples with different dip angles. Atα = 0°–60°, the fitting curves show nonlinear characteristics, and the degree of nonlinearity is different. At the same stress, the flow rate increases with an increase in the pressure gradient, indicating an increase in the permeability of the sample. In contrast, at α = 75°–90°, linear Darcy flow characteristics are observed.For the nonlinear flow relationship in fractures and porous media, Forchheimer proposed a zero-intercept quadratic equation to describe the nonlinear flow process from a macroscopic perspective [28, 29]: (1)−∇p=aQ+bQ2,(2)a=μkA,b=βρA2,where ▽p = dP/dL is the ratio of the pressure difference P between the inlet and outlet and the length L of the sample; Q is the overall volume velocity at the outlet; a and b are the model fitting coefficients, representing the specific gravity of the water pressure drop caused by the linear term and nonlinear term in the seepage test, respectively. When b = 0, the Forchheimer function degenerates to Darcy’s law; μ is the dynamic viscosity coefficient of the fluid; k is the intrinsic permeability of rock; A is the flow area of the rock; β is the nonlinear coefficient; ρ is the density of the fluid.Table1 lists the fitting coefficients a and b, the coefficients of determination R2, the values of the intrinsic permeability k, and the nonlinear coefficient β of the fractured rock obtained from equation (2). The R2 values of all test conditions are greater than 0.98, indicating a good agreement between the experimental and simulated values. Thus, the Forchheimer equation is suitable to describe the relationship between the volume velocity and pressure gradient of similar material samples of the columnar jointed rock mass.Table 1
The fitting results of the nonlinear Forchheimer equation and the calculated values ofk and β
Sample numberMaximum principal stressσ1 (MPa)Fitting coefficienta (Pa·s·m−4 × 1014)Fitting coefficientb (Pa·s2·m−7 × 1022)Coefficient of determinationR2Intrinsic permeabilityk (m2)Nonlinear coefficientβ (m−1)SS0°18.078E − 01−1.361E − 020.99931.579E − 15−6.765E + 1228.821E − 01−1.270E − 020.99851.446E − 15−6.312E + 1239.157E − 01−1.049E − 020.99911.393E − 15−5.214E + 1249.852E − 01−1.543E − 020.99911.294E − 15−7.669E + 1259.792E − 01−2.216E − 020.99891.302E − 15−1.101E + 1368.043E − 01−1.372E − 020.99761.585E − 15−6.819E + 12SS15°17.929E − 01−1.533E − 020.98831.608E − 15−7.619E + 1229.268E − 01−2.202E − 020.98771.376E − 15−1.094E + 1338.912E − 01−2.152E − 020.99331.431E − 15−1.070E + 1348.455E − 01−1.856E − 020.99401.508E − 15−9.225E + 1257.246E − 01−1.476E − 020.99461.760E − 15−7.336E + 1264.928E − 01−6.380E − 030.99352.588E − 15−3.171E + 12SS30°18.718E − 01−1.875E − 020.99821.463E − 15−9.319E + 1221.024E + 00−2.772E − 020.99901.246E − 15−1.378E + 1331.113E + 00−3.456E − 020.99601.146E − 15−1.718E + 1341.027E + 00−3.233E − 020.99631.242E − 15−1.607E + 1358.599E − 01−2.280E − 020.99841.483E − 15−1.133E + 1366.127E − 01−1.117E − 020.99872.081E − 15−5.552E + 12SS45°19.247E − 01−2.320E − 020.99331.379E − 15−1.153E + 1329.551E − 01−2.243E − 020.99521.335E − 15−1.115E + 1338.284E − 01−1.956E − 020.99731.539E − 15−9.722E + 1245.011E − 01−5.140E − 030.99902.545E − 15−2.555E + 1254.406E − 01−3.710E − 030.99962.894E − 15−1.844E + 1263.786E − 01−2.860E − 030.99833.368E − 15−1.421E + 12SS60°18.4222E − 01−1.8600E − 020.99651.514E − 15−9.245E + 1221.0840E + 00−2.3160E − 020.99791.176E − 15−1.151E + 1331.1346E + 00−1.9720E − 020.99611.124E − 15−9.801E + 1241.0895E + 00−2.1190E − 020.99631.170E − 15−1.053E + 1359.9456E − 01−2.4650E − 020.99911.282E − 15−1.225E + 1367.7904E − 01−1.3570E − 020.99921.637E − 15−6.745E + 12SS75°15.3832E − 011.9100E − 030.98672.369E − 159.493E + 1126.5698E − 012.1100E − 030.99461.941E − 151.049E + 1237.7073E − 01−1.0400E − 030.99501.655E − 15−5.169E + 1148.9612E − 01−6.2300E − 030.99161.423E − 15−3.096E + 1251.0872E + 00−9.0500E − 030.99631.173E − 15−4.498E + 1261.4753E + 00−2.6690E − 020.99358.643E − 16−1.327E + 13SS90°15.4456E − 011.6000E − 030.99072.342E − 157.952E + 1127.4407E − 01−6.7000E − 030.99621.714E − 15−3.330E + 1238.9512E − 01−8.2300E − 030.99751.425E − 15−4.091E + 1249.2963E − 01−1.6400E − 030.99751.372E − 15−8.151E + 1151.1876E + 00−2.8810E − 020.97791.074E − 15−1.432E + 1361.3268E + 00−3.7120E − 020.98029.611E − 16−1.845E + 13Table1 indicates that the variation trend of the linear term coefficient a in the fitting equation of the samples with different dip angles with an increase in the maximum principal stress σ1 is different. The trends are shown in Figure 9. It can be seen that at α = 0°–60°, the linear term coefficient a increases first and then decreases, at α = 75°–90°, and the linear term coefficient a increases monotonically.Figure 9
The relationship between the linear coefficienta and the stress for specimens with different dip angles. (a) SS0°, (b) SS15°, (c) SS30°, (d) SS45°, (e) SS60°, (f) SS75°, and (g) SS90°.
(a)(b)(c)(d)(e)(f)(g)Figure10 shows the relationship between the nonlinear coefficient β and the maximum principal stress σ1 for samples with different dip angles. It can be seen that only the nonlinear coefficient β of the samples with dip angles of 75° and 90° is positive when the stress has been initially loaded and the nonlinear coefficient β of other samples is negative. This phenomenon has also appeared in many scholars’ seepage tests; Liao et al. studied and discussed the possibility and causes of negative values of β [27, 30–32]. In this test, the author believes that the negative values of the nonlinear coefficients β reflect the test attributes of the similar material samples of the columnar jointed rock mass. There are two potential reasons: (1) A boundary layer exists on the surface of porous media. As the pressure gradient increases, the thickness of the boundary layer decreases, and the cross-sectional area of the fluid and the permeability increase, resulting in a negative value of β. (2) White cement was used as the filling medium between the columns of the samples. As the pressure gradient increases, the maximum stress the filling medium can bear is exceeded, resulting in the formation of new fractures or the connection of the original fractures in the medium. The pore structure of the medium is changed, increasing the porosity and permeability and resulting in a negative value of β. However, an in-depth investigation is required to reveal the underlying mechanism.Figure 10
The relationship between the nonlinear coefficientβ and the maximum principal stress σ1.


### 3.1.2. The Relationship between the Intrinsic Permeability and the Maximum Principal Stress

Figure11 shows the relationship between the intrinsic permeability and the maximum principal stress of the samples with different dip angles. The samples of the columnar jointed rock mass show obvious seepage anisotropy, and there are significant differences in the relationship between the natural permeability k and the maximum principal stress σ1 for samples with different dip angles. The relationship can be summarized into two types: In the Type I curves (α = 0°–60°), the intrinsic permeability decreases initially and then increases with an increase in the maximum principal stress, showing a U-shape. Define the stress when the permeability decreases and start to increase as the inflection point stress, and the inflection point stress is 5, 3, 4, 3, and 4 MPa for the samples with different dip angles. In the Type II curves (α = 75°–90°), the intrinsic permeability exhibits a monotonous decrease with an increase in the maximum principal stress, and the curve is L-shaped.Figure 11
The relationship between the intrinsic permeabilityk and the maximum principal stress σ1 for samples with different dip angles. (a) SS0°, (b) SS15°, (c) SS30°, (d) SS45°, (e) SS60°, (f) SS75°, and (g) SS90°.
(a)(b)(c)(d)(e)(f)(g)According to the variation characteristics of the permeability-stress curve, try to fit the test data with one-variable cubic polynomial function:(3)k=a+bσ1+cσ12+dσ13,where k is the intrinsic permeability of the rock; σ1 is the maximum principal stress; a, b, c, and d are the fitting parameters. The R2 values of all fitting curves are greater than 0.96, indicating an excellent goodness-of-fit of the cubic polynomial function describing the relationship between the intrinsic permeability and the maximum principal stress of similar material samples of a columnar jointed rock mass.In addition, the values of the intrinsic permeability of the samples with different dip angles also show obvious anisotropy. The initial permeabilityks, the minimum permeability kmin, the maximum permeability kmax, and the final permeability kf of the samples with different dip angles are plotted in Figure 12.Figure 12
Theks, kmin, kmax, and kf of samples with different dip angles.The maximum decrease in the minimum permeability relative to the samples’ initial permeability is 63.51% for the specimen with a 75° dip angle, and the minimum decrease is 3.18% for the specimen with a 45° dip angle, representing a 20-fold difference. The maximum permeability of the 45° angle specimen increases by 144.25% compared with the initial permeability. The permeability of the 75° and 90° dip angle samples decreases monotonically, and the maximum permeability is the same as the initial permeability. The maximum permeability of the 0° and 60° dip angle samples is basically equal to the initial permeability. The final permeability increases by 0.43%, 60.90%, 42.28%, 144.25%, and 8.11% for the 0°–60° dip angle samples, whereas decreases of 63.51% and 58.96% are observed for the 75° and 90° dip angle samples, respectively (compared with the initial permeability). The sensitivity of the permeability of the samples with different dip angles to the maximum principal stress isα = 45° > α = 75° > α = 15° > α = 90° > α = 30° > α = 60° > α = 0°.


## 3.1.1. The Relationship between Volume Velocity and Pressure Gradient

Figure8 shows the relationship between the volume velocity Q and the pressure gradient dP/dL of the similar material samples with different inclination angles under maximum principal stress loading.Figure 8
The relationship between the volume velocityQ and the pressure gradient dP/dL for samples with different dip angles. (a) SS0°, (b) SS15°, (c) SS30°, (d) SS45°, (e) SS60°, (f) SS75°, and (g) SS90°.
(a)(b)(c)(d)(e)(f)(g)Significant differences are observed in the relationship between the volume velocity and the pressure gradient for the samples with different dip angles. Atα = 0°–60°, the fitting curves show nonlinear characteristics, and the degree of nonlinearity is different. At the same stress, the flow rate increases with an increase in the pressure gradient, indicating an increase in the permeability of the sample. In contrast, at α = 75°–90°, linear Darcy flow characteristics are observed.For the nonlinear flow relationship in fractures and porous media, Forchheimer proposed a zero-intercept quadratic equation to describe the nonlinear flow process from a macroscopic perspective [28, 29]: (1)−∇p=aQ+bQ2,(2)a=μkA,b=βρA2,where ▽p = dP/dL is the ratio of the pressure difference P between the inlet and outlet and the length L of the sample; Q is the overall volume velocity at the outlet; a and b are the model fitting coefficients, representing the specific gravity of the water pressure drop caused by the linear term and nonlinear term in the seepage test, respectively. When b = 0, the Forchheimer function degenerates to Darcy’s law; μ is the dynamic viscosity coefficient of the fluid; k is the intrinsic permeability of rock; A is the flow area of the rock; β is the nonlinear coefficient; ρ is the density of the fluid.Table1 lists the fitting coefficients a and b, the coefficients of determination R2, the values of the intrinsic permeability k, and the nonlinear coefficient β of the fractured rock obtained from equation (2). The R2 values of all test conditions are greater than 0.98, indicating a good agreement between the experimental and simulated values. Thus, the Forchheimer equation is suitable to describe the relationship between the volume velocity and pressure gradient of similar material samples of the columnar jointed rock mass.Table 1
The fitting results of the nonlinear Forchheimer equation and the calculated values ofk and β
Sample numberMaximum principal stressσ1 (MPa)Fitting coefficienta (Pa·s·m−4 × 1014)Fitting coefficientb (Pa·s2·m−7 × 1022)Coefficient of determinationR2Intrinsic permeabilityk (m2)Nonlinear coefficientβ (m−1)SS0°18.078E − 01−1.361E − 020.99931.579E − 15−6.765E + 1228.821E − 01−1.270E − 020.99851.446E − 15−6.312E + 1239.157E − 01−1.049E − 020.99911.393E − 15−5.214E + 1249.852E − 01−1.543E − 020.99911.294E − 15−7.669E + 1259.792E − 01−2.216E − 020.99891.302E − 15−1.101E + 1368.043E − 01−1.372E − 020.99761.585E − 15−6.819E + 12SS15°17.929E − 01−1.533E − 020.98831.608E − 15−7.619E + 1229.268E − 01−2.202E − 020.98771.376E − 15−1.094E + 1338.912E − 01−2.152E − 020.99331.431E − 15−1.070E + 1348.455E − 01−1.856E − 020.99401.508E − 15−9.225E + 1257.246E − 01−1.476E − 020.99461.760E − 15−7.336E + 1264.928E − 01−6.380E − 030.99352.588E − 15−3.171E + 12SS30°18.718E − 01−1.875E − 020.99821.463E − 15−9.319E + 1221.024E + 00−2.772E − 020.99901.246E − 15−1.378E + 1331.113E + 00−3.456E − 020.99601.146E − 15−1.718E + 1341.027E + 00−3.233E − 020.99631.242E − 15−1.607E + 1358.599E − 01−2.280E − 020.99841.483E − 15−1.133E + 1366.127E − 01−1.117E − 020.99872.081E − 15−5.552E + 12SS45°19.247E − 01−2.320E − 020.99331.379E − 15−1.153E + 1329.551E − 01−2.243E − 020.99521.335E − 15−1.115E + 1338.284E − 01−1.956E − 020.99731.539E − 15−9.722E + 1245.011E − 01−5.140E − 030.99902.545E − 15−2.555E + 1254.406E − 01−3.710E − 030.99962.894E − 15−1.844E + 1263.786E − 01−2.860E − 030.99833.368E − 15−1.421E + 12SS60°18.4222E − 01−1.8600E − 020.99651.514E − 15−9.245E + 1221.0840E + 00−2.3160E − 020.99791.176E − 15−1.151E + 1331.1346E + 00−1.9720E − 020.99611.124E − 15−9.801E + 1241.0895E + 00−2.1190E − 020.99631.170E − 15−1.053E + 1359.9456E − 01−2.4650E − 020.99911.282E − 15−1.225E + 1367.7904E − 01−1.3570E − 020.99921.637E − 15−6.745E + 12SS75°15.3832E − 011.9100E − 030.98672.369E − 159.493E + 1126.5698E − 012.1100E − 030.99461.941E − 151.049E + 1237.7073E − 01−1.0400E − 030.99501.655E − 15−5.169E + 1148.9612E − 01−6.2300E − 030.99161.423E − 15−3.096E + 1251.0872E + 00−9.0500E − 030.99631.173E − 15−4.498E + 1261.4753E + 00−2.6690E − 020.99358.643E − 16−1.327E + 13SS90°15.4456E − 011.6000E − 030.99072.342E − 157.952E + 1127.4407E − 01−6.7000E − 030.99621.714E − 15−3.330E + 1238.9512E − 01−8.2300E − 030.99751.425E − 15−4.091E + 1249.2963E − 01−1.6400E − 030.99751.372E − 15−8.151E + 1151.1876E + 00−2.8810E − 020.97791.074E − 15−1.432E + 1361.3268E + 00−3.7120E − 020.98029.611E − 16−1.845E + 13Table1 indicates that the variation trend of the linear term coefficient a in the fitting equation of the samples with different dip angles with an increase in the maximum principal stress σ1 is different. The trends are shown in Figure 9. It can be seen that at α = 0°–60°, the linear term coefficient a increases first and then decreases, at α = 75°–90°, and the linear term coefficient a increases monotonically.Figure 9
The relationship between the linear coefficienta and the stress for specimens with different dip angles. (a) SS0°, (b) SS15°, (c) SS30°, (d) SS45°, (e) SS60°, (f) SS75°, and (g) SS90°.
(a)(b)(c)(d)(e)(f)(g)Figure10 shows the relationship between the nonlinear coefficient β and the maximum principal stress σ1 for samples with different dip angles. It can be seen that only the nonlinear coefficient β of the samples with dip angles of 75° and 90° is positive when the stress has been initially loaded and the nonlinear coefficient β of other samples is negative. This phenomenon has also appeared in many scholars’ seepage tests; Liao et al. studied and discussed the possibility and causes of negative values of β [27, 30–32]. In this test, the author believes that the negative values of the nonlinear coefficients β reflect the test attributes of the similar material samples of the columnar jointed rock mass. There are two potential reasons: (1) A boundary layer exists on the surface of porous media. As the pressure gradient increases, the thickness of the boundary layer decreases, and the cross-sectional area of the fluid and the permeability increase, resulting in a negative value of β. (2) White cement was used as the filling medium between the columns of the samples. As the pressure gradient increases, the maximum stress the filling medium can bear is exceeded, resulting in the formation of new fractures or the connection of the original fractures in the medium. The pore structure of the medium is changed, increasing the porosity and permeability and resulting in a negative value of β. However, an in-depth investigation is required to reveal the underlying mechanism.Figure 10
The relationship between the nonlinear coefficientβ and the maximum principal stress σ1.


## 3.1.2. The Relationship between the Intrinsic Permeability and the Maximum Principal Stress

Figure11 shows the relationship between the intrinsic permeability and the maximum principal stress of the samples with different dip angles. The samples of the columnar jointed rock mass show obvious seepage anisotropy, and there are significant differences in the relationship between the natural permeability k and the maximum principal stress σ1 for samples with different dip angles. The relationship can be summarized into two types: In the Type I curves (α = 0°–60°), the intrinsic permeability decreases initially and then increases with an increase in the maximum principal stress, showing a U-shape. Define the stress when the permeability decreases and start to increase as the inflection point stress, and the inflection point stress is 5, 3, 4, 3, and 4 MPa for the samples with different dip angles. In the Type II curves (α = 75°–90°), the intrinsic permeability exhibits a monotonous decrease with an increase in the maximum principal stress, and the curve is L-shaped.Figure 11
The relationship between the intrinsic permeabilityk and the maximum principal stress σ1 for samples with different dip angles. (a) SS0°, (b) SS15°, (c) SS30°, (d) SS45°, (e) SS60°, (f) SS75°, and (g) SS90°.
(a)(b)(c)(d)(e)(f)(g)According to the variation characteristics of the permeability-stress curve, try to fit the test data with one-variable cubic polynomial function:(3)k=a+bσ1+cσ12+dσ13,where k is the intrinsic permeability of the rock; σ1 is the maximum principal stress; a, b, c, and d are the fitting parameters. The R2 values of all fitting curves are greater than 0.96, indicating an excellent goodness-of-fit of the cubic polynomial function describing the relationship between the intrinsic permeability and the maximum principal stress of similar material samples of a columnar jointed rock mass.In addition, the values of the intrinsic permeability of the samples with different dip angles also show obvious anisotropy. The initial permeabilityks, the minimum permeability kmin, the maximum permeability kmax, and the final permeability kf of the samples with different dip angles are plotted in Figure 12.Figure 12
Theks, kmin, kmax, and kf of samples with different dip angles.The maximum decrease in the minimum permeability relative to the samples’ initial permeability is 63.51% for the specimen with a 75° dip angle, and the minimum decrease is 3.18% for the specimen with a 45° dip angle, representing a 20-fold difference. The maximum permeability of the 45° angle specimen increases by 144.25% compared with the initial permeability. The permeability of the 75° and 90° dip angle samples decreases monotonically, and the maximum permeability is the same as the initial permeability. The maximum permeability of the 0° and 60° dip angle samples is basically equal to the initial permeability. The final permeability increases by 0.43%, 60.90%, 42.28%, 144.25%, and 8.11% for the 0°–60° dip angle samples, whereas decreases of 63.51% and 58.96% are observed for the 75° and 90° dip angle samples, respectively (compared with the initial permeability). The sensitivity of the permeability of the samples with different dip angles to the maximum principal stress isα = 45° > α = 75° > α = 15° > α = 90° > α = 30° > α = 60° > α = 0°.


## 3.2. Analysis of the Seepage Characteristics Using the Izbash Equation

The Izbash equation, which is a power function, has also been widely used to characterize the nonlinear seepage characteristics of fluids [33]:(4)−∇p=λQm,where λ and m are empirical coefficients. When m = 1, the Izbash equation degenerates to the linear Darcy’s law. When m = 1-2, equation (4) represents the nonlinear seepage characteristics caused by a significant inertia effect. Values of m = 0-1 are often used to characterize nonlinear seepage in low-permeability rock media caused by a significant solid-liquid interface effect.The Izbash equation is used to fit the test data. The fitting coefficientsλ and m and the R2 values for different maximum principal stress are listed in Table 2. The results show that the Izbash equation is suitable for characterizing the nonlinear seepage characteristics of similar material samples of the columnar jointed rock mass well. The R2 values are all greater than 0.98. The order of magnitude of the fitting coefficient λ is 1013–1014, and the range of the fitting coefficient m is between 0.68 and 1.02.Table 2
The fitting parameters of the Izabsh equation.
Sample numberMaximum principal stressσ1 (MPa)Fitting coefficientλFitting coefficientmCoefficient of determinationR2SS0°19.326E + 130.850.997529.744E + 130.890.999239.798E + 130.920.999541.072E + 140.890.999551.140E + 140.820.998469.529E + 130.840.9993SS15°19.895E + 130.800.989221.142E + 140.790.989831.139E + 140.770.991441.067E + 140.780.994551.006E + 140.750.996667.177E + 130.770.9988SS30°11.043E + 140.810.993021.229E + 140.780.992831.353E + 140.760.988841.395E + 140.690.983351.266E + 140.680.984769.560E + 130.710.9921SS45°11.140E + 140.780.986421.117E + 140.820.990031.100E + 140.750.990846.255E + 130.850.999255.451E + 130.870.999364.901E + 130.860.9993SS60°11.069E + 140.780.995421.202E + 140.860.998331.226E + 140.890.997441.207E + 140.870.998351.188E + 140.800.998969.357E + 130.830.9989SS7515.369E + 131.020.986326.519E + 131.020.994437.757E + 130.990.995049.339E + 130.950.991951.100E + 140.960.995961.496E + 140.930.9928SS90°15.618E + 131.000.990328.247E + 130.910.997839.540E + 130.930.998249.515E + 130.980.997751.364E + 140.830.985761.497E + 140.830.9877The trends of the coefficientsλ with the maximum principal stress σ1 obtained from the Izbash equation for samples with different dip angles are plotted in Figure 13. It can be seen that the trend of the coefficients λ is basically consistent with that of the linear coefficient a obtained from the Forchheimer equation; there is a good corresponding relationship between the two. In addition, there is also a corresponding relationship between the two equations in describing the nonlinear seepage characteristics of similar material samples of columnar jointed rock mass; that is, the fitting coefficient m < 1 of the Izbash equation corresponds to the nonlinear coefficient β < 0 in the Forchheimer equation.Figure 13
The relationship between the fitting coefficientλ of the Izbash equation and the stress for samples with different dip angles. (a) SS0°, (b) SS15°, (c) SS30°, (d) SS45°, (e) SS60°, (f) SS75°, and (g) SS90°.
(a)(b)(c)(d)(e)(f)(g)Figure14 shows the R2 values of the Forchheimer equation and Izbash equations for fitting the test data. Both the zero-intercept quadratic function (Forchheimer equation) and the power function (Izbash equation) well describe the nonlinear flow of the similar material samples of the columnar jointed rock mass, and the difference between them is not significant. However, the Forchheimer equation has a clear theoretical basis and can be simplified by the Navier-Stokes equation if the flow channel does not change significantly. In contrast, the Izbash equation is purely empirical without clear physical meaning.Figure 14
TheR2 values of the Forchheimer and Izbash equations under all working conditions.


## 4. Conclusions

In this paper, the seepage characteristics of a columnar jointed rock mass under maximum principal stress were investigated using similar material samples with different dip angles and conducting seepage-stress coupling tests. The main conclusions are as follows:(1)
Under maximum principal stress, the seepage characteristics of the similar material samples of the columnar jointed rock mass showed anisotropic characteristics. The permeability of samples with different dip angles has different sensitivity to stress, and the law and degree of variation with stress are also different.(2)
The relationship between the volume velocityQ and the pressure gradient dP/dL of the samples with different dip angles showed nonlinear seepage characteristics. The Forchheimer equation accurately described the relationship, and the R2 values were greater than 0.98. In most cases, the nonlinear coefficient β was negative, which was attributed to the test properties and the unique geological structure of the similar material samples.(3)
The relationship between the intrinsic permeabilityk and the maximum principal stress σ1 of the samples with different dip angles was categorized into two types. In the Type I curves, the permeability initially decreased and then increased with an increase in the maximum principal stress, showing a U-shape. In the Type II curves, the permeability showed a monotonous decrease with an increase in the maximum principal stress, and the curve was L-shaped. The relationship can be well described by the cubic polynomial function equation k=a+bσ1+cσ12+dσ13, and the R2 values are greater than 0.96.(4)
The final permeability of the samples with 0°–60° dip angles increased compared with the initial permeability, and that of the 45° dip angle sample increased by 144.25%, exhibiting the highest sensitivity to stress. In contrast, the decreases were 63.51% and 58.96% for the 75° and 90° dip angle samples, respectively.(5)
The Izbash equation, a power function, also accurately described the nonlinear flow characteristics of the similar material samples of the columnar jointed rock mass, withR2 values greater than 0.98. The fitting coefficients λ and m were in good correspondence with the linear coefficient a and the nonlinear coefficient β in the Forchheimer equation.



---

*Source: 6661741-2021-01-28.xml*