Effect of shearing on hydraulic properties of rough‐walled fractures under different boundary conditions

This study experimentally investigated nonlinear flow characteristics of rough‐walled fractures during shear process under different boundary conditions, using the servo‐controlled shear‐flow apparatus. A series of shear‐flow tests were conducted on two types of rough fractures with different shear displacements under various boundary conditions. The effects of shear process, fracture surface roughness, boundary normal stiffness, and initial normal stress on nonlinear flow behaviors are analyzed. The results show that Forchheimer equation provides a good description of the nonlinear relationship between flow rate and pressure gradient in rough‐walled fractures. The linear and nonlinear coefficients decrease by 1‐2 and 2‐3 orders of magnitude during shear, respectively. Both the mechanical aperture and the hydraulic aperture increase with increasing shear displacement. The mechanical aperture is always larger than the hydraulic aperture, and their difference gradually increases with the increment of shear displacement. The contact area ratio exhibits a fast decrease as shear displacement ranges between 1 and 3 mm and then decreases gradually with the decreasing rate generally diminishes. With the increase in shear displacement, the critical Reynolds number shows an increasing trend. The critical Reynolds number is in range 3.22‐35.41 for fracture G1 and 1.23‐12.35 for fracture G3. As the boundary normal stiffness and initial normal stress increase, the critical Reynolds number decreases.


| INTRODUCTION
Understanding the hydraulic properties of rock fractures is of great importance in many geophysical processes and underground engineering applications, such as oil natural gas production, geothermal energy development, nuclear waste disposal, and CO 2 sequestration. [1][2][3][4][5] The rock masses commonly comprise matrix blocks encased with numerous individual fractures. The fractures provide the primary pathway of fluid flow due to the permeability of fracture is much larger than that of rock matrix. 6,7 The fracture permeability is sensitive to the stress conditions. Many previous studies have estimated the effects of normal stress and shear stress/ displacement on hydraulic properties of a single fracture. The results show that the permeability of a single fracture decreases due to the normal stress-induced fracture closure, [8][9][10][11][12][13] and many of empirical equations have been proposed to describe the relationships between normal stress and permeability. 14,15 When fracture undergoes shear, the permeability may decrease due to shear contraction or increase as a result of significant dilation. [16][17][18] In most of the previous studies, the fluid flow in a single fracture is usually described by classical cubic law. It assumes that fracture surface can be approximated by the smooth parallel plate model, and the flow rate is linearly proportional to the pressure drop and the cube of fracture aperture. 17,[19][20][21] The cubic law is valid for laminar flow with a sufficiently low Reynolds number (Re). However, the inertial effects are not negligible as the Re increases to a relatively large value. The rough-walled fracture exhibits a complexity of the geometrical characteristics on their surfaces, which gives rise to channeling flow and/or eddies resulting from the variations in fracture aperture. In this situation, the flow deviates from a linear flow regime, and using cubic law to calculate fluid flow overestimates the transmissivity of rock fractures. [22][23][24][25] Substantial efforts have been focusing on nonlinear flow characteristics in rough-walled fractures. Two reasons trigger the inertial losses during the fluid flow through the fractures. One reason is due to eddy flow as the Re increases. [26][27][28] However, many studies observed that the nonlinear flow generally occurs before the fully developed turbulence. They found that weak inertia regimes exist even when Re has a small value. 11,25,29,30 This phenomenon is mainly attributed to the roughness, contact areas, and tortuosity-caused non-negligible inertial losses. 18,31,32 To evaluate the nonlinear flow characteristics of rock fractures, most early studies have focused on the influences of normal load and surface roughness. The results show that nonlinearity in flow through fractures shows a decreasing trend with increasing fracture roughness [33][34][35][36][37] and normal loads. 11,25,30,38,39 However, these studies are focused on the hydraulic behaviors of mated fractures. The effect of shear stress/displacement on nonlinear flow characteristics was not taken into account. Rock fractures are often subject to shear resulting from various natural and human activities. [40][41][42][43][44][45][46] During fracture shearing, the shear-induced dilation and mismatching of fracture surfaces cause contact spot changes and aperture distributions, which influence the flow and transport behaviors in rough-walled fractures. Therefore, studying the effect of shear process on nonlinear flow behavior in rough-walled fractures is highly important. Recently, shear impacts on flow regimes and nonlinear flow behavior have been investigated by numerical methods and experimental investigations. The numerical simulation studies proposed that shear displacement can enhance the degree of flow nonlinearity and induce nonlinear laminar flow in rough-walled fractures. 47,48 However, the shear step is created by artificially moving one half of wall with respect to the other wall, which was not received from physically realistic laboratory test. This leads to that the opposite response was observed for the experimental investigations. The experimental efforts observed that flow nonlinearity decreases with increasing shear displacement. 23,49,50 However, these shear-flow tests were conducted under constant normal load (CNL) boundary conditions, in which the normal load applied on the specimen was constant during the shear process. The CNL condition is only appropriate for the nonreinforced rock slope or planar fractures with no dilation during shearing. [51][52][53][54] For deep underground scenarios, the normal load applied perpendicular to the shear direction is not a constant. The rock joints are subjected to a variable normal load and the stiffness of the surrounding rock mass could significantly affect the shear behavior. Thus, the constant normal stiffness (CNS) condition is more applicable than the CNL conditions. 53,55 Additionally, the evolution of fracture aperture under CNS conditions is more complicated than that under CNL conditions. Two competitive effects are involved under CNS conditions during shearing. On the one hand, shear-induced dilation increases fracture aperture. On the other hand, shear dilation-induced increase in normal stress results in asperity degradation, which restrains the dilation of rock fracture. This decreases dilation during shear process, and the variation in hydraulic behavior was not well understood. Therefore, the effect of shear process on nonlinear flow characteristics under CNS conditions should be studied, which has not been reported before, if any.
In this study, we experimentally investigated the effect of shear process on nonlinear flow behaviors of rough-walled fractures under different boundary conditions, using a servo-controlled shear-flow testing apparatus. The shear-flow tests were conducted on plaster specimens with two different rough fractures. For each fracture specimen, high precision hydraulic tests were performed at different shear displacements with various pressure gradients. Finally, the shear impacts on nonlinear flow characteristics through rough-walled fractures under different boundary conditions, as well as the mechanical and hydraulic aperture, contact area ratio, normalized transmissivity, and critical Reynolds number, were analyzed.

| THEORETICAL BACKGROUND
For viscous flow, the fluid flow in a single fracture can be simplified by the classical cubic law, which can be written as 56-59 : where Q is the flow rate, w is the width of a fracture, e h is the hydraulic aperture, ΔP is the pressure drop in the flow direction, and L is the fracture length over the pressure drop.
Based on the cubic law, the transmissivity (T) of the fracture can be written as follows 17 : However, the cubic law is only anticipated only for describing the fluid flow behavior in rock fractures under a sufficiently low Re. As the Re increases to a relatively large value, the fluid flow enters the nonlinear flow zone in which the increase in discharge is smaller than the proportional increase in pressure gradient (∇P). To describe the nonlinear relation between ∇P and Q, Forchheimer equation is the most extensively adopted method, [11][12][13]27,29,30 in which the ∇P is a quadratic function of the Q, as follows: where a and b are coefficients that represent the pressure drop caused by linear and nonlinear effects, respectively. is the non-Darcy coefficient.
The Re is defined as the ratio of inertial force to viscous force, which is commonly used to predict the termination of linear flow and onset of nonlinear flow behavior, written as 25,27,38 To quantitatively evaluate nonlinear flow effect in roughwalled fractures, a factor E was introduced to determine the fluid flow regime 60 : For engineering purposes, the critical condition for the flow regimes from linear to nonlinear has been defined at E = 0.1. 11,29 Therefore, when the nonlinear term contributes to 10% of the total pressure drop, the corresponding Re is defined as the Re c 23,25,35 : The mechanical aperture, e m , is defined as the mean pointto-point between the walls of a fracture. During the shear process, the e m can be calculated based on the following equation 16 : where e 0 is the initial mechanical aperture, Δe n is the change in mechanical aperture by normal stress, and the Δe s is the change in mechanical aperture by shearing. The e 0 under a certain normal stress can be calculated by using a normal loading-unloading test. 8 In some circumstances, the e 0 can be assumed to equal to the initial hydraulic aperture to keep a constant initial status. 22

| Specimen preparation and roughness measurement
Two granite blocks were prepared from the field of an underground power station in Japan. The blocks were cut and polished to cuboid specimen with a dimension of length × width × height = 200 × 100 × 100 mm. Then, an artificial tensile fracture was created using the Brazilian test. Two artificially created rough-walled fractures (Labelled G1 and G3) with different surface morphologies are shown in Figure 1. To ensure the specimens under different test conditions having the same roughness, plaster replicas were used to conduct shear-flow tests under various boundary stress conditions. The replicas were made of mixtures of high-strength plaster, water, and retardant with a weight ratio of 1:0.2:0.005. 52 A total of 8 specimens were prepared, which were maintained at a constant temperature of 25°C and placed in a humid box with a relative humidity of 95% for 28 days after the specimens were demoulded. The physic-mechanical properties of these rock-like specimens are referred to Li et al. 18 The surfaces of the granite fractures were measured using a 3D laser scanning profilometer system, which has an accuracy of ±20 μm in X-Y plane and ±10 μm in the Z axis. The joint roughness coefficient (JRC) of each surface is quantitatively described using the method proposed by Tse and Cruden 61 : where Z 2 is the root mean square slope of the profiles based on the extracted data, z i represents the coordinates of the fracture surface profile, n is the number of data points, and Δx is the interval between the data points. The mean JRC values of fractures G1 and G3 are 3.21 and 7.36, which are evaluated by cutting the surface using 10 profiles parallel to the shear direction.
The variation in frequency vs asperity height of two fracture surface is showed in Figure 1C,D. The results show that the asperity height follows a Gaussian distribution, which is consistent with many observations reported by Shapiro and Wilk, 62 Sharifzadeh et al, 63 and Adler et al. 64 Compared with the fracture G1, the rougher fracture G3 trends to generate a wide range of the asperity distribution, and the mean asperity height and standard deviation are larger.

| Experimental apparatus and procedure
The shear-flow tests under various boundary conditions were conducted using the servo-controlled shear-flow test apparatus, as shown in Figure 2. It mainly consists of a hydraulic-servo actuator unit, a hydraulic testing unit, and a visualization unit. A horizontal loading system and a stiffness controlling in vertical direction are used to apply the shear and normal loads with a maximum force of 100 kN. The shear displacements are measured by one linear variable differential transformer (LVDT) with an accuracy of 0.001 mm, which is attached on the load jacks in shear direction. Two LVDTs are fixed on the two sides of upper part of shear box to measure the normal displacements. The capacity of LVDTs in shear (horizontal) and normal (vertical) directions is 20 mm and 10 mm, respectively. The constant inlet water pressure is applied on the one side of the fracture, and the pressure is supplied by an air pump. The lateral sides of the fracture were sealed using soft and elastic gel sheets, which were very flexible and strong with a good sealing effect. The detailed description of the apparatus was reported by Li et al. 18 One of the important features of this apparatus is that it can automatically reproduce CNL and CNS boundary conditions with a good precision. The CNS boundary condition is achieved by a closed loop in the system control program, with electrical and hydraulic-servo controls. We employed a nonlinear feedback of control and measurement on PC windows through a multifunctional board. The digital control program is designed using LabVIEW programming language for building data acquisition and instrumentation systems.
Here, the change in the normal stress due to the application of normal stiffness (k n ) was calculated as follows: where Δ n and Δ n are changes in normal load and normal displacement, respectively. The effect of shear process on nonlinear flow behavior of rough-walled fractures under CNL and CNS conditions was investigated. First, the well mated two halves of plaster specimens were placed in the shear box. The lower part of the specimen was fixed, and the upper part could move vertically and horizontal without rotation during shear. Then, the CNL/ CNS condition was applied in the shear box through a normal load jack. Before the hydraulic tests, the inlet and outlet of the specimen were sealed using rubber sheets and two lateral sides were sealed using soft and elastic gel sheets. This is a simple and effective sealing method for coupled shear-flow tests. Many coupled shear-flow tests have been performed using this apparatus, which shows a good sealing effect (ie, 18,22,65,66 ). The flow tests were carried out at six different shear displacements (d) of 1, 2, 3, 5, 7, and 10 mm. The shear velocity was 0.5 mm/min. At each d, we performed 8-10 constant water head flow tests. The water head difference ranges from 0.1 to 1.0 m. The pressure drop was measured by a differential pressure transducer with a rated capacity of 10 kPa and a precision of 0.01 kPa. The water flowing out of the fracture was collected in a glass container, and the weight of the container was recorded using an electrical balance that has a precision of 0.01 g. In all hydraulic tests, the flow direction is parallel to the shear direction. The test cases and their corresponding boundary conditions are listed in Table 1. The entire flow test was conducted at a room temperature of about 20°C. The density and the dynamic viscosity of water are ρ = 0.998 × 10 3 kg/m 3 and μ = 1 × 10 −3 Pa s, respectively.

| Mechanical behavior
The tested shear behaviors of the rough fractures under various boundary conditions are presented in Figure 3 and listed in Table 2. As shown in Figure 3A,B, the curve of shear stress vs shear displacement under CNL conditions exhibits a similar shear mode. The shear stress rapidly increases at the beginning of shear until reaching the peak shear strength (τ p ), and then, stress-softening behavior is displayed in the postpeak stage. Finally, the shear stress becomes stably at the residual shear strength. For CNS conditions, the shear behavior is dependent on the JRC, initial normal stress (σ n0 ), and normal stiffness. At the prepeak stage, both the τ p and the d p increase with the increment of σ n0 . For CNS conditions at the same σ n0 , the prepeak stage of shear behavior is similar to CNL conditions. This means the k n in prepeak stage is considered to have no direct effect on the τ p . Similar conclusions are also reported by Lee et al. 54 Additionally, it is found that τ p is related to the fracture surface roughness. The rougher fracture shows a larger of τ p . At the postpeak stage, the shear stress increases with increasing k n and σ n0 . For fracture G1, the shear stress slowly decreases after τ p , and the stress-softening behavior is less significant with increasing k n . However, for fracture G3, the stress-hardening behavior can be observed, and the phenomenon of stress hardening becomes more pronounced with increasing k n . This indicates that rougher fracture increases the tendency toward that stress-hardening behavior. Figure 3C,D shows the variations in normal displacement vs the shear displacement. For all test cases, the normal displacement exhibits nonlinear descending and ascending behavior with increasing d. As shearing begins,  shear displacement. As shown in Figure 3E,F, the normal stress is kept constant for CNL condition. However, the stress paths for CNS conditions are slightly decreased at first and then trended toward increasing normal stress as d increases. The normal stress shows an increasing trend with increasing k n and σ n0 , restraining the dilation of the fractures. The restrained dilation under high normal stress is mainly because of the significant asperity damage. Additionally, the evolutions of normal stress and dilation are related to the fracture surface roughness. The rougher the fracture surface, the larger dilation and normal stress that could be obtained.  Figure 4 shows the relationship between ∇P and Q for fractures G1 and G3 corresponding to different shear displacements subjected to various boundary conditions. The correction between ∇P and Q shows obvious nonlinear characteristics, which can be fitted by a zero intercept quadratic polynomial. The residual squared R 2 for all cases is larger than 0.99, indicating that the Forchheimer equation fits the flow data very well. As d increases, the slope of ∇P -Q curves becomes flatter as a result of shear dilation-induced fracture aperture increase. Similarly, compared with fracture G1, the slope of ∇P-Q curves of fracture G3 becomes more flatter at the same d, indicating a smaller flow resistance as a result of larger dilation for rough fractures. Hence, less energy is needed to achieve the same flow rate. However, for fracture specimen with a certain d, as σ n0 and k n increase, the slope of ∇P-Q curves becomes steeper, implying a larger flow resistance due to the restrained dilation of the fracture.

F I G U R E 4 Relationships
Based on Equation (3), the linear coefficient a and nonlinear coefficient b are calculated, as listed in Table 3. Both coefficients a and b in Forchheimer equation show decreasing trends with increasing d. The decrease extent at d from 1 to 3 mm is more significant than that in the range between 3 and 10 mm. The a and b show a reduction of 1-2 and 2-3 orders of magnitude as d increases from 1 to 10 mm, respectively. According to Equation (4), the decrease in a is mainly due to the shear-induced dilation of rock fracture, and the decrease in b is ascribed to the dilation of fracture or/and flow path tortuosity decreases as the contact area decreases. However, for a certain d, as σ n0 and k n increase, both a and b decrease. This is due to the higher normal stress is more prone to restrain the dilation of the fracture. Meanwhile, the production of gouge material under higher normal stress in CNS conditions may impede the fluid flow through the fracture.

| Evolutions of aperture and contact area ratio during shear process
During the shear process, the mismatching of fracture surfaces significantly changes the contact area and aperture distributions, which is an important issue affecting the hydraulic behavior in rough-walled fractures. Based on the measured topographical data and Equation (8), the evolutions of the mechanical aperture distributions under different shear displacements are obtained. Figure 5 shows the aperture distributions of fractures G1 and G3 for three different shear displacements under CNL condition. For the two fracture surfaces, at d = 1 mm, the fracture exists a large number of contact areas (aperture is equal to 0). The mean aperture and the range of aperture distributions are very small. As d increases, the contact area ratio (C) decreases, but the mean aperture and the aperture range increase significantly. For two fracture surfaces, the rougher fracture trends to generate a wider range of the aperture distribution. Figure 6 shows the evolutions of e m and e h during shear process of two rough-walled fractures under various boundary conditions. The e m is calculated based on Equation (8), and the e h is calculated by the cubic law at a constant water head of 0.1 m. As shown in Figure 6, both the e m and the e h increase with increasing d, due to the shear-induced dilation, the fracture aperture. However, the e m is always larger than the e h , and their difference gradually increases with the increment of d. The degree of e h deviating from e m is larger for a rougher fracture surface. These phenomena may be caused by the tortuosity of the streamlines due to the existence of contact areas and unevenly fracture aperture distributions. 22 Additionally, for a certain d, the e m and e h decrease gradually with increasing k n and σ n0 . The increases in normal stiffness and initial normal stress caused a more significant asperity damage, which decreases the dilation of fracture aperture. These results are consistent with the observations of Olsson and Barton, 17 Li et al, 18 and Koyama et al. 66 Figure 7 shows the evolutions of C during shear process for two rough-walled fractures under various boundary conditions. The C is defined as the percentage of the number   18 As shown in Figure 7, the C changes inversely with the change in e h during shear. As d increases from 1 to 3 mm, the C decreases significantly and then gradually decreases. When the d exceeds 5 mm, the C gradually reaches a constant value. For the same test condition at the same d, the C shows an increasing trend with increasing k n and σ n0 . This is mainly because the asperity deformation increases induced by the high normal stress.

| Critical conditions for nonlinear flow
To illustrate the nonlinear flow characteristics, some previous studies adopted normalized transmissivity (T/T 0 ) to reflect the flow capability in rough-walled fractures. The normalized transmissivity is defined as the ratio of the apparent transmissivity (T) to the intrinsic transmissivity (T 0 ), which can be described as a function of Re and Forchheimer coefficient , as follows 29 : where T 0 is determined using the best fitted values of coefficient a (ie, T 0 = /a). Figure 8 shows the relationship between T/T 0 and Re for rough-walled fractures under various boundary conditions. When Re is sufficiently small (ie, less than 1), T/T 0 decreases slowly and approximately equals to 1 due to the viscosity. As Re increases, the T/T 0 decreases and the reduction rate increases because of inertia effect. The results indicate that variations in T/T 0 are sensitive to the Re. Additionally, as d increases, the T/T 0 vs Re curves generally shift upward. The shift degree of T/T 0 increases more significantly for d from 1 mm to 3 mm than that from 3 mm to 10 mm. The small value of d can easily lead to a transition from linear to nonlinear flow.
Note that the T/T 0 = 0.9 have the same physical meaning with E = 0.1, both of which indicate that the nonlinear term contributes to 10% of the total pressure drop. The relationship between the T/T 0 and E can be calculated as follows 39,67 : Based on Equation (13), when T/T 0 = 0.9, Re c is calculated. Figure 9 shows the variation in Re c during shear process for rough-walled fractures under different boundary conditions. The results show that Re c increases significantly for d in the range of 1-3 mm, and then, the rate of increase gradually decreases. At the beginning of shear, the Re c is very small, in the range of 3.22-5.32 for fracture G1 and 1.23-2.76 for fracture G3. This is mainly because the slightly mismatched fracture surface formed a small void space with many contact areas, which make the flow paths more tortuous. As d increases from 1 to 3 mm, the Re c increases significantly due to the shear dilation and decrease in the contact areas. With continuously increasing d, the dilation of the fracture aperture decreases the influence of fracture roughness on nonlinear flow and results in a slight increase in Re c . In addition, the roughness of the fracture surface plays a dominate role in controlling the fluid flow, and the rougher fracture results in a more dramatic increase in the degree of nonlinearity of fluid flow and a smaller Re c . The Re c is in range between 3.22 and 35.41 for fracture G1 and range of 1.23-13.01 for fracture G3. This is mainly caused by more complexity of the flow velocity distributions or direction along the flow channel in rougher fractures. As shown in Figure 1C,D, the rougher fracture trends to generate a wider range of the asperity distribution, which indicates that the aperture distribution becomes more anisotropic and heterogeneous for rougher fractures during shear process Figure 5. When shear occurs, the sharp changes in fracture aperture along the flow path will promote the variations in the plane velocity, due to the acceleration and deceleration of flow to maintain the conservation of mass. This causes the fluid flow more easily to become nonlinear. Additionally, the Re c shows a decreasing (13)   trend with increasing k n and σ n0 . The possible reason could be that the fracture closure increases the C for high normal stress under CNS conditions Figure 7. Moreover, the production of gouge materials under high normal stress conditions may make the flow paths more tortuous and result in a lower Re c . To further analyze the influences of aperture and contact area on nonlinear flow behavior, the relationship between e h and C with Re c is plotted in Figure 10, in which the Re c increases with increasing e h . A power function is suitable for fitting Re c as a function of e h . The rougher fracture corresponds to a flatter slope. This indicates that the rougher fracture is more easily to enter nonlinear flow regimes under the same e h . As shown in Figure 10B, the Re c decreases with increasing C. The Re c decreases exponentially with increasing C. The relationship can be well described by a power function. The Re c shows a decreasing trend with increasing JRC. Their difference becomes small with increasing C.

| CONCLUSIONS
This study experimentally investigates the effect of shearing on nonlinear flow characteristics of fluid through the roughwalled fractures under different boundary conditions. We performed shear-flow tests on plaster replicas of two rough granite fractures with different shear displacements under both CNL and CNS conditions. At each shear displacement, a series of water flow tests were conducted with respect to different water pressures (1-10 kPa). The influences of shear, fracture surface roughness, and boundary stress conditions on nonlinear flow characteristics through rough-walled fracture, as well as the mechanical and hydraulic aperture, contact area ratio, normalized transmissivity, and critical Reynolds number, were analyzed.
The results show that Forchheimer equation provides a good description of the nonlinear flow process during shearing. The linear coefficient a and nonlinear coefficient b decrease by 1-2 and 2-3 orders of magnitude as shear displacement increases from 1 to 10 mm, respectively. Both the mechanical aperture and the hydraulic aperture increase with increasing shear displacement, due to the shear-induced dilation. The mechanical aperture is always larger than the hydraulic aperture, and their difference gradually increases with the increment of shear displacement. For a certain shear displacement, mechanical aperture and hydraulic aperture show increasing trends with increasing fracture surface roughness, but decrease with increasing normal stiffness and initial normal stress. As shear displacement increases from 1 to 3 mm, the contact area ratio decreases significantly and then gradually decreases. When the shear displacement exceeds 5 mm, the contact area ratio gradually reaches a constant value. The contact area ratio shows an increasing trend with increasing normal stiffness and initial normal stress due to the asperity deformation increases induced by the high normal stress. The normalized transmissivity decreases with increasing Reynolds number. As shear displacement increases, the normalized transmissivity curves generally shift upward. With the increase in shear displacement, the  In the present study, we characterized the evolutions of aperture and contact condition during shear process and analyzed the nonlinear flow behaviors of two rough-walled fractures under various boundary conditions. However, the mechanical aperture and contact area ratio are measured from the shear deformation data, which ignores the effects of plastic deformation and production of gouge materials on real aperture changes and contact areas. Therefore, new visualization techniques should be explored to quantitatively investigate the effects of aperture field, contact condition, and gouge material distribution on nonlinear flow behavior. Additionally, the fluid flow in rough-walled fracture is anisotropic. In the further, the works are needed to consider the shear-induced aperture distributions on anisotropic flow behavior in rough-walled fractures.